Zhiqiang Wu1,2†
Jianmin Liu
Zijian Kang
Wangjun Yan- 1Department of Musculoskeletal Surgery, Shanghai Cancer Center, Fudan University, Shanghai, China
- 2Department of Oncology, Shanghai Medical College, Fudan University, Shanghai, China
- 3Neurovascular Center, Changhai Hospital, Naval Medical University, Shanghai, China
- 4Department of Rheumatology and Immunology, Second Affiliated Hospital of Naval Medical University, Shanghai, China
Background: Glioma is one of the most typical tumors in the central nervous system with a poor prognosis, and the optimal management strategy remains controversial. Lactate in the tumor microenvironment is known to promote cancer progression, but its impact on clinical outcomes of glioma is largely unknown.
Methods: Glioma RNA-seq data were obtained from TCGA and GCGA databases. Lactate metabolism genes (LMGs) were then evaluated to construct an LMG model in glioma using Cox and LASSO regression. Immune cell infiltration, immune checkpoint gene expression, enriched pathways, genetic alteration, and drug sensitivity were compared within the risk subgroups. Based on the risk score and clinicopathological features, a nomogram was developed to predict prognosis in patients with glioma.
Results: Five genes (LDHA, LDHB, MRS2, SL16A1, and SL25A12) showed a good prognostic value and were used to construct an LMG-based risk score. This risk score was shown as an independent prognostic factor with good predictive power in both training and validation cohorts (p < 0.001). The LMG signature was found to be correlated with the expression of immune checkpoint genes and immune infiltration and could shape the tumor microenvironment. Genetic alteration, dysregulated metabolism, and tumorigenesis pathways could be the underlying contributing factors that affect LMG risk stratification. The patients with glioma in the LMG high-risk group showed high sensitivity to EGFR inhibitors. In addition, our nomogram model could effectively predict overall survival with an area under the curve value of 0.894.
Conclusion: We explored the characteristics of LMGs in glioma and proposed an LMG-based signature. This prognostic model could predict the survival of patients with glioma and help clinical oncologists plan more individualized and effective therapeutic regimens.
1. Introduction
Glioma is one of the most typical tumors in the central nervous system, comprising approximately 80% of primary brain tumors (1). Glioma can be classified into four grades: grades II and III are defined as diffuse lower-grade gliomas, and grade IV glioma is also termed glioblastoma (2). Globally, glioma is mostly treated by surgery followed by postoperative radiotherapy and chemotherapy (3). Despite considerable advances in the development of treatments, the prognosis of patients with glioma remains poor, and the optimal management strategy remains controversial. The 5-year survival rate for patients with high-grade glioma is only ~5% (1). Although low-grade gliomas have a better prognosis than glioblastomas, they are more likely to recur and transit into high-grade gliomas (4, 5). Thus, identifying new and reliable biomarkers of glioma development and progression is urgently needed to guide clinical management and find potential targets for patients with glioma.
The presence of lactate in human tumors has long been neglected. It is now rediscovered as an important carbon source for cellular metabolism and as a signaling molecule in cancerous tissues (6). Lactate in the tumor microenvironment (TME) promotes cancer progression by creating an active niche that can shape tumor pathogenesis and evolution (7). It also reserves the acidic phenotype and increases tumor progression by modulating the TME, including cell invasion, angiogenesis, survival signaling, metastasis development, and immune surveillance escape (6). Extracellular acidosis suppresses T cell-mediated immunity, thus reducing cytolytic activity and cytokine production. Numerous studies have demonstrated that neutralization of tumor acidity in immunotherapy can improve antitumor responses (8). Extracellular lactate levels can be sensed by several cell types, including cancer cells, T cells, NK cells, dendritic cells, and macrophages, triggering intracellular signaling and strongly impacting cell behaviors and functions in the TME (9–12). Therefore, the identification of how lactate metabolism regulators mediate the TME may help improve the survival prognosis of patients with glioma.
In this study, we analyzed the expression and characterized the genomic dysregulation of lactate metabolism genes (LMGs) in glioma and determined that LMGs were associated with the tumorigenesis and prognosis of glioma. In addition, we further constructed a prediction model and revealed a high efficacy for prognosis prediction. We also explored the association between LMGs and the immune microenvironment using the constructed signature. Our study aimed to comprehensively assess the correlation of LMGs with the prognosis and immune microenvironment in glioma.
2. Methods
2.1. Data collection
The gene expression data and clinical characteristics of patients with glioma in the training set were obtained from the Chinese Glioma Genome Altas (CGGA)-RNA-seq dataset (693) (http://www.cgga.org.cn/). The Cancer Genome Atlas (TCGA) dataset (https://portal.gdc.cancer.gov) and the CGGA RNA-seq dataset (325) were used as two independent datasets. The included patients and detailed metabolism-related genes are listed in Supplementary Table 1. Genetic mutation profiles were obtained from the cBioPortal dataset (http://www.cbioportal.org/). Notably, 26 LMGs were downloaded from the Molecular Signature Database version 7.0 (MSigDB) (http://www.broad.mit.edu/gsea/msigdb/). The representative immunohistochemical (IHC) staining images of LMGs were obtained from the Human Protein Atlas (HPA) database (https://www.proteinatlas.org/).
2.2. Clinical sample collection
To determine the expression of LMGs between tumor tissues and paracancerous tissues of glioma, six samples from Changhai Hospital were collected for quantitative real-time polymerase chain reaction (qRT-PCR). The primer of the six LMGs is shown in Supplementary Table 2. Ethical approval was granted by the Ethics Committee of Shanghai Changhai Hospital (Ethics approval number: CHEC2020-164). All participants were fully informed of the research and provided the informed consent forms.
2.3. Protein–protein interaction network analysis
The STRING functional protein association network (https://string-db.org) (13) was used to construct the interaction network. The input genes consisted of LMGs to a high confidence (0.4) of active interaction. Moreover, functional enrichment analysis was carried out on the LMGs.
2.4. Construction and evaluation of the LMG-based prognostic risk score model
Univariate cox regression analysis on the LMGs in the training cohort was performed to identify the association between the expression levels of the genes and overall survival (OS) time using the survival package. Significant genes with a p-value of < 0.05 identified by univariate cox regression were further screened by least absolute shrinkage and selection operator (LASSO) cox regression. Then, six optimal LMGs were used to construct a prognostic risk score model using the following formula: risk score = ∑Coefi·Expi (Coef is the regression coefficient calculated by the LASSO). According to the calculated prognostic risk score, all patients were divided into a high-risk and a low-risk group. Kaplan-Meier survival curve and time-dependent receiver operating characteristic (ROC) curve were used to evaluate the prognostic performance of the constructed risk model.
2.5. Immune cell infiltration analysis
The abundance of immune cells in tumor tissue was calculated using the Cibersort algorithm (http://CIBERSORT.stanford.edu/), which transformed the normalized gene expression matrix into the composition of infiltrating immune cells. The LM22 signature matrix defined 22 types of immune cell components and was used as a reference expression signature with 1,000 permutations.
2.6. Gene set variation analysis and gene set enrichment analysis
Hallmark gene sets and curated gene sets were downloaded from the MsigDB dataset. Gene set variation analysis (GSVA) was utilized to evaluate the relative enrichment of the 50 hallmark pathways across samples using a non-parametric approach (14). The correlation between the LMG score and the 50 cancer Hallmark Pathway score was calculated by Spearman analysis and was visualized in a heatmap. Gene set enrichment analysis (GSEA) for the curated gene sets was performed using the R package ClusterProfiler (15) based on the differentially expressed genes. Significant gene sets were defined by false discovery rate (FDR) < 0.05 and normalized enrichment score (NES) > 2. The results were visualized by R package enrichplot.
2.7. Nomogram construction and verification
The clinical features and LMG signature were used to develop the nomogram for glioma prognostic prediction using the “rms” package (16). Calibration curves and ROC curves for 1, 2, 3, and 5 years were used to assess the accuracy and discrimination ability of the nomogram.
2.8. Single nucleotide variation and copy number variation analysis
The single nucleotide variation (SNV) and copy number variation (CNV) data of patients with glioma were obtained from the TCGA dataset. In brief, the TCGAbiolinks R package was used to download the SNV and CNV data from the TCGA database. The maftools R package was utilized to analyze SNV data and visualize the mutation waterfall plots. CNV data from the two different risk groups were analyzed using the GISTIC2.0 module on GenePattern (https://cloud.genepattern.org/gp/pages/index.jsf) with default parameters.
2.9. Drug sensitivity analysis
Genomics of Drug Sensitivity in Cancer (GDSC; https://www.cancerrxgene.org/) (17) is the largest public pharmacogenomics database, which is used for exploring molecular cancer therapy and mutation. The R package pRRophetic (18) was used to screen chemotherapeutic agents in the high- and low-risk groups.
2.10. RNA isolation, reverse transcription, and qRT-PCR
First, RNA from glioma and paracancerous tissues was extracted using Trizol (Invitrogen, USA) extraction reagents. Isolated RNA was converted into cDNA via a cDNA reverse transcription kit (TaKaRa, Cat: RR047A, Japan). Following this, the relative expression of cDNA was visualized using the SYBR qPCR Master Mix (Takara, Cat: RR820A, Japan) on a Real-Time PCR system (Applied Biosystems, Foster City, USA). β-Actin was used as a reference gene, and relative gene expression was quantified using the formula: 2–ΔΔCt.
2.11. Statistical analysis
Statistical analyses were mainly performed using R version 3.6.2. Kaplan–Meier and log-rank analyses were performed to evaluate the survival differences between different groups of patients. Student's t-test and one-way ANOVA analysis were used to estimate the differences between two groups and more than two groups. The correlation analysis was calculated using the “Spearman” method. Two-sided p < 0.05 was regarded as statistically significant.
3. Results
3.1. Correlation between LMG expression and the clinical characteristics in glioma
To identify the significance of lactate metabolism in glioma, we collected 26 LMGs from the MSigDB dataset and analyzed their expression in glioma with different clinical characteristics. It was found that the expression levels of LDHB, LDHD, SLC16A7, SLC25A12, PER2, and TP53 in 1p/19q Codel glioma were higher than those in non-Codel glioma. The expression levels of EMB, LDHA, SLC16A1, SLC16A3, SLC16A8, PARK7, and PFKFB2 were lower in 1p/19q Codel glioma than those in non-Codel glioma (Figure 1A). Furthermore, the expressions of HAGH, LDHB, LDHD, SLC16A7, SLC25A12, and TP53 were upregulated, and the expressions of EMB, LDHA, SLC16A3, and SLC16A8 were downregulated in isocitrate dehydrogenase (IDH) wild-type glioma compared with mutant glioma (Figure 1B). To investigate the protein interactions between LMGs and their biological functions, we constructed the interactions between genes related to lactate metabolism and the enriched their gene ontology (GO) functions by these interactions from the STRING database. The results showed that there were interactions between LMGs and that LDHA, LDHB, and SLC16A3 had more interactions relative to other genes (Figure 1C). These genes were enriched in the lactate transmembrane transport, lactate metabolic process, and biosynthetic process (Figure 1C). Correlation analysis showed a strong positive correlation between some of the LMGs. For example, TP53, HIF1A, MRS2, DRAS2, and SLC16A1 were positively correlated, and LDHB, PARK7, HAGH, and PNKD were also positively correlated, suggesting that they may play a synergistic role in lactate metabolism (Figure 1D). We also examined the mutational panorama of LMGs and observed that the frequencies of mutation change with the LMGs were quite low in the TCGA cohort, except for TP53, which may imply that their transcriptional levels were relatively stable in performing biological functions (Figure 1E).
Figure 1. Transcriptome and genetic alteration profiles of lactic acid metabolism-related genes (LMGs) in glioma. (A) Box plot showing the expression of LMGs in patients with X1p9q codel and non-Codel glioma. (B) Box plot showing the expression of LMGs in IDH mutate and patients with wild glioma. (C) Protein–protein interaction (PPI) network and gene ontology enrichment analysis of the LMGs. (D) Heatmap showing the correlation of the LMGs in glioma. (E) Genetic alteration profiles of the LMGs from patients with glioma inn TCGA dataset. *p < 0.05, **p < 0.01, ***p < 0.001, ****p < 0.0001.
3.2. Prognostic significance of LMGs in glioma
To explore the prognostic significance of LMGs in glioma, we performed univariate cox proportional hazards regression analysis and evaluated the relationship between the gene expression value of LMGs and the OS status of patients with glioma. We found that 11 genes were significantly associated with OS time (p < 0.05, Figure 2A). Among the 11 genes, 6 genes (HIF1A, LDHA, SLC16A1, MRS2, PFKB2, and TP53) were identified as high-risk factors (hazard ratio >1) and 5 genes (HAGH, LDHB, SLE16A1, SLC25A12, and PER2) were identified as protective factors (hazard ratio < 1). These genes were subjected to Kaplan-Meier analysis, and the representative genes are shown in Figure 2B. Moreover, the IHC profiles generated from the HPA datasets also provided information about the location of LMGs and their expression status (Supplementary Figures S1A–L). We used qPCR to further validate the LMG expression. The results showed that glioma tissues have higher expression of LDHA, HLF1A, SLC16A1, and MRS2 and lower expression of LDHB and SLC25A12 compared with paracancerous tissues (Figure 2C). These results suggested that LMGs were correlated with outcomes in patients with glioma and might play a role in glioma progression.
Figure 2. Prognostic significance of LMGs in patients with glioma. (A) Cox regression analysis showing the hazard ratios (HRs) for LMGs with 95% confidence intervals (CIs). (B) Kaplan–Meier curves for survival states of the representative genes identified by univariate cox regression. (C) Expression of six LMGs at the mRNA level by qRT-PCR. *p < 0.05, **p < 0.01.
3.3. A prognostic risk model constructed from the LMGs
Considering the prognostic value of some LMGs in glioma, we intended to construct a risk score model to evaluate the prognosis status of glioma samples more accurately. The 11 prognosis-related genes were further performed using the LASSO regression analysis to screen the most valuable predictive genes (Figures 3A, B). Five genes were selected with a p < 0.05 and the risk signature included MRS2, LDHA, SLC16A1, LDHB, and SLC25A12 (Figure 3C). Based on the risk score, we divided the patients with glioma into high-risk and low-risk groups (Figure 3D). It was found that the LMG signature could effectively differentiate cancer prognosis. Patients with glioma with high-risk scores had a shorter survival time and poorer prognosis as compared with those with low-risk scores (Figures 3E, F). In the high-risk group, the expression of MRS2, LDHA, and SLC6A1 was higher, and the expression of LDHB and SLC25A12 was lower than that in the low-risk group (Figure 3G). The risk score exhibited a high prognostic validity, with an area under curve (AUC) value of 0.729 (Figure 3H). To further validate the robustness of the LMG model, we utilized two independent sources of glioma to test the robustness of our model. The result of the Kaplan-Meier survival analysis showed that the LMG model was effective in both datasets (Figures 3I, J). The above results suggest that our LMG-based prognostic model could effectively predict glioma prognosis.
Figure 3. Construction of LMG-based risk model in patients with glioma. (A) LASSO coefficient profiling of the five LMGs. (B) Cross-validation for tuning parameter selection in the proportional hazards model. (C) The five genes' coefficients screened by lasso regression. (D) LMG signature distribution of patients with glioma in the CGGA-693 cohort. (E) Kaplan–Meier plot showing the overall survival (OS) for low- and high-LMG groups using the log-rank test. (F) Survival state distribution of patients with glioma in the CGGA-693 cohort. (G) Heatmap showing the expression level of five LMGs in the CGGA-693 cohort. (H) Time-dependent ROC analysis for the risk score in the CGGA-693 cohort. (I, J) Kaplan–Meier plot showing the OS for low- and high-risk LMG signature groups in the TCGA cohort and the CGGA-325 cohort.
3.4. Prognostic value of the LMG signature for subgroups of clinical classifications
We further compared the risk score between patients with different clinical characteristics and found that the LMG signature was associated with age, tumor grade, IDH mutation, radiotherapy, and the O-methylguanine-DNA-methyltransferase (MGMT) status (Figure 4A). However, the LMG signature was not associated with gender (Figure 4A). To further identify the role of the LMG signature for clinical subgroups, patients in the training cohort were divided into subgroups according to the clinical features, including the age, sex, IDH status, 1p/19q status, and WHO classification of the CGGA-693 cohort to compare the survival curves between the subgroup based on the median risk score. The results showed that the LMG signature also had prognostic value in each subgroup, except for patients with grade II in the training cohort (Figures 4B–H).
Figure 4. Systematic dissection of LMG signature and clinical parameters in patients with glioma. (A) The boxplots illustrating the correlation between LMG signature and different clinicopathological characteristics of patients with glioma in the CGGA-693 cohort. (B–H) Kaplan–Meier analysis for OS in different clinical subtypes including age (B), gender (C), chemotherapy (D), IDH mutation (E), MGMTp methylation (F), 1p19q (G), and tumor grade in the CGGA-693 cohort.
3.5. Construction and validation of the nomogram
To test whether the LMG signature was an independent predictor for glioma, multivariate analysis was performed using a cox proportional hazards model with the predictive risk score and clinical factors. The results showed that IDH mutation, 1p19q codel, and the LMG-based risk score were independent prognostic predictors of OS (Figure 5A). To evaluate the clinical application of the LMG signature, we constructed a nomogram model by combining the clinical features (age, sex, and neoadjuvant) with the risk score of patients with glioma (Figure 5B). The nomogram showed a favorable predictive ability for OS rates, with a high AUC value of 0.894 (Figure 5C). The calibration curve showed that the predictions of our nomogram model were in good agreement with the 1-, 2-, 3-, and 5-year actual observations (Figure 5D). Thus, our nomogram was more accurate and clinically valuable than any single independent recurrence factor.
Figure 5. Construction of the prognostic nomogram containing LMG signature in the training set. (A) Multivariate cox regression analysis of the clinical parameters and prognostic model for the overall status. (B) Nomogram for predicting the survival of patients with glioma. (C) Time-dependent ROC curves comparing prognostic accuracy of the nomogram model in the CGGA-693 cohort. (D) The calibration curves for predicting the 1-, 2-, 3-, and 5-year OS of the nomogram.
3.6. Relationship between the LMG signature and TME
To evaluate the relationship between these TME indicators and the LMG signature, we used the CIBERSORT algorithm (19) to calculate the proportion of the tumor-infiltrating immune cells in the TME in patients with glioma. The proportion of M0 macrophages, monocytes, CD4+ T memory resting cells, and Treg cells was high and that of memory B cells, NK cells, and CD8+ T cells was low in patients of the high-risk group (Figure 6A), suggesting a lack of the cell killing activity within tumors of patients in the high-risk group, most probably due to the reduced number of NK cells and less CD8+T infiltration. Meanwhile, we compared the expression level of immune-checkpoint genes between the low- and high-risk LMG signature groups and found that CD80, CD86, PDL1, PDL2, HAVCR2, and TIGIT were highly expressed in the high-risk groups (Figure 6B). The upregulation of PDL1, PDL2, HAVCR2, and TIGIT could be in an immunosuppressive status, thus facilitating tumor progression.
Figure 6. Potential implications for immune infiltration and immune checkpoint gene expression between high- and low-risk LMG signature groups. (A) Box plots to display the proportions of 22 immune infiltrating cells in patients with glioma. (B) Box plots to show the expression levels of immune checkpoint genes in high- and low-risk LMG signature groups. (C) Heat map showing the correlation between the risk score and hallmark gene sets. (D) Correlation analysis of immune scores in patients with glioma. (E, F) GSEA analysis showing the enriched pathways in the high- (E) and low-risk (F) LMG signature groups. *p < 0.05, ***p < 0.001.
3.7. Relationship between the LMG signature and cancer hallmark pathways
To evaluate the relationship between the LMG signature and cancer hallmark pathways, we collected the hallmark gene set from the MsigDB dataset and calculated the pathway activity in each sample using GSVA (14). The GSVA scores of these samples were then subjected to Spearman correlation analysis with the risk coefficients of the samples. It was found that the risk score was highly correlated with the “IL6 JAK STAT3 pathway,” “epithelial-mesenchymal response,” and “TNF-α signaling via NFKB” pathway (Figure 6C). We further used ESTIMATE (20) to calculate the immune score of each patient and then analyzed the correlation between the risk score and the immune score. It was found that there was a significant positive correlation between the risk score and the immune score, and the r-value calculated by Pearson was 0.35 (Figure 6D). GSEA was further performed to identify the involvement of pathways regulating tumorigenesis in the high-risk group. We found that multiple classic tumor-related pathways were enriched in the high-risk group, such as “hypoxia,” “MYC targets,” and “cancer poor survival” pathways (Figure 6E). In contrast, metabolism-related pathways were enriched in the low-risk group, including “Oxidative phosphorylation,” “selenoamino acid metabolism,” and “pyruvate metabolism” pathways (Figure 6F). These results indicate that the LMG signature was associated with the immune response and tumorigenesis pathway. Lactic acid metabolism may affect immune cell function to affect glioma progression.
3.8. The underlying genetic alteration of LMGs
To investigate the underlying genetic alteration of the risk level defined by LMG signature in glioma, we obtained TCGA-LGG/GBM somatic mutation profiles and analyzed the mutation landscape for patients with high- and low-risk (Figures 7A, B). We explored the top 20 mutated genes in two groups, respectively. The gene with the most mutation frequency was IDH1 (80%) in the low-risk group and that in the high-risk group was TP53 (51%). We also calculated the significant differentially mutated genes between the two risk groups. As shown in Figures 7C, D, IHD1, CIC, FUBP1, and NOTCH1 were found with a much higher mutation rate in the low-risk group than in the high-risk group. In contrast, PTEN, EGFR, TTN, and RB1 showed a much lower mutation rate in the low-risk group than in the high-risk group. Additionally, the CNV alteration landscapes of the high- and low-risk groups were further compared. The two groups had distinct CNV alteration profiles (Figures 7E, F). The above results suggested that genetic alteration may affect the expression of LMG, which further affect tumor progression.
Figure 7. Genetic alterations between high- and low-risk LMG signature groups. (A, B) Waterfall plots showing the top 20 mutation landscapes of the low- (A) and high-risk (B) LMG groups in the TCGA dataset. (C) Forest plot showing the significantly different mutated genes between the low- and high-risk LMG groups. (D) Waterfall plots showing the significantly different mutated genes between the low- and high-risk LMG signature groups. (E, F) The distribution of copy number variation (CNV) features across all chromosomes for the low- (E) and high-risk (F) LMG signature groups.
3.9. Prediction of the chemotherapy response by LMGs
To explore potential therapeutic agents for patients with glioma, GDSC analysis was used to predict the chemotherapy response between the two groups (Figures 8A–L). Several chemotherapeutic agents were significantly different between high- and low-risk groups with regard to their sensitivity. EGFR inhibitor Lapatinib and Gefitinib had higher IC50 in the low-risk group than that in the high-risk group (Figures 8A, B), suggesting EGFR inhibitors could be potential pharmacological therapies for patients with high-risk glioma. Other chemotherapeutic agents, such as Vinblastine, Nilotinib, Vorinostat, Axitinib, and Sorafenib also had distinct sensitivity between the two risk groups. The IC50 values of these agents were significantly higher in patients with high-risk glioma than the patients with low-risk glioma (Figures 8C–L), suggesting that patients with high-risk glioma may have tolerance for these drugs and ultimately be able to resist.
Figure 8. Potential candidate drugs using GDSC analysis. (A–L) Boxplots showing the IC50 values of chemotherapy agents for high- and low-risk LMG signature groups. ***p < 0.001, ****p < 0.0001.
4. Discussion
Despite considerable advances in the diagnosis and treatment of glioma, it remains cancer with high morbidity and mortality. Increasing evidence suggests that metabolic changes in tumors can modify their microenvironment, and then the newly remodeled microenvironment confers an advantage to tumor cells (21). However, no research has been published on a signature as a prognostic indicator in glioma. In the present study, we generated an LMG-based signature on a training cohort and two independent validation cohorts to predict the survival of patients with glioma. We found that patients in the high-risk group had shorter OS with a high average AUC >0.729. The lactate metabolism signature was also correlated with clinical features, immune infiltration, the expression of checkpoint molecules genetic alteration, activated pathways, and drug sensitivity. These findings illustrate that the novel LMG signature had strong predictive power in OS prediction, suggesting that lactate metabolism may affect tumor progression.
Based on LMGs, we constructed a prognostic risk model to predict the prognosis of patients with glioma. The result showed that our model had good performance in predicting survival state in both the training and validation cohorts. The LMG signature was composed of five crucial genes, including LDHA, LDHB, MRS2, SLC16A1, and SLC25A12 among which LDHA is responsible for the conversion of pyruvate into lactate and NAD+, while LDHB is responsible for converting lactate into pyruvate, fuelling oxidative metabolism (7). The inhibition of MAPK by suppressing LDHA activity decreases the production of pro-inflammatory cytokines in macrophages (22). The LDH inhibitor also reduces the ATP levels and induces oxidative stress and cell death in tumor cells (23). Several studies have shown that lactate and LDHA contribute to tumor progression (24). The downregulation of LDHB is an early event in the development of breast, prostate, and pancreatic cancer (25–27) and is correlated with high proliferation, increased tumor cell invasion, and unfavorable survival outcomes (28). MRS2 is involved in mMg 2+ uptake machinery via lactate-mediated conductance (29). Increasing MRS2 levels suppress gastric cancer cell apoptosis and reduce levels of MRS2, leading to abnormal mitochondrial function and decreased mitochondrial ΔΨ (30). SLC16A1 is responsible for transporting monocarboxylic acid metabolites. SLC16A1 plays a crucial role in cancer metabolism (31). Lactate maintains glycolysis efficiency by regulating the pH level in cells and interstitium and upregulating SLC16A1 expression. Cancer cells can increase carcinogenicity and invasion by overexpressing SLC16A1 and transporting lactic acid inward. Targeting SLC16A1 may, therefore, prove to be a promising therapeutic strategy for some cancers (32, 33). The SLC25A12 gene encodes a calcium-binding mitochondrial carrier protein involved in the exchange of aspartate for glutamate across the inner mitochondrial membrane. Loss of SLC25A12 impairs the cytosolic aspartate levels, NAD+/NADH ratio, mitochondrial respiration, and tumor growth (34). Our study showed that LDHA, MRS2, and SLC16A1 were indicators of poor survival, while LDHB and SLC25A12 were indicators of favorable prognosis, suggesting that lactate metabolism is dysregulated in glioma and this dysregulated expression is closely related to tumor progression.
Hypoxia and abnormal metabolite levels, specifically lactate, often characterize disordered metabolic states that contribute to the immunosuppression of the TME. Excessive lactate secreted by metabolism-reprogrammed cancer cells regulates immune responses via causing extracellular acidification, acting as an energy source by shuttling between different cell populations, and inhibiting the mTOR pathway in immune cells (35). We found that there were great differences in the immune microenvironment between high- and low-risk groups, showing that the abundance of CD4+ T cells, CD8+ T cells, and B cells was relatively low in the high-risk group. It was reported that CD8+ T cells activate cytolytic activity and effector function, which in turn enhances antitumor responses (36). CD4+ T cell and B cell also help and innate signals to DC functions and CD8+ T cell priming, enhancing the antitumor immunity (37). Furthermore, in addition, we found that the expression of PDL1, PDL2, HAVCR2, and TIGIT was elevated in the high-risk group. There is evidence that lactate-induced expression of GPR81 induces tumors to express immune checkpoint ligand PD-L1 (38), indicating that lactate-mediated immune dysregulation can distort host immune checkpoints in various ways to escape immune responses and promote the development and progression of glioma. Our findings are consistent with other studies showing that increased LDH activity triggers tumor immune escape by inhibiting immune function (39, 40). Although tumor immunotherapy has demonstrated good therapeutic outcomes in some cancers, few studies have reported its efficacy in gliomas (41). Our LMG-based model may be useful for screening patients for glioma immunotherapy.
Our study illustrated that patients with glioma in the two risk groups had distinct genetic alterations, activated pathways, and potentially sensitive drugs. A common characteristic of gliomas is the presence of markers for IDH1 and TP53 mutations, which influence the fate of the cells (42). The two risk groups were IDH1- and TP53-dominant mutations, respectively, suggesting discrepant oncogenesis mechanisms. These two mutations could be the cause of metabolism disorders. Recent studies reported that the IDH1 mutant glioma tissue displayed massive alterations in glycolysis and lipid metabolism compared with IDH1 wild-type glioma tissue. Both groups showed similar levels of tricarboxylic acid (TCA) cycle intermediates, but IDH1 mutant gliomas accumulated more pyruvate (43). TP53 can affect glycolysis and mitochondrial oxidative phosphorylation pathways and promote tumor development (44). Furthermore, the high-risk group had increased mutations of PTEN and EGFR. PTEN occurs in ~20% of glioblastomas but is rare in lower-grade gliomas, suggesting PTEN suppressor genes involved in the development of glioblastomas (45). There have been several EGFR gene alterations identified in gliomas, especially in glioblastomas, which acted as a prognostic factor and a predictor of treatment response in patients with glioma (46). Thus, EGFR inhibitors were predicted to be potential pharmacological therapies for high-risk groups in patients with glioma. Our study showed that the LMGs may have cross-talks with pro-oncogenic signaling pathways or other metabolism pathways. MYC could act as a glycolysis regulator and modulate both glucose transport and glucose breakdown into lactate by targeting genes (47). The activated MYC targets in high-risk groups could be a potential way to promote tumorigenesis. The other metabolism pathways, such as oxidative phosphorylation, were also reported to affect invasion ability, drug sensitivity, and prognosis (4, 48). Through drug sensitivity analysis, we found that the two LMG risk groups had distinct sensitivity to chemotherapy drugs. Our model may help to develop individualized treatment for patients with glioma.
In addition, our multivariate cox model showed that our LMG model remained independent of other factors such as IDH mutation and X1p19q deletion for patients with glioma. Using routine clinical factors associated with OS, we developed a nomogram model for clinical applications. The calibration plot showed good consistency between the prediction by nomogram and actual observation of the survival in glioma, suggesting that the proposed nomogram model could be used as a supportive tool to help clinicians distinguish, assess, and evaluate the risk and prognosis of patients with glioma.
In conclusion, we comprehensively explored the characteristics of LMGs in glioma and established a novel prognostic model to stratify different risk groups of patients with glioma, which may help plan individualized treatment and improve clinical survival outcomes of patients with glioma.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary material, further inquiries can be directed to the corresponding authors.
Ethics statement
Ethical approval was granted by the Ethics Committee of Shanghai Changhai Hospital (Ethics approval number: CHEC2020-164). The patients/participants provided their written informed consent to participate in this study. Written informed consent was obtained from the individual(s) for the publication of any potentially identifiable images or data included in this article.
Author contributions
WY, ZK, and JL designed the overall research strategy. ZW, JW, ZK, and YL performed the bioinformatics analysis and drafted and revised the manuscript. JL and WY participated in the data discussion. All authors contributed to this manuscript and approved the submitted version.
Funding
This study was supported by the National Natural Science Foundation of China (Grant 82071278).
Acknowledgments
We thank Dr. Jianming Zeng (University of Macau) and all the members of his bioinformatics team and biotrainee for generously sharing their experiences and codes.
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.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fneur.2022.1064349/full#supplementary-material
Supplementary Figure 1. Immunohistochemistry staining of lactic acid metabolism-related genes (LMGs). (A–L) Immunohistochemical (IHC) staining of LDHA (A, B), LDHB (C, D), SLC25A12 (E, F), HIA1A (G, H), MRS2 (I, J), and SLC16A (K, L) retrieved from Human Protein Atlas.
Supplementary Table 1. Clinical information in the included datasets.
Supplementary Table 2. The primer sequence of the six LMG-genes for the qRT-PCR.
References
1. Ostrom QT, Cioffi G, Waite K, Kruchko C, Barnholtz-Sloan JS. CBTRUS statistical report: primary brain and other central nervous system tumors diagnosed in the United States in 2014–2018. Neuro Oncol. (2021) 23:iii1–iii105. doi: 10.1093/neuonc/noab200
2. Louis DN, Perry A, Wesseling P, Brat DJ, Cree IA, Figarella-Branger D, et al. The 2021 WHO classification of tumors of the central nervous system: a summary. Neuro Oncol. (2021) 23:1231–51. doi: 10.1093/neuonc/noab106
3. Kieran MW, Goumnerova L, Manley P, Chi SN, Marcus KJ, Manzanera AG, et al. Phase I study of gene-mediated cytotoxic immunotherapy with AdV-tk as adjuvant to surgery and radiation for pediatric malignant glioma and recurrent ependymoma. Neuro Oncol. (2019) 21:537–46. doi: 10.1093/neuonc/noy202
4. Ceccarelli M, Barthel FP, Malta TM, Sabedot TS, Salama SR, Murray BA, et al. Molecular profiling reveals biologically discrete subsets and pathways of progression in diffuse glioma. Cell. (2016) 164:550–63. doi: 10.1016/j.cell.2015.12.028
5. Jansen E, Hamisch C, Ruess D, Heiland DH, Goldbrunner R, Ruge MI, et al. Observation after surgery for low grade glioma: long-term outcome in the light of the 2016 WHO classification. J Neurooncol. (2019) 145:501–7. doi: 10.1007/s11060-019-03316-7
6. Ippolito L, Morandi A, Giannoni E, Chiarugi P. Lactate: a metabolic driver in the tumour landscape. Trends Biochem Sci. (2019) 44:153–66. doi: 10.1016/j.tibs.2018.10.011
7. Certo M, Tsai C-H, Pucino V, Ho P-C, Mauro C. Lactate modulation of immune responses in inflammatory versus tumour microenvironments. Nat Rev Immunol. (2021) 21:151–61. doi: 10.1038/s41577-020-0406-2
8. Pilon-Thomas S, Kodumudi KN, El-Kenawi AE, Russell S, Weber AM, Luddy K, et al. Neutralization of tumor acidity improves antitumor responses to immunotherapy. Cancer Res. (2016) 76:1381–90. doi: 10.1158/0008-5472.CAN-15-1743
9. Goetze K, Walenta S, Ksiazkiewicz M, Kunz-Schughart LA, Mueller-Klieser W. Lactate enhances motility of tumor cells and inhibits monocyte migration and cytokine release. Int J Oncol. (2011) 39:453–63. doi: 10.3892/ijo.2011.1055
10. Haas R, Smith J, Rocher-Ros V, Nadkarni S, Montero-Melendez T, D'Acquisto F, et al. Lactate regulates metabolic and pro-inflammatory circuits in control of T cell migration and effector functions. PLoS Biol. (2015) 13:e1002202. doi: 10.1371/journal.pbio.1002202
11. Husain Z, Seth P, Sukhatme VP. Tumor-derived lactate and myeloid-derived suppressor cells: linking metabolism to cancer immunology. Oncoimmunology. (2013) 2:e26383. doi: 10.4161/onci.26383
12. Chen P, Zuo H, Xiong H, Kolar MJ, Chu Q, Saghatelian A, et al. Gpr132 sensing of lactate mediates tumor-macrophage interplay to promote breast cancer metastasis. Proc Natl Acad Sci USA. (2017) 114:580–5. doi: 10.1073/pnas.1614035114
13. Jensen LJ, Kuhn M, Stark M, Chaffron S, Creevey C, Muller J, et al. STRING 8—a global view on proteins and their functional interactions in 630 organisms. Nucleic Acids Res. (2009) 37:D412–6. doi: 10.1093/nar/gkn760
14. Hänzelmann S, Castelo R, Guinney J. GSVA: Gene set variation analysis for microarray and RNA-Seq data. BMC Bioinform. (2013) 14:1–7. doi: 10.1186/1471-2105-14-7
15. Wu T, Hu E, Xu S, Chen M, Guo P, Dai Z, et al. clusterProfiler 40: a universal enrichment tool for interpreting omics data. Innovation (New York, NY). (2021) 2:100141. doi: 10.1016/j.xinn.2021.100141
16. Iasonos A, Schrag D, Raj GV, Panageas KS. How to build and interpret a nomogram for cancer prognosis. J Clin Oncol Off J Am Soc Clin Oncol. (2008) 26:1364–70. doi: 10.1200/JCO.2007.12.9791
17. Yang W, Soares J, Greninger P, Edelman EJ, Lightfoot H, Forbes S, et al. Genomics of drug sensitivity in cancer (GDSC): a resource for therapeutic biomarker discovery in cancer cells. Nucleic Acids Res. (2013) 41:D955–61. doi: 10.1093/nar/gks1111
18. Geeleher P, Cox N, Huang RS. pRRophetic: an R package for prediction of clinical chemotherapeutic response from tumor gene expression levels. PLoS ONE. (2014) 9:e107468. doi: 10.1371/journal.pone.0107468
19. Newman AM, Liu CL, Green MR, Gentles AJ, Feng W, Xu Y, et al. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods. (2015) 12:453–7. doi: 10.1038/nmeth.3337
20. Yoshihara K, Shahmoradgoli M, Martínez E, Vegesna R, Kim H, Torres-Garcia W, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. (2013) 4:2612. doi: 10.1038/ncomms3612
21. Liu Z, Liu L, Jiao D, Guo C, Wang L, Li Z, et al. Association of RYR2 mutation with tumor mutation burden, prognosis, and antitumor immunity in patients with esophageal adenocarcinoma. Front Genet. (2021) 12:669694. doi: 10.3389/fgene.2021.669694
22. Song Y-J, Kim A, Kim G-T, Yu HY, Lee E-S, Park MJ, et al. Inhibition of lactate dehydrogenase A suppresses inflammatory response in RAW 2647 macrophages. Mol Med Rep. (2019) 19:629–37. doi: 10.3892/mmr.2018.9678
23. Le A, Cooper CR, Gouw AM, Dinavahi R, Maitra A, Deck LM, et al. Inhibition of lactate dehydrogenase A induces oxidative stress and inhibits tumor progression. Proc Natl Acad Sci USA. (2010) 107:2037–42. doi: 10.1073/pnas.0914433107
24. Certo M, Marone G, de Paulis A, Mauro C, Pucino V. Lactate: fueling the fire starter. Wiley Interdiscip Rev Syst Biol Med. (2020) 12:e1474. doi: 10.1002/wsbm.1474
25. Leiblich A, Cross SS, Catto JWF, Phillips JT, Leung HY, Hamdy FC, et al. Lactate dehydrogenase-B is silenced by promoter hypermethylation in human prostate cancer. Oncogene. (2006) 25:2953–60. doi: 10.1038/sj.onc.1209262
26. Maekawa M, Taniguchi T, Ishikawa J, Sugimura H, Sugano K, Kanno T. Promoter hypermethylation in cancer silences LDHB, eliminating lactate dehydrogenase isoenzymes 1–4. Clin Chem. (2003) 49:1518–20. doi: 10.1373/49.9.1518
27. Cui J, Quan M, Jiang W, Hu H, Jiao F, Li N, et al. Suppressed expression of LDHB promotes pancreatic cancer progression via inducing glycolytic phenotype. Med Oncol. (2015) 32:143. doi: 10.1007/s12032-015-0589-8
28. Yang Z, Fujii H, Mohan SV, Goronzy JJ, Weyand CM. Phosphofructokinase deficiency impairs ATP generation, autophagy, and redox balance in rheumatoid arthritis T cells. J Exp Med. (2013) 210:2119–34. doi: 10.1084/jem.20130252
29. Daw CC, Ramachandran K, Enslow BT, Maity S, Bursic B, Novello MJ, et al. Lactate elicits ER-mitochondrial Mg(2+) dynamics to integrate cellular metabolism. Cell. (2020) 183:474–89.e17. doi: 10.1016/j.cell.2020.08.049
30. Yamanaka R, Shindo Y, Oka K. Magnesium is a key player in neuronal maturation and neuropathology. Int J Mol Sci. (2019) 20:3439. doi: 10.3390/ijms20143439
31. Xue M, Song M, Yan D, Sun S, Wang Y, Fu T, et al. Effect of SLC16A1 on hepatic glucose metabolism in newborn and post-weaned Holstein bulls. Front Genet. (2022) 13:1–12. doi: 10.3389/fgene.2022.811849
32. Rassl DM, Marciniak SJ, Silvestri R, Landi S, Gemignani F. Tissue expression of lactate transporters (MCT1 and MCT4) and prognosis of malignant pleural mesothelioma (brief report). J Transl Med. (2020) 18:341. doi: 10.1186/s12967-020-02487-6
33. Hashimoto T, Hussien R, Oommen S, Gohil K, Brooks GA. Lactate sensitive transcription factor network in L6 cells: activation of MCT1 and mitochondrial biogenesis. FASEB J Off Publ Fed Am Soc Exp Biol. (2007) 21:2602–12. doi: 10.1096/fj.07-8174com
34. Alkan HF, Vesely PW, Hackl H, Foßelteder J, Schmidt DR, Vander Heiden MG, et al. Deficiency of malate-aspartate shuttle component SLC25A12 induces pulmonary metastasis. Cancer Metab. (2020) 8:1–13. doi: 10.1186/s40170-020-00232-7
35. Wang Z-H, Peng W-B, Zhang P, Yang X-P, Zhou Q. Lactate in the tumour microenvironment: from immune modulation to therapy. EBioMedicine. (2021) 73:103627. doi: 10.1016/j.ebiom.2021.103627
36. Loo Yau H, Bell E, Ettayebi I, de Almeida FC, Boukhaled GM, Shen SY, et al. DNA hypomethylating agents increase activation and cytolytic activity of CD8(+) T cells. Mol Cell. (2021) 81:1469–1483.e8. doi: 10.1016/j.molcel.2021.01.038
37. Borst J, Ahrends T, Babała N, Melief CJM, Kastenmüller W. CD4(+) T cell help in cancer immunology and immunotherapy. Nat Rev Immunol. (2018) 18:635–47. doi: 10.1038/s41577-018-0044-0
38. Feng J, Yang H, Zhang Y, Wei H, Zhu Z, Zhu B, et al. Tumor cell-derived lactate induces TAZ-dependent upregulation of PD-L1 through GPR81 in human lung cancer cells. Oncogene. (2017) 36:5829–39. doi: 10.1038/onc.2017.188
39. Husain Z, Huang Y, Seth P, Sukhatme VP. Tumor-derived lactate modifies antitumor immune response: effect on myeloid-derived suppressor cells and NK cells. J Immunol. (2013) 191:1486–95. doi: 10.4049/jimmunol.1202702
40. 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:657–71. doi: 10.1016/j.cmet.2016.08.011
41. Xu S, Tang L, Li X, Fan F, Liu Z. Immunotherapy for glioma: current management and future application. Cancer Lett. (2020) 476:1–12. doi: 10.1016/j.canlet.2020.02.002
42. Gulaia V, Shmelev M, Romanishin A, Shved N, Farniev V, Goncharov N, et al. Single-nucleus transcriptomics of IDH1- and TP53-mutant glioma stem cells displays diversified commitment on invasive cancer progenitors. Sci Rep. (2022) 12:18975. doi: 10.1038/s41598-022-23646-3
43. Zhou L, Wang Z, Hu C, Zhang C, Kovatcheva-Datchary P, Yu D, et al. Integrated metabolomics and lipidomics analyses reveal metabolic reprogramming in human glioma with IDH1 mutation. J Proteome Res. (2019) 18:960–9. doi: 10.1021/acs.jproteome.8b00663
44. Eriksson M, Ambroise G, Ouchida AT, Lima Queiroz A, Smith D, Gimenez-Cassina A, et al. Effect of mutant p53 proteins on glycolysis and mitochondrial metabolism. Mol Cell Biol. (2017) 37:e00328-17. doi: 10.1128/MCB.00328-17
45. Rasheed BK, Stenzel TT, McLendon RE, Parsons R, Friedman AH, Friedman HS, et al. gene mutations are seen in high-grade but not in low-grade gliomas. Cancer Res. (1997) 57: 4187–90.
46. Saadeh FS, Mahfouz R, Assi HI. EGFR: as a clinical marker in glioblastomas and other gliomas. Int J Biol Markers. (2018) 33:22–32. doi: 10.5301/ijbm.5000301
47. Sabnis HS, Somasagara RR, Bunting KD. Targeting MYC dependence by metabolic inhibitors in cancer. Genes (Basel). (2017) 8:114. doi: 10.3390/genes8040114
Keywords: glioma, lactate metabolism, immune infiltration, prognostic model, CGGA
Citation: Wu Z, Wang J, Li Y, Liu J, Kang Z and Yan W (2023) Characterization of a lactate metabolism-related signature for evaluation of immune features and prediction prognosis in glioma. Front. Neurol. 13:1064349. doi: 10.3389/fneur.2022.1064349
Received: 08 October 2022; Accepted: 12 December 2022;
Published: 09 January 2023.
Edited by:
Quan Cheng, Central South University, ChinaReviewed by:
Qing Yuan, Chinese Academy of Medical Sciences and Peking Union Medical College, ChinaJing Luo, Tsinghua University, China
Copyright © 2023 Wu, Wang, Li, Liu, Kang and Yan. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Jianmin Liu, 
Y2hzdHJva2UmI3gwMDA0MDsxNjMuY29t; Zijian Kang, a2FuZ3ppamlhbjE5OTImI3gwMDA0MDsxNjMuY29t; Wangjun Yan, eWFud2omI3gwMDA0MDtmdWRhbi5lZHUuY24=
†These authors have contributed equally to this work and share first authorship