ORIGINAL RESEARCH article

Front. Oncol., 13 April 2021

Sec. Neuro-Oncology and Neurosurgical Oncology

Volume 11 - 2021 | https://doi.org/10.3389/fonc.2021.592211

An Immune-Related Gene Pairs Signature Predicts Prognosis and Immune Heterogeneity in Glioblastoma

  • NZ

    Nijia Zhang 1

  • MG

    Ming Ge 1*

  • TJ

    Tao Jiang 2,3

  • XP

    Xiaoxia Peng 4

  • HS

    Hailang Sun 1

  • XQ

    Xiang Qi 1

  • ZZ

    Zhewei Zou 1

  • DL

    Dapeng Li 1

  • 1. Department of Pediatric Neurosurgery, Beijing Children’s Hospital, Capital Medical University, National Center for Children’s Health, Beijing, China

  • 2. Department of Neurosurgery, Beijing Tiantan Hospital, Capital Medical University, Beijing, China

  • 3. Department of Neurosurgery, Beijing Neurosurgical Institute, Capital Medical University, Beijing, China

  • 4. Clinical Epidemiology and Evidence-based Medicine Center, Beijing Children’s Hospital, Capital Medical University, National Center for Children’s Health, Beijing, China

Abstract

Purpose:

Glioblastoma is one of the most aggressive nervous system neoplasms. Immunotherapy represents a hot spot and has not been included in standard treatments of glioblastoma. So in this study, we aim to filtrate an immune-related gene pairs (IRGPs) signature for predicting survival and immune heterogeneity.

Methods:

We used gene expression profiles and clinical information of glioblastoma patients in the TCGA and CGGA datasets, dividing into discovery and validation cohorts. IRGPs significantly correlative with prognosis were selected to conduct an IRGPs signature. Low and high risk groups were separated by this IRGPs signature. Univariate and multivariate cox analysis were adopted to check whether risk can be a independent prognostic factor. Immune heterogeneity between different risk groups was analyzed via immune infiltration and gene set enrichment analysis (GSEA). Some different expressed genes between groups were selected to determine their relationship with immune cells and immune checkpoints.

Results:

We found an IRGPs signature consisting of 5 IRGPs. Different risk based on IRGPs signature is a independent prognostic factor both in the discovery and validation cohorts. High risk group has some immune positive cells and more immune repressive cells than low risk group by means of immune infiltration. We discovered some pathways are more active in the high risk group, leading to immune suppression, drug resistance and tumor evasion. In two specific signaling, some genes are over expressed in high risk group and positive related to immune repressive cells and immune checkpoints, which indicate aggression and immunotherapy resistance.

Conclusion:

We identified a robust IRGPs signature to predict prognosis and immune heterogeneity in glioblastoma patients. Some potential targets and pathways need to be further researched to make different patients benefit from personalized immunotherapy.

Introduction

Glioblastoma multiforme (GBM) is the most aggressive and malignant tumor in the central nervous system. The current standard therapy involving tumor resection, radiotherapy and chemotherapy implemented in 2005 and have yet to be modified (1). Despite this conventional treatments, the median survival time for GBM is despondingly 12-18 months (2). The evolvement in genomics and proteomics has made researchers acquire prominent molecular biomarkers, while few lead to a robust and innovative signature on GBM therapy (3). Immune system serves as a defensive mechanism against the formation and progression of tumors. Cancer immunotherapy has attracted attention worldwide owing to remarkable success treating advanced cancers (4). GBM immunotherapy is a research hot spot in recent years. However, majority of patients response to immunotherapy ineffectively (5). So there is an urgent need to develop new biomarkers for guiding individual immunotherapy in the treatment of GBM.

Recently, several studies have developed prognostic signatures to dividing GBM patients into different risk groups according to gene expression profiles (68). However, gene expression values measured by different platforms might generate sampling bias (9). Effective analysis of large-scale gene expression inflicts a great computational challenge that requires the use of appropriate methods. In order to eliminate the defects of data normalization and scale in gene expression data processing, some researches have invented a new method, which based on relative sequences of gene expression level, and obtained reliable results (10, 11). To date there are no studies using the new method to differentiate prognosis and immune heterogeneity of GBM patients on account of immune-related genes. So in this study, we aimed to explore prognositic immune-related gene pairs (IRGPs) in GBMs and find potential targets could be used to develop new immunotherapeutic agents.

Materials and Methods

Data Resources and Processing

This study was a retrospective study using public data. We obtained gene expression profiles and corresponding intact clinical information on 149 GBM samples in the open-source database named TCGA (https://portal.gdc.cancer.gov/) (12). Via another independent database, the Chinese Glioma Genome Atlas (CGGA) (http://www.cgga.org.cn/), we acquired molecular and clinical information of 374 GBM patients from different platforms (13, 14). For TCGA, the gene expression profile on probe level was transformed into gene symbol level. When multiple same gene symbols exist, the highest expression was selected. All expression data in both datasets were not further standardized during establishment of signature.

Identification and Verification of Immune-Related Gene Pairs (IRGPs) Signature

In this study, we selected GBM patients in the TCGA dataset as the training group, correspondingly the CGGA cohort as the validation group. Previous articles have involved how to construct IRGPs (11). We obtained 2498 immune-related genes (IRGs) in the ImmPort database (https://immport.niaid.nih.gov) (15). These immune-related genes include 17 immune classifications according to different molecular function, such as antimicrobials, antigen processing and presentation, cytokines, interleukins, natural killer cell cytotoxicity, TNF family receptors. The IRGs owning relatively high variable quantity (measured by median absolute deviation >0.5) were selected across all different platforms (16). Pairwise comparison was undertaken to obtain IRGPs using the gene expression level in one particular sample. In simple terms, a score of IRGP was 1 if IRG 1 was higher than IRG 2; conversely the IRGP score was 0. We abandoned IRGPs with low variations and the rest of IRGPs were optimized to conduct prognostic IRGPs by means of cox regression, log-rank test and multiple lasso regression. Immune-related gene pairs index (IRGPI) was produced in the training cohort by means of lasso penalized cox regression with 10-fold internal cross-validation in the glmnet package (version 3.0-2). A time-dependent receiver operating characteristic (ROC) curve (survivalROC, version 1.0.3) was generated to ensure the optimal cut-off value of IRGPI using overall survival in the TCGA dataset for distinguishing high risk from low risk patients. Between different risk groups, we used log-rank test to evaluate the established model in both the two datasets. Then we assessed whether risk based on this IRGPs signature could be an independent prognostic factor compared with other clinical factors using the univariate and multivariate cox proportional-hazards analysis.

Estimation and Comparison of the Immune Infiltration Pattern Between Different Risk Groups

We used the RNA-seq data, which include 149 GBM patients from the TCGA database and 139 patients from the same platform in the CGGA cohort, via a sample-level enrichment method named single sample gene set enrichment analysis (ssGSEA) in the GSVA package (version 1.34.0) to calculate the relative abundance of 30 immune cells of each patients in two distinct risk groups (17). From previous publications (18, 19), distinct genes which are highly expressed in each cell type were selected to represent immune populations. We used heatmap in the package ComplexHeatmap (20) to present the overall immune infiltration associated with some important mutations and clinical information in two groups. Then we detailedly compared whether average normalized enrichment scores (NES) of immune cells are significantly different between high risk and low risk groups.

Gene Set Enrichment Analysis (GSEA)

In order to determine the different expression of the same gene between high and low risk groups, we used log2 fold change of average gene expression from two groups. We conducted gene set enrichment analysis with 1000 permutations from clusterProfiler package (version 3.14.3) (21). Hallmark gene sets concerned with this study were downloaded from the Molecular Signature Database (version 7.1) (http://www.broadinstitute.org/gsea/msigdb/collections.jsp) (22). The particular gene sets whose adjusted p-value less than 0.05 were kept as statistically significant pathways.

Initial Estimate Potential Target Genes

Some genes expressed differently between two groups were extracted in specific pathways which are significant in both cohorts. To determine whether these genes could be potential targets for immunotherapy, we analyzed the relationship between these genes and immune infiltration. Immune checkpoints including PD-L1/CD274 and CTLA4 have been targets of immunotherapy in some solid cancers (23). So correlation among differential expression genes, PD-L1 and CTLA4 were also evaluated.

Statistical Analysis

All the statistical analysis were conducted using R (version 3.6.3, www.r-project.org). The log-rank test from survival package (version: 3.1-11) was adopted to analysis the survival differences in the TCGA and CGGA datasets. We perform univariate and multivariate analyses to construct the Cox proportional hazards regression model. We compared the differences of immune infiltration between groups using the Mann-Whitney test. Spearman method was adopted to assess correlations.

Results

Establishment of Prognostic Immune-Related Gene Pairs Signature

The gene expression profiles of GBM patients in the TCGA dataset (n= 149) were used as the discovery cohort. According to the evaluation criterion that median absolute deviation (MAD)>0.5, the high variation genes were retained as candidate genes. The filtered data used immune-related genes (IRGs) downloaded from the ImmPort database to establish immune-related gene pairs (IRGPs). 56554 IRGPs were conducted from 724 immune-related genes. Screening IRGPs by means of the log-rank test, cox regression and multiple lasso regression, we finally selected 5 IRGPs and calculated immune-related gene pairs index (IRGPI). These 5 IRGPs are composed of 9 unique IRGs, most of which relate to antigen processing and presentation, antimicrobials and cytokines (Table 1). Then each GBM patient’s risk score in the training group was operated. Through one year time-dependent ROC curve analysis, the best cut-off value of the IRGPI was 0.197 for criterion to distinguish different risk groups (Figure 1). High risk group have a worse prognosis than low risk groups (Figure 2A), indicating IRGPI dividing patients significantly. Risk and other clinical factors including age, gender, radiotherapy and chemotherapy were analyzed in the univariate and multivariate Cox analysis (Table 2). Risk based on IRGPI signature showed statistical significant in both Cox analysis. Supplementary Table 1 showed detailed information of each patient including survival time, status, riskscore and risk group in the TCGA dataset.

Table 1

IRG1CategoryIRG2CategoryCoef
HSPA6Antigen Processing and PresentationBMP2TGFb Family Member0.464
PSMC3Antigen Processing and PresentationMDKCytokines-0.358
FGF2AntimicrobialsOSMRCytokine Receptors-0.460
PPP4CAntimicrobialsMDKCytokines-0.293
LEFTY2CytokinesMSTNCytokines0.555

Immune-related gene pairs information and coefficients.

Figure 1

Figure 2

Table 2

DatasetsFactorUnivariate cox analysisMultivariate cox analysis
HR(95%CI)P-valueHR(95%CI)P-value
TCGAAge1.02(1.00-1.04)0.011.01(0.99-1.02)0.44
Gender0.92(0.62-1.35)0.650.77(0.49-1.21)0.25
Radiotherapy0.16(0.09-0.27)3.22E-120.14(0.07-0.29)1.45E-07
Chemotherapy0.45(0.30-0.67)8.72E-050.71(0.42-1.22)0.22
Risk2.79(1.87-4.16)5.27E-073.09(2.01-4.74)2.40E-07
CGGAAge1.01(1.00-1.02)0.031.01(1.00-1.02)0.04
Gender0.92(0.74-1.16)0.490.91(0.71-1.15)0.42
Radiotherapy0.69(0.52-0.92)0.010.66(0.49-0.90)0.01
Chemotherapy0.46(0.34-0.61)9.04E-080.49(0.37-0.66)3.01E-06
Risk1.72(1.37-2.16)2.84E-061.70(1.33-2.16)1.66E-05

Univariate and multivariate analysis of prognostic factors in the TCGA and CGGA dataset.

Validation of Prognostic Prognostic Gene Pairs Signature

The gene expression data of GBM patients in the CGGA dataset (n= 374) were used to verify whether IRGPI signature has the same role as that in the TCGA dataset. Different risk groups according to the IRGP signature have distinct prognosis, as similar with the result of training group (Figure 2B). Through the univariate and multivariate Cox analysis, risk groups also show significantly independent factor for survival (Table 2). Similar items of each patient in the CGGA database were displayed (Supplementary Table 2).

Immune Infiltration Pattern Between Different Risk Groups

Previous publication has elaborated immune infiltration is related to cancer prognosis (24). And single sample gene set enrichment analysis (ssGSEA) has been adopted to assess immune cell infiltration (25). For each patients from different risk groups in the TCGA dataset and a subset from the same platform in the CGGA database, we made use of ssGSEA to calculate the relative abundance of 30 immune cells from 16 immune populations. Based on overall immune infiltration (Figure 3), we found that majority of patients in high risk group had higher immune infiltration and belonged to mesenchymal subtype based on Verhaak subtype (26). While in the low risk group, more patients harbored MGMT promoter methylation and IDH mutation and proneural subtype increased, which indicated better prognosis (2629). Comparison and significant differences of specific immune cells between high and low risk groups have been exhibited (Figure 4). Although some immune positive cells including activated CD4 T cell (30), activated dendritic cell, effector memory CD8 T cell, natural killer cell and type 1 T helper cell are higher in the high risk group, more immune repressive cells involving gamma delta T cell (31), macrophage M2 (32), MDSC (myeloid derived suppressor cells) (33), neutrophil (34), regulatory T cell (35) contribute to worse prognosis, similar results from the TCGA and CGGA datasets.

Figure 3

Figure 4

Functional Pathways Assessment of Different Groups

In order to find out which pathways change more in high risk group than low risk group, we adopted gene set enrichment analysis (GSEA). The bubble plot showed top twelve pathways enriched in the high risk group from two cohorts (Figure 5). Statistic value of these pathways were presented (Supplementary Tables 3 and 4). Epithelial-Mesenchymal Transition (EMT), interferon gamma response, IL2/STAT5 signaling and some other pathways are activated in the high risk group, which manifest immune suppression (36), drug resistance (37), tumor evasion (38). Particularly in the epithelial-mesenchymal transition pathway (Figures 6A, C), we can see some genes including COL6A3, COL1A1, COL1A2 and LRRC15 express higher in the high risk group (Figure 6B), which also can be seen in high risk group from the CGGA dataset (Figure 6D). From the response to IFNγ pathway (Figures 7A, C), IDO1, IL6 and PTGS2 expressed higher in the high risk group between two cohorts (Figures 7B, D).

Figure 5

Figure 6

Figure 7

Relationship Between Selected Genes and Immune Reaction

Immune checkpoints such as PD-L1/CD274 and CTLA4 are often overexpressed in tumor cells to devitalize effector T cells and suppress immune responses, which could be targets for immunotherapy (3941). To determine whether these selected genes have an impact on immune inhibition, they were estimated through their correlation with immune cells and two immune checkpoints. The selected genes were mostly correlated to gamma delta T cell, regulatory T cell, central memory T cell (42), which are immune negative or noneffective cells (Figure 8 and Supplementary Figures 1 and 2). Moreover, all seven genes had positive correlation with PD-L1/CD274 and CTLA4 (Figure 9 and Supplementary Tables 5 and 6).

Figure 8

Figure 9

Discussion

Glioblastoma multiforme (GBM) is one of the lethal tumor in the central nervous system. The standard treatment includes surgery, radiotherapy, chemotherapy. However, no novel regimen has been introduced into clinical practice to obviously improve survival of GBM patients in recent years. Immune system is one of defensive lines against tumor and immunotherapy has become a research focus among different types of cancer to prolong patients’ survival. The current immunotherapy comprises cell therapy, peptide vaccine, and immune checkpoint inhibitors (33). As for GBM, only less than 10% patients respond to checkpoint inhibition (43). This phenomenon implies advanced malignancies have complex interactions with the immune system. It is urgent to find out new immune related biomarkers to predict prognosis and give each patient personalized treatment.

As for exploring prognostic signatures, it’s difficult to analyze the gene expression profiles of every tumor specimen. Data need to be standardized properly due to different sequencing platforms and tumor samples. In order to avoid technical bias from different platforms, we adopted an new method not require scaling and normalization (44). In this study, we found out an immune-related gene pairs (IRGPs) signature independently predict prognosis. Proneual subtype, IDH mutant and MGMT promoter methylation are well known prognostic factors for longer survival of GBM patients (27, 28, 45, 46). A substantial proportion of the high risk group patients belong to mesenchymal subtype, while proneural subtype increases in the low risk group. Compare with high risk group, more patients in the low risk group have IDH mutation and MGMT promoter methylation. Other factors including EGFR amplification, TERT promoter mutation, TP53 mutation don’t show difference between two groups because limited mutation information of GBM patients in CGGA dataset. Meanwhile, different risk group based on this signature has distinct immune cells infiltration. High risk group has some immune positive cells and more immune suppressive cells, indicating battle between immune cells and tumor cells is a continuous process of elimination, equilibrium, and escape (47). Previous study has shown that macrophage dominated and high M2 macrophage polarization consistent with an immunosuppressed tumor microenvironment, which foreboded a poor outcome (48). In our study, high risk group patients have more M2 macrophage than low risk group patients. It manifests there is a discrepancy of immunosuppressed microenvironment between different groups. Furthermore, we discover epithelial-mesenchymal transition signaling (49), TNFαsignaling via NF-κB (50), IL-2/STAT5 signaling (36), interferonγ response signaling and some other pathways are activated in high risk group, most of which trigger tumor cells malignant progression, immune evasion, metastasis and poor prognosis. In two specific pathways (epithelial-mesenchymal transition and interferonγ response), we obtain some genes express higher in high risk group than that in low risk group. Previous article has verified silence of COL6A3 and COL1A2 can inhibit tumor cell proliferation, migration, and invasion in the gastric cancer (51). COL6A3 also has similar role in the colorectal cancer (52). Another studies have reported that COL1A1 is related to facilitate cell invasion in glioma (53). LRRC15 has been confirmed as a immunotherapy resistant target in the single-cell RNA sequencing experiment (54). Experiments showed the combination of IL-6 and PD-1/PD-L1 inhibitors promotes antitumor immunity (55) and PTGS2 deletion sensitized tumors to immunotherapy (56). One research has demonstrated increased levels of IDO1 in the glioblastoma cell had positive correlation with human-infiltrating T cells leading to poor prognosis (57). In murine GBM model, IDO1 inhibition combine with radiotherapy and PD-1 blockade increased survival (58). Most important of all, IDO1 inhibitor could benefit a subset of patients with recurrent malignant glioma in a phase 1 study (59). In our study, we found that these selected genes were most correlated to immune repressive cells and noneffective memory T cells. They also had a positive relationship with CD274 and CTLA4. All these findings indicated that high risk group might be more aggressive and immunosuppressive than low risk group. Immune heterogeneity existed between different risk groups. The mechanism of tumor invasion and immune resistance involving COL6A3, COL1A1, COL1A2, LRRC15, IDO1, IL-6 and PTGS2 need to be further researched aiming to improve prognosis of aggressive glioblastoma.

Our study also have some limitations. First, this signature is based on gene expression profiles, which are not widespread applied owing to expenditure and high requirement of bioinformatics knowledge. In addition, the mechanism of seven genes in different groups has not been explored, though they were potential targets for immunotherapy.

Conclusion

In conclusion, we identified an immune related gene pairs signature. Different groups based on this signature have distinct prognosis and immune heterogeneity. Some biological processes and genes have indicated the poor prognosis is related to tumor immune evasion, malignant progression, metastasis. The detailed function of these targets need to be explored to correct immune dysfunction and make all patients benefit from personalized immunotherapy.

Statements

Data availability statement

Publicly available datasets were analyzed in this study. The TCGA-GBM dataset can be found from TCGA database (https://cancergenome.nih.gov/). The GBM data from CGGA database can be gained in the website (http://www.cgga.org.cn/).

Ethics statement

Written informed consent for participation was not required for this study in accordance with the national legislation and the institutional requirements.

Author contributions

NZ and MG contributed to the conception of the study. HS, XQ, ZZ, and DL contributed to manuscript preparation. NZ performed the data analyses and wrote the manuscript. MG, TJ, and XP helped perform the analysis with constructive discussions. All authors contributed to the article and approved the submitted version.

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.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2021.592211/full#supplementary-material

References

  • 1

    StuppRMasonWPvan den BentMJWellerMFisherBTaphoornMJet al. Radiotherapy plus concomitant and adjuvant temozolomide for glioblastoma. N Engl J Med (2005) 352(10):987–96. doi: 10.1056/NEJMoa043330

  • 2

    PollackIFHamiltonRLJamesCDFinkelsteinSDBurnhamJYatesAJet al. Rarity of PTEN deletions and EGFR amplification in malignantgliomas of childhood: Results from the Children’s Cancer Group 945 cohort. J Neurosurg (2006) 105:418–24. doi: 10.3171/ped.2006.105.5.418

  • 3

    SasmitaAOWongYPLingA. Biomarkers and therapeutic advances in glioblastoma multiforme. Asia Pac J Clin Oncol (2018) 14(1):4051. doi: 10.1111/ajco.12756

  • 4

    YangY. Cancer immunotherapy: harnessing the immune system to battle cancer. J Clin Invest (2015) 125(9):3335–7. doi: 10.1172/JCI83871

  • 5

    JacksonCMLimM. Immunotherapy for glioblastoma: playing chess, not checkers. Clin Cancer Res (2018) 24:4059–61. doi: 10.1158/1078-0432.CCR-18-0491

  • 6

    ChengWRenXZhangCCaiJLiuYHanSet al. Bioinformatic profiling identifies an immune-related risk signature for glioblastoma. Neurology (2016) 86(24):2226–34. doi: 10.1212/WNL.0000000000002770

  • 7

    JiaDLiSLiDXueHYangDLiuY. Mining TCGA database for genes of prognostic value in glioblastoma microenvironment. Aging (Albany NY) (2018) 10(4):592605. doi: 10.18632/aging.101415

  • 8

    YinAALuNEtcheverryAAubryMBarnholtz-SloanJZhangLHet al. A novel prognostic six-CpG signature in glioblastomas. CNS Neurosci Ther (2018) 24(3):167–77. doi: 10.1111/cns.12786

  • 9

    LeekJTScharpfRBBravoHCSimchaDLangmeadBJohnsonWEet al. Tackling the widespread and critical impact of batch effects in high-throughput data. Nat Rev Genet (2010) 11(10):733–9. doi: 10.1038/nrg2825

  • 10

    HeinäniemiMNykterMKramerRWienecke-BaldacchinoASinkkonenLZhouJXet al. Gene-pair expression signatures reveal lineage control. Nat Methods (2013) 10(6):577–83. doi: 10.1038/nmeth.2445

  • 11

    LiBCuiYDiehnMLiR. Development and validation of an individualized immune prognostic signature in early-stage nonsqua-mous non-small cell lung cancer. JAMA Oncol (2017) 3(11):1529–37. doi: 10.1001/jamaoncol.2017.1609

  • 12

    Cancer GARN. Comprehensive genomic characterization defines human glioblastoma genes and core pathways. Nature (2008) 455(7216):1061–8. doi: 10.1038/nature07385

  • 13

    WangYQianTYouGPengXChenCYouYet al. Localizing seizure-susceptible brain regions associated with low-grade gliomas using voxel-based lesion-symptom mapping. Neuro Oncol (2015) 17(2):282–8. doi: 10.1093/neuonc/nou130

  • 14

    ZhaoZMengFWangWWangZZhangCJiangT. Comprehensive RNA-seq transcriptomic profiling in the malignant progression of gliomas. Sci Data (2017) 4:170024. doi: 10.1038/sdata.2017.24

  • 15

    BhattacharyaSAndorfSGomesLDunnPSchaeferHPontiusJet al. ImmPort: disseminating data to the public for the future of immunology. Immunol Res (2014) 58(2-3):234–9. doi: 10.1007/s12026-014-8516-1

  • 16

    GuinneyJDienstmannRWangXde ReynièsASchlickerASonesonCet al. The consensus molecular subtypes of colorectal cancer. Nat Med (2015) 21(11):1350–6. doi: 10.1038/nm.3967

  • 17

    HanzelmannSCasteloRGuinneyJ. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinf (2013) 14:7. doi: 10.1186/1471-2105-14-7

  • 18

    AngelovaMCharoentongPHacklHFischerMLSnajderRKrogsdamAMet al. Characterization of the immunophenotypes and antigenomes of colorectal cancers reveals distinct tumor escape mechanisms and novel targets for immunotherapy. Genome Biol (2015) 15(54). doi: 10.1186/s13059-015-0620-6

  • 19

    BindeaGMlecnikBTosoliniMKirilovskyAWaldnerMObenaufACet al. Spatiotemporal dynamics of intratumoral immune cells reveal the immune landscape in human cancer. Immunity (2013) 39:782–95. doi: 10.1016/j.immuni.2013.10.003

  • 20

    Gu ZEilsRSchlesnerM. Complex heatmaps reveal patterns and correlations in multidimensional genomic data. Bioinformatics (2016) 32:2847–9. doi: 10.1093/bioinformatics/btw313

  • 21

    YuGWangLHanYQing-YuH. ClusterProfiler: an R package for comparing biological themes among gene clusters. OMICS: A J Integr Biol (2012) 16(5):284–7. doi: 10.1089/omi.2011.0118

  • 22

    LiberzonABirgerCThorvaldsdóttirHGhandiMMesirovJPTamayoP. The Molecular Signatures Database (MSigDB) hallmark gene set collection. Cell Syst (2015) 1(6):417–25. doi: 10.1016/j.cels.2015.12.004

  • 23

    Gong JChehrazi-RaffleAReddiSSalgiaR. Development of PD-1 and PD-L1 inhibitors as a form of cancer immunotherapy: a comprehensive review of registration trials and future considerations. J Immunother Cancer (2018) 6:8. doi: 10.1186/s40425-018-0316-z

  • 24

    ZhouRZhangJZengDSunHRongXShiMet al. Immune cell infiltration as a biomarker for the diagnosis and prognosis of stage I-III colon cancer. Cancer Immunol Immunother (2019) 68(3):433–42. doi: 10.1007/s00262-018-2289-7

  • 25

    YeLZhangTKangZGuoGSunYLinKet al. Tumor-Infiltrating Immune Cells Act as a Marker for Prognosis in Colorectal Cancer. Front Immunol (2019) 10(2368). doi: 10.3389/fimmu.2019.02368

  • 26

    VerhaakRGHoadleyKAPurdomEWangVQiYWilkersonMDet al. Integrated genomic analysis identifies clinically relevant subtypes of glioblastoma characterized by abnormalities in PDGFRA, IDH1, EGFR, and NF1. Cancer Cell (2010) 17:98110. doi: 10.1016/j.ccr.2009.12.020

  • 27

    Labussière MBBMokhtari KMBoisselierBMokhtariKDi StefanoA-LRahimianARossettoMet al. Combined analysis of TERT, EGFR, and IDH status defines distinct prognostic glioblastoma classes. Neurology (2014) 83:1200–6. doi: 10.1212/WNL.0000000000000814

  • 28

    BinabajMMBahramiAShahidSalesSJoodiMMashhadMJHassanianSMet al. The prognostic value of MGMT promoter methylation in glioblastoma: A meta-analysis of clinical trials. J Cell Physiol (2018) 233:378–86. doi: 10.1002/jcp.25896

  • 29

    BarbagalloGMParatoreSCaltabianoRPalmucciSParraHSPriviteraGet al. Long-term therapy with temozolomide is a feasible option for newly diagnosed glioblastoma: a single-institution experience with as many as 101 temozolomide cycles. Neurosurg Focus (2014) 37:E4. doi: 10.3171/2014.9.FOCUS14502

  • 30

    BorstJAhrendsTBąbałaNMeliefCKastenmüllerW. CD4+ T cell help in cancer immunology and immunotherapy. Nat Rev Immunol (2018) 18(10):635–47. doi: 10.1038/s41577-018-0044-0

  • 31

    FlemingCMorrisseySCaiYYanJ. γδ T Cells: Unexpected Regulators of Cancer Development and Progression. Trends Cancer (2017) 3(8):561–70. doi: 10.1016/j.trecan.2017.06.003

  • 32

    MantovaniASozzaniSLocatiMAllavenaP. Macrophage polarization: tumor-associated macrophages as a paradigm for polarized M2 mononuclear phagocytes. Trends Immunol (2002) 23:549–55. doi: 10.1016/S1471-4906(02)02302-5

  • 33

    MiyauchiJTTsirkaSE. Advances in immunotherapeutic research for glioma therapy. J Neurol (2018) 265(4):741–56. doi: 10.1007/s00415-017-8695-5

  • 34

    SwierczakAMouchemoreKAHamiltonJAAndersonRL. Neutrophils: important contributors to tumor progression and metastasis. Cancer Metastasis Rev (2015) 34(4):735–51. doi: 10.1007/s10555-015-9594-9

  • 35

    TanakaASakaguchiS. Regulatory T cells in cancer immunotherapy. Cell Res (2017) 27(1):109–18. doi: 10.1038/cr.2016.151

  • 36

    RaniAMurphyJJ. STAT5 in Cancer and Immunity. J Interferon Cytokine Res (2016) 36(4):226–37. doi: 10.1089/jir.2015.0054

  • 37

    DuBShimJS. Targeting Epithelial-Mesenchymal Transition (EMT) to Overcome Drug Resistance in Cancer. Molecules (2016) 21(7):965. doi: 10.3390/molecules21070965

  • 38

    CastroFCardosoAPGonçalvesRMSerreKOliveiraMJ. Interferon-Gamma at the Crossroads of Tumor Immune Surveillance or Evasion. Front Immunol (2018) 9:847. doi: 10.3389/fimmu.2018.00847

  • 39

    ZhangJDangFRenJWeiW. Biochemical Aspects of PD-L1 Regulation in Cancer Immunotherapy. Trends Biochem Sci (2018) 43:1014–32. doi: 10.1016/j.tibs.2018.09.004

  • 40

    RibasAWolchokJD. Cancer immunotherapy using checkpoint blockade. Science (2018) 359:1350–5. doi: 10.1126/science.aar4060

  • 41

    HavelJJChowellD. The evolving landscape of biomarkers for checkpoint inhibitor immunotherapy. Nat Rev Cancer (2019) 19:133–50. doi: 10.1038/s41568-019-0116-x

  • 42

    BuschDHFräßleSPSommermeyerDBuchholzVRRiddellSR. Role of memory T cell subsets for adoptive immunotherapy. Semin Immunol (2016) 28:2834. doi: 10.1016/j.smim.2016.02.001

  • 43

    OmuroAVlahovicGLimMSahebjamSBaehringJCloughesyTet al. Nivolumab with or without ipilimumab in patients with recurrent glioblastoma: results from exploratory phase I cohorts of CheckMate 143. Neuro Oncol (2018) 20(5):674–86. doi: 10.1093/neuonc/nox208

  • 44

    EddyJASungJGemanDPriceND. Relative expression analysis for molecular cancer diagnosis and prognosis. Technol Cancer Res Treat (2010) 9(2):149–59. doi: 10.1177/153303461000900204

  • 45

    AldapeKZadehGMansouriSReifenbergerGvon DeimlingA. Glioblastoma: pathology, molecular mechanisms and markers. Acta Neuropathol (2015) 129:829–48. doi: 10.1007/s00401-015-1432-1

  • 46

    GuanXVengoecheaJZhengSSloanAEChenYBratDJet al. Molecular subtypes of glioblastoma are relevant to lower grade glioma. PloS One (2014) 9:e91216. doi: 10.1371/journal.pone.0091216

  • 47

    DunnGPOldLJSchreiberRD. The three Es of cancer immunoediting. Annu Rev Immunol (2004) 22:329–60. doi: 10.1146/annurev.immunol.22.012703.104803

  • 48

    ThorssonVGibbsDLBrownSDWolfDBortoneDSYangT-HOet al. The Immune Landscape of Cancer. Immunity (2018) 48:812–30. doi: 10.1016/j.immuni.2018.03.023

  • 49

    DiepenbruckMChristoforiG. Epithelial-mesenchymal transition (EMT) and metastasis: yes, no, maybe. Curr Opin Cell Biol (2016) 43:713. doi: 10.1016/j.ceb.2016.06.002

  • 50

    LimSOLiCWXiaWChaJHChanLCWuYet al. Deubiquitination and Stabilization of PD-L1 by CSN5. Cancer Cell (2016) 30(6):925–39. doi: 10.1016/j.ccell.2016.10.010

  • 51

    AoRGuanLWangYWangJ-N. Silencing of COL1A2, COL6A3, and THBS2 inhibits gastric cancer cell proliferation, migration, and invasion while promoting apoptosis through the PI3k-Akt signaling pathway. J Cell Biochem (2018) 119:4420–34. doi: 10.1002/jcb.26524

  • 52

    LiuWLiLYeHTaoH. Role of COL6A3 in colorectal cancer. Oncol Rep (2018) 39:2527–36. doi: 10.3892/or.2018.6331

  • 53

    SunSWangYWuYGaoYLiQAbdulrahmanAAet al. Identification of COL1A1 as an invasion−related gene in malignant astrocytoma. Int J Oncol (2018) 53:2542–54. doi: 10.3892/ijo.2018.4568

  • 54

    DominguezCXMüllerSKeerthivasanSKoeppenHHungJGierkeSet al. Single-Cell RNA Sequencing Reveals Stromal Evolution into LRRC15+ Myofibroblasts as a Determinant of Patient Response to Cancer Immunotherapy. Cancer Discovery (2020) 10(2):232–53. doi: 10.1158/2159-8290.CD-19-0644

  • 55

    TsukamotoHFujiedaKMiyashitaAFukushimaSIkedaTKuboYet al. Combined Blockade of IL6 and PD-1/PD-L1 Signaling Abrogates Mutual Regulation of Their Immunosuppressive Effects in the Tumor Microenvironment. Cancer Res (2018) 78:5011–22. doi: 10.1158/0008-5472

  • 56

    MarkosyanNLiJSunYHRichmanLPLinJHYanFet al. Tumor cell-intrinsic EPHA2 suppresses anti-tumor immunity by regulating PTGS2 (COX-2). J Clin Invest (2019) 129:3594–609. doi: 10.1172/JCI127755

  • 57

    ZhaiLLadomerskyELauingKLWuMGenetMGritsinaGet al. Infiltrating T Cells Increase IDO1 Expression in Glioblastoma and Contribute to Decreased Patient Survival. Clin Cancer Res (2017) 23:6650–60. doi: 10.1158/1078-0432.CCR-17-0120

  • 58

    LadomerskyEZhaiLLenzenALauingKLQianJScholtensDMet al. IDO1 Inhibition Synergizes with Radiation and PD-1 Blockade to Durably Increase Survival Against Advanced Glioblastoma. Clin Cancer Res (2018) 24:2559–73. doi: 10.1158/1078-0432.CCR-17-3573

  • 59

    ReardonDADesjardinsARixeOCloughesyTAlekarSWilliamsJHet al. A phase 1 study of PF-06840003, an oral indoleamine 2,3-dioxygenase 1 (IDO1) inhibitor in patients with recurrent malignant glioma. Invest New Drugs (2020) 38:1784–95. doi: 10.1007/s10637-020-00950-1

Summary

Keywords

prognosis, immune heterogeneity, immune-related gene pairs signature, glioblastoma, immunotherapy

Citation

Zhang N, Ge M, Jiang T, Peng X, Sun H, Qi X, Zou Z and Li D (2021) An Immune-Related Gene Pairs Signature Predicts Prognosis and Immune Heterogeneity in Glioblastoma. Front. Oncol. 11:592211. doi: 10.3389/fonc.2021.592211

Received

06 August 2020

Accepted

16 March 2021

Published

13 April 2021

Volume

11 - 2021

Edited by

Leonora Balaj, Massachusetts General Hospital and Harvard Medical School, United States

Reviewed by

Juan Manuel Sepulveda Sanchez, University Hospital October 12, Spain; Sara Grazia Maria Piccirillo, University of New Mexico Health Sciences Center, United States

Updates

Copyright

*Correspondence: Ming Ge,

This article was submitted to Neuro-Oncology and Neurosurgical Oncology, a section of the journal Frontiers in Oncology

Disclaimer

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.

Outline

Figures

Cite article

Copy to clipboard


Export citation file


Share article

Article metrics