Development and Validation of a Combined Glycolysis and Immune Prognostic Model for Melanoma

Background: Glycolytic effects and immune microenvironments play important roles in the development of melanoma. However, reliable biomarkers for prognostic prediction of melanoma as based on glycolysis and immune status remain to be identified.

Methods: Glycolysis-related genes (GRGs) were obtained from the Molecular Signatures database and immune-related genes (IRGs) were downloaded from the ImmPort dataset. Prognostic GRGs and IRGs in the TCGA (The Cancer Genome Atlas) and GSE65904 datasets were identified. Least absolute shrinkage and selection operator (LASSO) Cox regression and multivariate Cox regression were used for model construction. Glycolysis expression profiles and the infiltration of immune cells were analyzed and compared. Finally, in vitro experiments were performed to assess the expression and function of these CIGI genes.

Results: Four prognostic glycolysis- and immune-related signatures (SEMA4D, IFITM1, KIF20A and GPR87) were identified for use in constructing a comprehensive glycolysis and immune (CIGI) model. CIGI proved to be a stable, predictive method as determined from different datasets and subgroups of patients and served as an independent prognostic factor for melanoma patients. In addition, patients in the high-CIGI group showed increased levels of glycolytic gene expressions and exhibited immune-suppressive features. Finally, SEMA4D and IFITM1 may function as tumor suppressor genes, while KIF20A and GPR87 may function as oncogenes in melanoma as revealed from results of in vitro experiments.

Conclusion: In this report we present our findings on the development and validation of a novel prognostic classifier for use in patients with melanoma as based on glycolysis and immune expression profiles.

Introduction

Skin cutaneous melanoma (SKCM) is a cancer resulting from malignant transformations of melanocytes in the basal layer of the skin. This condition represents one of the most aggressive types of skin cancer worldwide, with factors such as ultraviolet light (UV) and malignant transformations of moles being high-risk factors for SKCM (1). At present, most treatments for SKCM involve surgical resection removal, radiotherapy and chemotherapy (2). Although significant progress has been achieved with these treatments, more than 95% of metastatic melanoma patients die within one year (3). Therefore, there is an urgent need to identify prognostic biomarkers, which would provide clinicians the ability to promptly and accurately predict the clinical outcome of melanoma as well as initiate a protocol for a personalized treatment regimen.

Abnormal immune microenvironments and tumor metabolic reprogramming represent two significant characteristics of tumors (4). Tumor metabolic reprogramming is one of the critical mechanisms involved with regulating the immune microenvironment (5, 6). Tumor metabolic reprogramming includes enhanced aerobic glycolysis (Warburg effect), increased glucose uptake and consumption, enhanced lipid and protein synthesis and enhanced glutamine uptake and catabolism, with glycolysis as the main feature (7–9). Even with sufficient oxygen availability, tumor cells continue to metabolize glucose mainly through glycolysis to produce lactic acid (10, 11). The capacity for this glycolysis within tumor cells along with its production of lactic acid to participate in the regulation of immune microenvironments have become notable areas of research interest in recent years. It has been reported that lactic acid, as secreted from tumor cells, can inhibit the cytolytic ability of CD8+ effector T cells, but it may also be used by Treg cells to support cell metabolism (12). The lactic acid produced by melanoma cells can reduce the immune surveillance ability of T and NK cells by inhibiting the nuclear factor produced by the activated T cell (NFAT) dependent IFN-γ (13). In addition, the findings that glycolysis is the main metabolic pathway required after T cell activation, and enhances the eradicating effects of T cells (14, 15), indicates that glycolytic effects play different roles in immune versus tumor cells. Although a clear connection exists between glycolytic effects and immune microenvironments, scant attention has been directed toward examining this relationship in any detail.

In this study, we constructed and validated a comprehensive index of glycolysis and immune (CIGI) model as based on glycolysis- (GRGs) and immune- (IRGs) related genes. This CIGI model shows stable, prognostic prediction performance in different data sets and in different subgroups of melanoma patients. In additional, we demonstrate that CIGI was correlated with glycolysis and immune status. Finally, the expression and function of CIGI genes in melanoma were evaluated in in vitro experiments.

Materials and Method

Data Collection

Gene expression profile and clinical follow-up information of RNA sequencing samples from patients with SKCM were downloaded from the TCGA database, and the GSE65904 cohort was downloaded from the Gene Expression Omnibus (GEO) database. The RNA sequencing data from TCGA-SKCM were preprocessed as follows: 1) the samples without follow-up information were removed; 2) ENSEMBL ID were converted to gene symbols; 3) Remove the samples with OS <30 days; 4) Genes with expression level is lower than 1 and the proportion is higher than 50% in all samples were eliminated. The GEO cohort was processed via the following steps: 1) the samples without clinical follow-up information were removed; 2) ENSEMBL IDs were converted to gene symbols; 3) the probes corresponding to multiple genes were removed; 4) the median of multiple gene expression values was used. TCGA database and GSE65904 dataset were selected as they contained the largest sample set in the same platform with detailed follow-up information on melanoma. Glycolysis-related genes were identified in the Molecular Signatures database, while immune-related genes were downloaded from the ImmPort dataset.

Construction and Validation of the CIGI

The independent prognostic predictors among glycolysis- and immune-related genes were selected via univariate analysis based on the “survival” package. LASSO analysis and the stepwise Cox proportional hazards regression model were used to construct CIGI. High- and low-CIGI groups were differentiated as based on the optimized risk value. Kaplan-Meier survival analyses were used to analyze differences in overall survival between the high- and low-CIGI groups. Time‐dependent ROC curve analysis was used to evaluate the predictive value of the CIGI. Univariate and multivariate Cox regression analyses were performed to explore the independent prognostic value of the CIGI.

Potential Regulatory Pathways Analysis

Single-sample gene set enrichment analysis (ssGSEA) based on the “GSVA” package was used to quantify the scores of pathways in each sample. And, the ssGSEA score from each sample was used to analyze the relationship between CIGI and potentially regulatory pathways.

Immune Infiltration Analysis

The infiltration levels of 28 immune cell types in the TCGA dataset were assessed using ssGSEA analysis. The gene markers of 28 immune cells were obtained from a previous study (16).

Univariate and Multivariate Cox Analysis

Univariate and multivariate COX regressions were used to analyze the relationship between CIGI and other variables and the clinical prognosis of patients.

IHC (Immunohistochemistry) Analysis

Melanoma samples and matched nontumorous tissue were obtained from the First Hospital of China Medical University. Tissues were fixed in 10% formalin, embedded in paraffin, and processed as 4-µm continuous sections. IHC staining was performed according to the manufacturers’ instructions (UltraSensitiveTM SP; MXB, China). The antibodies used included: SEMA4D (1:200; ab134128; Abcam), IFITM1 (1:2000; ab233545; Abcam), KIF20A (1:100; ab70791; Abcam) and GPR87 (1:100; ab188901; Abcam). Each sample was independently assessed by two pathologists and scored using a semiquantitative scoring system with histoscores ranging from 0 (minimum) to 300 (maximum).

This, as well as any other portions of this study involving human tissue samples, were approved by the Ethics Committee of the First Hospital of China Medical University.

Cell Culture

The human melanocytes PIG1 cell line, melanoma cell line, A375, A875, and MeWo were obtained from the China Infrastructure of Cell Line Resource and was cultured in DMEM (10% FBS) at 37°C in a humidified 5% CO₂ incubator.

Cell Transfection

Small interfering RNA (siRNA) transfection was performed using Lipofectamine 3000 (Invitrogen, Shanghai, China). The sequences were as follows: 5’-GGCCTGAGGACCTTGCAGAAGA-3’ for SEMA4D-specific siRNA, 5’- CCTAGATACAGCAGTTTATAC-3’ for IFITM1-specific siRNA, 5’- GGCCAGGUUUCUGCCAAAATT-3’ for KIF20A-specific siRNA and 5’- UCUUAAUCGCGUAUAAGAGTT-3’ for GPR87-specific siRNA. The sequence for the negative control (NC) was 5′- UUCUCCGAACGUGUCACGUTT-3′.

CCK8 Assay

A375 cells (1500/well) were cultured in 96-well plates and transfected with NC-siRNA or siRNAs (SEMA4D-specific siRNA, IFITM1-specific siRNA, KIF20A-specific siRNA, GPR87-specific siRNA). After culture for 0, 24, 48, or 72 h, cells were cultured with the CCK8 solution (C0038, Beyotime, Shanghai, China) for an additional 2 h. Cell viability was expressed as an optical density (OD) value at 450 nm.

Colony-Forming Experiments

In order to examine the effects of SEMA4D, IFITM1, KIF20A and GPR87 expression on human melanoma cell proliferation, A375 cells (500/well) transfected with NC-siRNA or siRNAs were added to the 12-well plates. After two weeks, the number of colonies were counted.

Detection of Lactate, ATP and Glucose Uptake Levels

A375 cells (5×10⁵/well) were cultured in 96-well plates and transfected with NC-siRNA or siRNAs (KIF20A-specific siRNA, GPR87-specific siRNA). The culture medium and cells were collected after 48 h. Lactate levels in the medium were determined with use of the lactate assay kit (ab65331, Abcam), ATP levels with the ATP assay kit (ab83355, Abcam) and glucose uptake levels with the glucose uptake assay kit (ab136955, Abcam). All determinations were normalized with cell numbers.

Extracellular Acidification Rate (ECAR)

Measurement of ECAR was performed using the Seahorse XF96 Flux Analyser (Seahorse Bioscience) following instructions provided by the manufacturer. In short, after transfection with NC-siRNA or siRNAs (KIF20A-specific siRNA, GPR87-specific siRNA), A375 cells (1.5 × 10⁴/well) were cultured in 96-well plates. The culture medium was replaced with the test buffer prior to the ECAR test. Cells were then incubated with buffered medium containing 10mM glucose, 1mM oligomycin and 50mM 2-DG.

Statistical Analysis

The Kaplan-Meier analysis was used to compare differences in OS between the high- and low-CIGI groups. Statistical comparisons between two groups were performed using the Student’s two-tailed t-test. A p-value of < 0.05 was required for results to be considered as statistically significant.

Results

Development of CIGI in Melanoma

TCGA and GSE65904 datasets were included in this study. Univariate Cox regression analysis was conducted to select genes significantly related with prognosis based on these TCGA and GSE65904 datasets. A total of 6155 genes were significantly related to the overall survival of melanoma patients in the TCGA dataset, with hazard ratios of the top 10 genes(sorted by p value) presented in Figure 1A. A total of 2835 genes were significantly related to the overall survival of melanoma patients in the GSE65904 dataset, with hazard ratios of the top 10 genes(sorted by p value) presented in Figure 1B. A Venn diagram indicating that 182 prognostic GRGs and IRGs were identified in the GSE14520 and TCGA cohorts as shown in Figure 1C. As an excessively large number of genes is not conducive to clinical detection, we further narrowed the range of genes. We employed the Lasso regression analysis and the resulting change trajectory of each independent variable is shown in Figure 1D. With the gradual increase of lambda, the number of independent variable coefficients gradually increased to zero. Five-fold cross-validation was used to build the model, and the confidence interval under each lambda. Thus, we selected 17 genes at lambda= 0.0748 as the candidate genes. 17 genes were as follow: SEMA4D, PDGFRB, C5, PSMC6, CCL8, CNTFR, IL27RA, PIK3R2, IFITM1, KIF20A, GPR87, LEP, MAP2K1, COL5A1, SSTR2, SEMA6A and KLRD1. To optimize this model and identify only the most predictive genes, a stepwise Cox proportional hazards regression model was used, which resulted in a final set of 4 genes (Figure 1E). As a result of these final analyses the CIGI was constructed: CIGI = (-0.164 × SEMA4D expression) + (-0.15 × IFITM1 expression) + (0.278 × KIF20A expression) + (0.137 × GPR87 expression).

FIGURE 1

Figure 1 Identification of prognostic GRGs and IRGs in melanoma. Forest plots showing hazard ratios of the top 10 genes(sorted by p value) from the (A) TCGA and (B) GSE65904 datasets. (C) Venn diagram indicating that 182 prognostic GRGs and IRGs were identified in the GSE14520 and TCGA cohorts. (D) Cross-validation (100-fold) for tuning parameter selection in the LASSO model (upper panel) and LASSO coefficient profiles of the most relevant prognostic genes (lower panel). (E) Results of the cox proportional hazards regression model based on 17 genes.

Prognostic Analysis of CIGI in the TCGA Dataset

The risk score was calculated and distribution of the TCGA cohort is shown in Figure 2A. As based on the optimized risk value, patients were assigned to either a high- or low-risk CIGI group. Results from the Kaplan-Meier survival analyses indicated that the overall survival of patients in the high-CIGI group was significantly lower than that of patients in the low-risk group (Figure 2B; p < 0.0001). The ROC curve analyses showed that the AUC values for the 1-, 3-, and 5-year survival rates were 0.738, 0.655 and 0.702, respectively (Figure 2C). Additionally, patients with low-risk scores experienced a significantly longer disease-free survival (DSS, p <0.001; Figure 2D) and progression-free interval (PFI, p = 0.0015; Figure 2E).

FIGURE 2

Figure 2 Prognostic analysis of CIGI in the TCGA dataset. (A) Risk scores, survival times, survival status and expressions of the four genes in CIGI. (B) Kaplan–Meier analysis of the OS in the high- versus low-CIGI groups. (C) Time-dependent ROC analysis of CIGI for OS and survival status. (D) Kaplan–Meier analysis of the DSS in the high- versus low-CIGI group. (E) Kaplan–Meier analysis of the PFI in the high- versus low-CIGI group.

Verification of CIGI Based on the GSE65904 Dataset

Next, the stability and reliability of CIGI in the GSE65904 dataset were assessed. The CIGI distribution in the GSE65904 cohort is shown in Figure 3A. Results from further prognostic analysis using Kaplan-Meier revealed that the overall survival of melanoma patients in the high-CIGI group was significantly lower than that in the low-CIGI group (p=0.00054; Figure 3B). ROC curve analyses showed that AUC values for the 1-, 3- and 5-year survival rates were 0.622, 0.631 and 0.58, respectively (Figure 3C).

FIGURE 3

Figure 3 Validation of CIGI in the GSE65904 dataset. (A) Risk scores, survival times, survival status and expressions of the six genes. (B) Kaplan–Meier analysis of the OS in the high- versus low-CIGI group of the TCGA cohort. (C) Time-dependent ROC analysis of CIGI for OS and survival status.

Prognostic Value of CIGI in Different Subgroups of Patients With Melanoma

To further evaluate the clinical application value of CIGI, the prognostic value of CIGI in melanoma patients with different clinical characteristics was analyzed. As summarized in Figure 4, among melanoma patients with different clinical characteristics the OS of the high-CIGI group was significantly lower than that of the low-CIGI group (Figures 4A–L). Overall, CIGI effectively distinguished the prognosis of patients within different subgroups, further demonstrating the accuracy of CIGI.

FIGURE 4

Figure 4 Prognostic significance of CIGI in melanoma patients with different clinical features. (A) Age < 60, (B) Age > 60, (C) Male (D) Female (E) T1+T2 (F) T3+T4 (G) N0 (H) N1 (I) M0 (J) M1 (K) Stage 1+ Stage 2 (L) Stage 3+ Stage 4.

Univariate and Multivariate Cox Analysis of CIGI

Univariate and multivariate COX regressions were used to analyze the relationships among CIGI, clinical features and prognosis of melanoma patients. Results from the Univariate COX analysis demonstrated that CIGI (HR=3.675, p<0.001), age (HR=1.025, p<0.001), M stage (HR=1.893, p=0.041), N stage (HR=1.749, p<0.001), T stage (HR=1.461, p<0.001), and tumor stage (HR=1.422, p<0.001) were all significantly correlated with the prognosis of melanoma patients (Figure 5A). Multivariate COX analysis showed that CIGI (HR=3.454, p<0.001), age (HR=1.014, p=0.014), M stage (HR=3.821, p=0.007), N stage (HR=3.019, p=0.001), and T stage (HR=1.499, p<0.001) were all significantly correlated with the prognosis of melanoma patients (Figure 5B). Taken together, these results suggest that a high CIGI value served as an independent factor for prognosis in melanoma patients.

FIGURE 5

Figure 5 Univariate and multivariate cox analysis of CIGI. (A) Univariate Cox analyses of OS in the TCGA dataset. (B) Multivariate Cox analyses of OS in the TCGA dataset.

Nomogram Construction

A nomogram was generated to provide a quantitative method for predicting the probability of 1-, 3-, and 5-year OS in patients with melanoma, which could then be used in clinical practice. The nomogram integrated clinicopathological features and CIGI as based on results of the multivariate Cox regression analysis (Figure 6A). The c-index value of the nomogram was 0.74, indicating a satisfactory overlap with actual observations. Calibration curves showing calibration points in 5-years demonstrate a high degree of coincidence with the standard curves, indicating that the model provides an effective level of prediction performance (Figure 6B). In additional, results from the decision curve analysis indicate that the nomogram can better predict OS than that obtained using single clinicopathological features (Figure 6C). Therefore, this nomogram, as based on CIGI, could be used to predict the prognosis of melanoma patients in clinical practice.

FIGURE 6

Figure 6 Construction and validation of the nomogram model. (A) Nomogram model combining CIGI and traditional clinical features. (B) Five-year calibration curves of the nomogram model. (C) DCA of the nomogram model.

Glycolysis Profile in the CIGI

Next, the correlation between CIGI and glycolysis was assessed. As reported previously (17), multiple genes are involved in encoding glucose transporters and critical kinases (glycolytic genes) that directly regulate cell glycolysis (Figure 7A). Here, we analyzed the changes in mRNA of these genes between the low- versus high-CIGI group in the TCGA and GSE65904 databases. With the exception of HK2, PFKL and SLC2A1, the high-CIGI group showed higher expressions of other glycolytic genes than that of the low-CIGI group as based on the TCGA dataset (Figure 7B). Similarly, within the GSE65904 dataset, the expression of glycolytic genes were again higher in the high-CIGI group, with the exception of PFKL, PGAM1 and SLC2A1 (Figure 7C). In additional, ssGSEA scores of the “HALLMARK GLYCOLYSIS” pathway via GSVA analysis, as based on the TCGA and GSE65904 datasets were also generated (Figures 7D, E). From this analysis, we found that activity of the glycolysis pathway increased as a function of increases in CIGI. Overall, these results suggest a glycolysis-overexpression status is present in the high-CIGI group.

FIGURE 7

Figure 7 Glycolysis Profile in CIGI. (A) Summary of the glycolytic genes. (B) Expressions of glycolytic genes between the high- versus low-CIGI group in the TCGA dataset. (C) Expressions of glycolytic genes between the high- versus low-CIGI group in the GSE65904 dataset. (D, E) Correlations of ssGSEA scores of “HALLMARK GLYCOLYSIS” pathway with CIGI in the TCGA (R = 0.31, p < 0.001) and GSE65904 (R=0.51, p < 0.001) datasets. *p < 0.05, **p < 0.01, ***p < 0.001, ****p < 0.0001. ns, no significance.

Immune Profile in the CIGI

In this experiment, the relationship between CIGI and the immune status of melanoma was evaluated. To accomplish this goal a correlation analysis was first performed between CIGI and immune infiltrate cells, with the result that infiltration levels of 28 immune cell types were obtained using ssGSEA. With the exclusion of CD56bright natural killer and Type 2 T helper cells, the remaining immune cells were all significantly down-regulated in the high-CIGI group (Figure 8A). Next, the relationship between CIGI and tumor immune microenvironments were assessed. Stromal and immune scores were used to estimate tumor immune microenvironments. As shown in Figures 8B, C, CIGI was negatively correlated with immune (R=−0.46, p<0.001) and stromal (R=−0.36, p<0.001) scores. Finally, the relationship between CIGI and immune-related pathways were determined. With this analysis, we found that the activity of the Interferon-α response, Interferon-γ response, TNF-α via NF-KB, Inflammatory response and the IL2-STAT5 signaling pathway were all increased with increases in CIGI (Figure 8D). Overall, the results of these analyses indicate that the high-CIGI group exhibited immune-suppressive features.

FIGURE 8

Figure 8 Immune Profile in CIGI. (A) Distribution level of 28 types of immune cells in the high- versus low-CIGI group. (B, C) Correlations between CIGI and Immune (R = -0.46, p < 0.001) and Stromal (R = -0.36, p < 0.001) scores. (D) Immune-related pathways with differential expressions between the high- versus low-CIGI group. *p < 0.05, **p < 0.01, ***p < 0.001, ****p < 0.0001. ns, no significance.

Protein Level Validation and Functional Analysis of CIGI Genes

In Figures 9A, B the mRNA and protein level of SEMA4D and IFITM1 were decreased in A375, A875, and MeWo cells compared with melanocytes PIG1 cells. And, the mRNA and protein level of KIF20A and GPR87 were upregulated in A375, A875, and MeWo cells compared with melanocytes PIG1 cells. To verify the protein expressions of SEMA4D, IFITM1, KIF20A and GPR87 in melanoma tissue, 20 melanoma and paired normal tissue samples were compared. Results from the immunohistochemistry assay revealed that SEMA4D and IFITM1 were down-regulated, while KIF20A and GPR87 were over-expressed in melanoma tissue (Figure 9C). When analyzing the potential biological functions of these proteins in melanoma, we found that silencing SEMA4D or IFITM1 promoted, while silencing KIF20A or GPR87 inhibited the proliferation of melanoma cells in vitro (Figures 9D, E). These results suggest that SEMA4D and IFITM1 may function as tumor suppressor genes while KIF20A and GPR87 may function as oncogenes in melanoma.

FIGURE 9

Figure 9 Protein level validation and functional analysis of SEMA4D, IFITM1, KIF20A and GPR87 in melanoma versus normal samples. (A) and (B) The mRNA and prorenin level expression of four genes based on melanocytes and melanoma cells (A375, A875, and MeWo). (C) IHC analysis of SEMA4D, IFITM1, KIF20A and GPR87 in normal and melanoma tissue. (D) CCK8 assay. (E) Colony formation assays. Data represent means ± SDs; *P < 0.05, **P < 0.01 and ***P < 0.001 (versus control group).

KIF20A and GPR87 Regulate the Glycolysis Ability of Melanoma Cells

As KIF20A and GPR87 are glycolysis-related genes with poor prognostic potential for patients with melanoma, we next examined the effects of KIF20A and GPR87 on glycolytic ability within melanoma cells. After transfection with KIF20A-siRNA or GPR87-siRNA, levels of lactate, ATP and glucose uptake as well as extracellular acidification rates were all significantly reduced in these melanoma cells (Figures 10A–D). Based on these results, it seems likely KIF20A and GPR87 may function as oncogenes through their capacity to regulate glycolysis within melanoma cells. A large number of studies have shown that the AKT/LDHA pathway plays a vital role in the process of cellular glycolysis (18–20). We tried to explore the regulatory effects of KIF20A and GPR87 on the AKT/LDHA pathway. The results showed that after knocking out KIF20A and GPR87, the AKT/LDHA pathway activity of A375 cells was significantly inhibited (Figure 10E). These results suggest that: KIF20A and GPR87 may regulate the glycolytic ability of melanoma cells by regulating the AKT/LDHA pathway.

FIGURE 10

Figure 10 KIF20A and GPR87 regulate glycolysis in melanoma cells. Levels of (A) lactate, (B) ATP, (C) glucose uptake and (D) extracellular acidification rates (ECAR) as measured and analyzed in A375 cells transfected with NC-siRNA, KIF20A-siRNA or GPR87-siRNA. (E) Western blot analysis in A375 cells transfected with NC-siRNA, KIF20A-siRNA or GPR87-siRNA. *p < 0.05.

Discussion

Melanoma is a highly aggressive cancer, with a poor prognosis largely due to limited treatments currently available (21). Therefore, prognosis signatures for melanoma patients are sorely needed. With the development of bioinformatics and sequencing technology, several potential prognostic evaluation programs for melanoma patients have been generated (22–24). However, most of the parameters analyzed in these programs originate from the genome or transcriptome, with no consideration of biological processes. Therefore, these models do not adequately reflect the disease characteristics of melanoma. Glycolysis and immune microenvironments represent two significant biological hallmarks of tumors that have been demonstrated to be of value in evaluating the prognosis of patients with melanoma (25, 26). In this study, we included GRGs and IRGs to construct a CIGI as based on gene expression data in public databases. As demonstrated within different data sets and subgroups of patients, CIGI exhibits an effective degree of predictive performance, and can serve as an independent prognostic factor for melanoma patients. In this way, the accuracy displayed by CIGI for the prognosis of melanoma indicates a great potential for clinical application.

Several studies have constructed prognostic models of melanoma based on gene expression profiles from public databases. Wang et al. identified an 8-genes risk signature in melanoma patients based on entire genome (24). Song et al. constructed a 12-gene signature for survival prediction in malignant melanoma patients (27). However, these studies were derived from the entire genome or transcriptome with no consideration of biological processes. Consequently, they were simply mathematic models that did not reflect the intrinsic character of the cancer itself. Recently, researchers try to build a prognostic risk model based on tumor cell characteristics such as immune, metabolism, hypoxia, etc. Xue et al. identified an immune-related signature for predicting prognosis in melanoma, which contained 23 immune-related gene pairs (28). Shou et al. determined of hypoxia signature (20 hypoxia-related genes) to predict prognosis in melanoma (29). Our research is the first to construct a prognostic model for the two major tumor characteristics of glycolysis and immunity, which can simultaneously reflect the tumor’s immune status and glycolysis changes. In addition, our model contains less genes, which is convenient for clinical application. Importantly, the function of genes in our model has been confirmed based on experiments.

Results from increasing numbers of studies have indicated that metabolic changes in tumor microenvironments can inhibit immune cell infiltrations and other anti-tumor immunity processes through the production of immunosuppressive metabolites (30, 31). Glycolysis represents the main feature of tumor metabolism. Glycolysis forms an acidic tumor microenvironment by increasing lactic acid efflux, inhibiting anti-tumor responses mediated by T cells and inhibiting the activity of tumor infiltrating myeloid cells (12, 13). Our current findings reveal that melanoma patients in the CIGI-high group show higher expressions of glycolytic genes and exhibit immune-suppressive features.

The CIGI we constructed contains four genes, SEMA4D, IFITM1, KIF20A and GPR87. Among these, SEMA4D and IFITM1 are immune-related genes, while KIF20A and GPR87 are glycolysis-related genes. SEMA4D functions as a cell surface receptor for PLXNB1 and PLXNB2 and plays an important role in cell-cell signaling (32). SEMA4D has the capacity of inducing B-cells to aggregate and improves their viability (33),and high expression levels of SEMA4D are significantly correlated with poor prognosis for a number of cancers including colon (34), ovarian epithelial (35) and cervical (36) as well as soft tissue sarcoma (37). With regard to IFITM1, the protein encoded by this gene is an interferon-induced antiviral protein, which mediates the innate immunity of cells to influenza A H1N1, West Nile and Dengue viruses by inhibiting their initial replication (38), and plays a key role in the antiproliferative action of IFN-gamma (39). IFITM1 has been reported to be abnormally expressed in tumor tissues and it is an independent prognostic biomarker for patients with acute myeloid leukemia (40), lung adenocarcinoma (41) and gallbladder cancer (42). KIF20A plays a crucial role in cell mitosis and cell migration. Findings from several studies (43–45) have indicated that KIF20A is abnormally expressed in tumor tissues and associated with poor prognosis for patients with soft tissue sarcoma, ovarian cancer, and breast cancer. Finally, GPR87, which belongs to the G protein-couple receptor family, has been shown to be overexpressed in cancers such as pancreatic (46), non-small-cell lung (47), bladder (48) and hepatocellular (49). In this study, we present the first evidence that SEMA4D, IFITM1, KIF20A and GPR87 may possess a prognostic value for melanoma. In specific, the results of our in vitro experiments suggest that SEMA4D and IFITM1 may function as tumor suppressor genes while KIF20A and GPR87 may function as oncogenes in melanoma. Moreover, our findings provide the first indication that KIF20A and GPR87 can regulate the glycolytic ability of melanoma cells.

Although the results of our current study suggest that CIGI can serve as an effective prognostic tool for use in melanoma patients, the limitations associated with this study indicate the need for additional analyses prior to a clinical application of this protocol. First, as all samples used in our study were obtained retrospectively, the inclusion/analysis of prospective samples will need to be included for verification of our findings. Second, we focused our analyses only as related to the prognostic value and clinical significance of CIGI. The other potential functions of SEMA4D, IFITM1, KIF20A and GPR87 in CIGI will require further investigation with use of additional in vivo and in vitro experiments

Conclusion

In conclusion, here we present the findings of our endeavors to develop and validate a novel prognostic classifier for use in patients with melanoma as based on glycolysis and immune expression profiles.

Data Availability Statement

Publicly available datasets were analyzed in this study. This data can be found here: TCGA GSE65904.

Ethics Statement

The studies involving human participants were reviewed and approved by The Ethics Committee of the First Hospital of China Medical University. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

YY and YL made substantial contributions to the conception and design, acquisition, analysis, and interpretation of data. RQ and LZ was involved in drafting the manuscript and revising it critically for important intellectual content. All authors contributed to the article and approved the submitted version.

Funding

This work was supported by the National Natural Science Foundation of China (No. 81803148) and Natural Science Foundation of Liaoning Province (No.2019-KF-01-07).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

1. Shain AH, Bastian BC. From Melanocytes to Melanomas. Nat Rev Cancer (2016) 16(6):345–58. doi: 10.1038/nrc.2016.37

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Schadendorf D, van Akkooi A, Berking C, Griewank KG, Gutzmer R, Hauschild A, et al. Melanoma. Lancet (2018) 392(10151):971–84. doi: 10.1016/S0140-6736(18)31559-9

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Little EG, Eide MJ. Update on the Current State of Melanoma Incidence. Dermatol Clin (2012) 30(3):355–61. doi: 10.1016/j.det.2012.04.001

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Hanahan D, Weinberg RA. Hallmarks of Cancer: The Next Generation. Cell (2011) 144(5):646–74. doi: 10.1016/j.cell.2011.02.013

PubMed Abstract | CrossRef Full Text | Google Scholar

5. McCarthy SA, Mufson RA, Pearce EJ, Rathmell JC, Howcroft TK. Metabolic Reprogramming of the Immune Response in the Tumor Microenvironment. Cancer Biol Ther (2013) 14(4):315–8. doi: 10.4161/cbt.23616

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Wu S, Kuang H, Ke J, Pi M, Yang DH. Metabolic Reprogramming Induces Immune Cell Dysfunction in the Tumor Microenvironment of Multiple Myeloma. Front Oncol (2020) 10:591342. doi: 10.3389/fonc.2020.591342

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Wu Z, Lee YF, Yeo XH, Loo SY, Tam WL. Shifting the Gears of Metabolic Plasticity to Drive Cell State Transitions in Cancer. Cancers (Basel) (2021) 13(6). doi: 10.3390/cancers13061316

CrossRef Full Text | Google Scholar

8. Faubert B, Solmonson A, DeBerardinis RJ. Metabolic Reprogramming and Cancer Progression. Science (2020) 368(6487). doi: 10.1126/science.aaw5473

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Renner K, Singer K, Koehl GE, Geissler EK, Peter K, Siska PJ, et al. Metabolic Hallmarks of Tumor and Immune Cells in the Tumor Microenvironment. Front Immunol (2017) 8:248. doi: 10.3389/fimmu.2017.00248

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Ward PS, Thompson CB. Metabolic Reprogramming: A Cancer Hallmark Even Warburg did Not Anticipate. Cancer Cell (2012) 21(3):297–308. doi: 10.1016/j.ccr.2012.02.014

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Corbet C, Feron O. Tumour Acidosis: From the Passenger to the Driver's Seat. Nat Rev Cancer (2017) 17(10):577–93. doi: 10.1038/nrc.2017.77

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Watson MJ, Vignali P, Mullett SJ, Overacre-Delgoffe AE, Peralta RM, Grebinoski S, et al. Metabolic Support of Tumour-Infiltrating Regulatory T Cells by Lactic Acid. Nature (2021) 591(7851):645–51. doi: 10.1038/s41586-020-03045-2

PubMed Abstract | CrossRef Full Text | Google Scholar

13. >Brand A, Singer K, Koehl GE, Kolitzus M, Schoenhammer G, Thiel A, et al. LDHA-Associated Lactic Acid Production Blunts Tumor Immunosurveillance by T and NK Cells. Cell Metab (2016) 24(5):657–71. doi: 10.1016/j.cmet.2016.08.011

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Kim J, DeBerardinis RJ. Mechanisms and Implications of Metabolic Heterogeneity in Cancer. Cell Metab (2019) 30(3):434–46. doi: 10.1016/j.cmet.2019.08.013

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Pollizzi KN, Powell JD. Integrating Canonical and Metabolic Signalling Programmes in the Regulation of T Cell Responses. Nat Rev Immunol (2014) 14(7):435–46. doi: 10.1038/nri3701

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Bindea G, Mlecnik B, Tosolini M, Kirilovsky A, Waldner M, Obenauf AC, et al. Spatiotemporal Dynamics of Intratumoral Immune Cells Reveal the Immune Landscape in Human Cancer. Immunity (2013) 39(4):782–95. doi: 10.1016/j.immuni.2013.10.003

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Li L, Liang Y, Kang L, Liu Y, Gao S, Chen S, et al. Transcriptional Regulation of the Warburg Effect in Cancer by SIX1. Cancer Cell (2018) 33(3):368–85.e7. doi: 10.1016/j.ccell.2018.01.010

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Liu Y, Zhao Y, Guo L. Effects of Orexin A on Glucose Metabolism in Human Hepatocellular Carcinoma In Vitro via PI3K/Akt/mTOR-Dependent and -Independent Mechanism. Mol Cell Endocrinol (2016) 420:208–16. doi: 10.1016/j.mce.2015.11.002

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Yu T, Li L, Liu W, Ya B, Cheng H, Xin Q. Silencing of NADPH Oxidase 4 Attenuates Hypoxia Resistance in Neuroblastoma Cells SH-SY5Y by Inhibiting PI3K/Akt-Dependent Glycolysis. Oncol Res (2019) 27(5):525–32. doi: 10.3727/096504018X15179668157803

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Ma P, Xing M, Han L, Gan S, Ma J, Wu F, et al. High PD−L1 Expression Drives Glycolysis via an Akt/mTOR/HIF−1α Axis in Acute Myeloid Leukemia. Oncol Rep (2020) 43(3):999–1009. doi: 10.3892/or.2020.7477

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Luke JJ, Flaherty KT, Ribas A, Long GV. Targeted Agents and Immunotherapies: Optimizing Outcomes in Melanoma. Nat Rev Clin Oncol (2017) 14(8):463–82. doi: 10.1038/nrclinonc.2017.43

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Tong X, Qu X, Wang M. A Four-Gene-Based Prognostic Model Predicts Overall Survival in Patients With Cutaneous Melanoma. Front Oncol (2021) 11:639874. doi: 10.3389/fonc.2021.639874

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Gao Y, Li Y, Niu X, Wu Y, Guan X, Hong Y, et al. Identification and Validation of Prognostically Relevant Gene Signature in Melanoma. BioMed Res Int (2020) 2020:5323614. doi: 10.1155/2020/5323614

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Wang J, Kong PF, Wang HY, Song D, Wu WQ, Zhou HC, et al. Identification of a Gene-Related Risk Signature in Melanoma Patients Using Bioinformatic Profiling. J Oncol (2020) 2020:7526204. doi: 10.1155/2020/7526204

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Liao M, Zeng F, Li Y, Gao Q, Yin M, Deng G, et al. A Novel Predictive Model Incorporating Immune-Related Gene Signatures for Overall Survival in Melanoma Patients. Sci Rep (2020) 10(1):12462. doi: 10.1038/s41598-020-69330-2

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Tian M, Yang J, Han J, He J, Liao W. A Novel Immune Checkpoint-Related Seven-Gene Signature for Predicting Prognosis and Immunotherapy Response in Melanoma. Int Immunopharmacol (2020) 87:106821. doi: 10.1016/j.intimp.2020.106821

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Song LB, Zhang QJ, Hou XY, Xiu YY, Chen L, Song NH, et al. A Twelve-Gene Signature for Survival Prediction in Malignant Melanoma Patients. Ann Transl Med (2020) 8(6):312. doi: 10.21037/atm.2020.02.132

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Xue YN, Xue YN, Wang ZC, Mo YZ, Wang PY, Tan WQ. A Novel Signature of 23 Immunity-Related Gene Pairs Is Prognostic of Cutaneous Melanoma. Front Immunol (2020) 11:576914. doi: 10.3389/fimmu.2020.576914

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Shou Y, Yang L, Yang Y, Zhu X, Li F, Xu J. Determination of Hypoxia Signature to Predict Prognosis and the Tumor Immune Microenvironment in Melanoma. Mol Omics (2021) 17(2):307–16. doi: 10.1039/d0mo00159g

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Kaymak I, Williams KS, Cantor JR, Jones RG. Immunometabolic Interplay in the Tumor Microenvironment. Cancer Cell (2021) 39(1):28–37. doi: 10.1016/j.ccell.2020.09.004

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Certo M, Tsai CH, Pucino V, Ho PC, Mauro C. Lactate Modulation of Immune Responses in Inflammatory Versus Tumour Microenvironments. Nat Rev Immunol (2021) 21(3):151–61. doi: 10.1038/s41577-020-0406-2

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Janssen BJ, Robinson RA, Pérez-Brangulí F, Bell CH, Mitchell KJ, Siebold C, et al. Structural Basis of Semaphorin-Plexin Signalling. Nature (2010) 467(7319):1118–22. doi: 10.1038/nature09468

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Hall KT, Boumsell L, Schultze JL, Boussiotis VA, Dorfman DM, Cardoso AA, et al. Human CD100, A Novel Leukocyte Semaphorin That Promotes B-Cell Aggregation and Differentiation. Proc Natl Acad Sci USA (1996) 93(21):11780–5. doi: 10.1073/pnas.93.21.11780

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Wang JS, Jing CQ, Shan KS, Chen YZ, Guo XB, Cao ZX, et al. Semaphorin 4D and Hypoxia-Inducible Factor-1α Overexpression Is Related to Prognosis in Colorectal Carcinoma. World J Gastroenterol (2015) 21(7):2191–8. doi: 10.3748/wjg.v21.i7.2191

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Chen Y, Zhang L, Pan Y, Ren X, Hao Q. Over-Expression of Semaphorin4d, Hypoxia-Inducible Factor-1α and Vascular Endothelial Growth Factor Is Related to Poor Prognosis in Ovarian Epithelial Cancer. Int J Mol Sci (2012) 13(10):13264–74. doi: 10.3390/ijms131013264

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Liu H, Yang Y, Xiao J, Yang S, Liu Y, Kang W, et al. Semaphorin 4D Expression Is Associated With a Poor Clinical Outcome in Cervical Cancer Patients. Microvasc Res (2014) 93:1–8. doi: 10.1016/j.mvr.2014.02.007

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Ch'ng E, Tomita Y, Zhang B, He J, Hoshida Y, Qiu Y, et al. Prognostic Significance of CD100 Expression in Soft Tissue Sarcoma. Cancer (2007) 110(1):164–72. doi: 10.1002/cncr.22764

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Yánez DC, Ross S, Crompton T. The IFITM Protein Family in Adaptive Immunity. Immunology (2020) 159(4):365–72. doi: 10.1111/imm.13163

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Yang G, Xu Y, Chen X, Hu G. IFITM1 Plays an Essential Role in the Antiproliferative Action of Interferon-Gamma. Oncogene (2007) 26(4):594–603. doi: 10.1038/sj.onc.1209807

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Liu Y, Lu R, Cui W, Pang Y, Liu C, Cui L, et al. High IFITM3 Expression Predicts Adverse Prognosis in Acute Myeloid Leukemia. Cancer Gene Ther (2020) 27(1-2):38–44. doi: 10.1038/s41417-019-0093-y

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Koh YW, Han JH, Jeong D, Kim CJ. Prognostic Significance of IFITM1 Expression and Correlation With Microvessel Density and Epithelial-Mesenchymal Transition Signature in Lung Adenocarcinoma. Pathol Res Pract (2019) 215(7):152444. doi: 10.1016/j.prp.2019.152444

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Li D, Yang Z, Liu Z, Zou Q, Yuan Y. DDR2 and IFITM1 Are Prognostic Markers in Gallbladder Squamous Cell/Adenosquamous Carcinomas and Adenocarcinomas. Pathol Oncol Res (2019) 25(1):157–67. doi: 10.1007/s12253-017-0314-3

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Zhu Z, Jin Z, Zhang H, Zhang M, Sun D. Knockdown of Kif20a Inhibits Growth of Tumors in Soft Tissue Sarcoma In Vitro and In Vivo. J Cancer (2020) 11(17):5088–98. doi: 10.7150/jca.44777

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Li Y, Guo H, Wang Z, Bu H, Wang S, Wang H, et al. Cyclin F and KIF20A, FOXM1 Target Genes, Increase Proliferation and Invasion of Ovarian Cancer Cells. Exp Cell Res (2020) 395(2):112212. doi: 10.1016/j.yexcr.2020.112212

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Nakamura M, Takano A, Thang PM, Tsevegjav B, Zhu M, Yokose T, et al. Characterization of KIF20A as a Prognostic Biomarker and Therapeutic Target for Different Subtypes of Breast Cancer. Int J Oncol (2020) 57(1):277–88. doi: 10.3892/ijo.2020.5060

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Wang L, Zhou W, Zhong Y, Huo Y, Fan P, Zhan S, et al. Overexpression of G Protein-Coupled Receptor GPR87 Promotes Pancreatic Cancer Aggressiveness and Activates NF-κb Signaling Pathway. Mol Cancer (2017) 16(1):61. doi: 10.1186/s12943-017-0627-6

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Nii K, Tokunaga Y, Liu D, Zhang X, Nakano J, Ishikawa S, et al. Overexpression of G Protein-Coupled Receptor 87 Correlates With Poorer Tumor Differentiation and Higher Tumor Proliferation in non-Small-Cell Lung Cancer. Mol Clin Oncol (2014) 2(4):539–44. doi: 10.3892/mco.2014.292

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Zhang X, Liu D, Hayashida Y, Okazoe H, Hashimoto T, Ueda N, et al. G Protein-Coupled Receptor 87 (GPR87) Promotes Cell Proliferation in Human Bladder Cancer Cells. Int J Mol Sci (2015) 16(10):24319–31. doi: 10.3390/ijms161024319

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Yan M, Li H, Zhu M, Zhao F, Zhang L, Chen T, et al. G Protein-Coupled Receptor 87 (GPR87) Promotes the Growth and Metastasis of CD133⁺ Cancer Stem-Like Cells in Hepatocellular Carcinoma. PloS One (2013) 8(4):e61056. doi: 10.1371/journal.pone.0061056

PubMed Abstract | CrossRef Full Text | Google Scholar