A Multi-Element Expression Score Is A Prognostic Factor In Glioblastoma Multiforme.

PURPOSE: Glioblastoma multiforme (GBM) is a highly malignant tumor of the central nervous system. Although primary GBM patients receive extensive therapies, tumors may recur within months, and there is no objective and scientific method to predict prognosis. Adoptive immunotherapy holds great promise for GBM treatment. However, the expression profiles of the tumor-associated antigens (TAAs) and tumor immune microenvironment (TME) genes used in immunotherapy of GBM patients have not been fully described. The present study aimed to develop a predictive tool to evaluate patient survival based on full analysis of the expression levels of TAAs and TME genes.
METHODS: Expression profiles of a panel of 87 TAAs and 8 TME genes significantly correlated with poor prognosis were evaluated in 44 GBM patients and 10 normal brain tissues using quantitative real-time polymerase chain reaction (qRT-PCR). A linear formula (the LASSO algorithm based in the R package) weighted by regression coefficients was used to develop a multi-element expression score to predict prognosis; this formula was cross-validated by the leave-one-out method in different GBM cohorts.
RESULTS: After analysis of gene expression, clinical features, and overall survival (OS), a total of 8 TAAs (CHI3L1, EZH2, TRIOBP, PCNA, PIK3R1, PRKDC, SART3 and EPCAM), 1 TME gene (FOXP3) and 4 clinical features (neutrophil-to-lymphocyte (NLR), number of basophils (BAS), age and treatment with standard radiotherapy and chemotherapy) were included in the formula. There were significant differences between high and low scoring groups identified using the formula in different GBM cohorts (TCGA (n=732) and GEO databases (n=84)), implying poor and good prognosis, respectively.
CONCLUSION: The multi-element expression score was significantly associated with OS of GBM patients. The improve understanding of TAAs and TMEs and well-defined formula could be implemented in immunotherapy for GBM to provide better care.

Introduction

Glioblastoma multiforme (GBM, also known as astrocytoma grade IV) is the most common and deadliest primary brain tumor, representing 30% of all central nervous system tumors.1 GBM exhibits various pathophysiological features of malignancy, including necrosis, vascular proliferation and pleomorphism.2 Owing to the blood-brain barrier, which restricts the infiltration of most antitumor drugs into the central nervous system, the standard treatment of GBM is limited to surgical resection followed by radiotherapy in combination with a chemotherapy (temozolomide, TMZ).3 However, surgical treatment for GBM is often compromised by the complexity of the intracranial operation and the dislodgement of tumor tissues. Residual tumor cells can lead to tumor recurrence within a relatively short time.4,5 Owing to the strong resistance of tumor cells to conventional therapies, including surgery, chemotherapy and radiotherapy, the median survival time of GBM patients with treatment is approximately only 12.5 months, and the two-year survival rate is less than 25%.6 Without any treatment, most GBM patients survive for only a few months.7

Recently, a variety of adoptive immunotherapies, including chimeric antigen receptor T cell immunotherapy (CAR-T), tumor-specific T cell receptors (TCR-T) and a multi-epitope-pulsed dendritic cell (DC) vaccine, have been used for the clinical treatment of GBM and showed major advantages.8–10 T cells can move through tissues, scan for MHC complexes and then activate their specific T cell receptors. In addition, tumor-specific T cells can be activated when encountering the TAAs of specialized antigen-presenting cells, including DCs.11 Moreover, a positive immune response to immunotherapy depends on dynamic interactions between tumor cells and immunomodulators inside the tumor microenvironment (TME).12 The TME is composed of tumor cells, stromal cells, inflammatory cell vasculature and extracellular matrices.13 Immunotherapies, which are capable of activating the immune system, expanding effector cells, infiltrating activated effector cells to the tumor tissue, and destroying tumor cells, exhibited successful tumor control.14 However, the TME usually prevents effective priming of lymphocyte, reduces there infiltration, and also suppresses infiltrating effector cells, leading to failure of immunotherapy.12 For example, PD-L1, the ligand of programmed cell death protein 1 (PD-1), can combine with T cells and inhibit their activation and then induce their exhaustion.15 Therefore, to improve the effectiveness of immunotherapies, especially those using specific TAAs as vaccine for personalized precision immunotherapies, a clear and more precise understanding of the expression of TAAs and TME gene in tumor cells is essential.9,10

With regard to the extreme short survival time of GBM patients, it is also very important to identify a method to accurately predict prognosis and to find an appropriate therapeutic scheme for patients. Currently, evaluation of prognosis for patients relies mainly on the clinical experience of doctors; a more comprehensive and exhaustive analysis is urgently needed. Recently, several gene signatures have been shown to predict prognostic outcomes. A previous study analyzed four data sets and identified a liver-specific, 7-gene signature that was correlated with a poor prognosis in Hepatocellular carcinomas (HCCs).16 Another study reported prognostic signatures derived from an optimized 5-gene platform to predict metastatic outcome independent of adjuvant chemotherapy use.17 Ng et al generated a 17-gene leukemia stem cell (LSC) score by extracting a list of genes differentially expressed in 78 acute myeloid leukemia (AML) patients and used this for analysis of five independent AML cohorts. The score was predictive of therapy resistance and patients with high LSC scores generally had a poor prognosis.18 However, the gene signatures for GBM remain to be elucidated. Therefore, it is imperative to develop a more objective and scientific evaluation method to predict the prognosis of GBM patients.

In this study, we designed specific primers for qRT-PCR to amplify 87 TAAs and 8 TME genes that were associated with tumorigenesis or had been used in clinical trials, and analyzed the expression levels of these genes in brain tumor tissues from 44 GBM patient tissues and 10 controls.19–28 We detected and quantified the mRNAs of these genes by qRT-PCR. All of the TAAs we selected induced immune responses, and some had already been used in the immunotherapy of other cancers.29–32 More importantly, we also analyzed the relationships between gene expression levels of these TAA/TME genes, clinical characteristics and OS using a linear regression method, and designed a system to predict prognosis of the patients. This may be helpful for designing clinical treatment and immunotherapy for GBM patients.

Materials And Methods

Patients And Tissues

Primary brain tumor tissues were obtained from 44 patients with stage IV GBM who underwent surgery at the Guangdong 999 Brain Hospital. None of these patients had received chemotherapy before surgery in this study. Normal brain tissues were obtained during surgical treatment of patients with non-malignant tumors or with trauma. Informed consent was obtained before collection of tissue samples in accordance with the Declaration of Helsinki and under the protocols approved by the Ethics Committee of the Guangdong 999 Brain Hospital. Clinic and pathological patient information is summarized in Table 1.

Table 1

Clinical Characteristics Of The GBM Patients In This Study

Variable	No.
No. of patients	44
WHO stage
Ⅳ	44
Gender
Male	31
Female	11
T & R*
Yes	23
No	21
Age
Mean	45.87
Range	6-67

Notes: *Standard chemotherapy with temozolomide (TMZ) & radiotherapy.

RNA Isolation And RT-PCR

Total RNA was isolated from 100 mg of tissue by using an RNeasy mini kit (Qiagen, Germany) according to the user’s manual. The cDNA was transcribed from 1 μg of RNA using a High-Capacity cDNA Reverse Transcription kit (Thermo Fisher, USA). Expression levels of the genes of interest were evaluated by gene amplification with specific primers (Table 2) on a CFX96 Real-Time system (Bio-Rad, USA) with human house gene glyceraldehyde-3-phosphate dehydrogenase (GAPDH) as a loading control, using SYBR kit (Thermo Fisher, USA). The expression value of GAPDH was used to normalize the quantification of the relative expression of the TAAs and TME genes. The relative expression level of the targeted genes in patients were compared with the average expression level of the targeted genes in normal brain tissue using the 2−ΔΔC formula.

Table 2

Primers Used For Amplifying TAAs And TME Genes

Primers For TAAs
NO.	Gene Symbol	Forward Primer 5ʹ-3ʹ	Reverse Primer 5ʹ-3ʹ
1	AIM2	GCCTCACGTGTGTTAGATGC	ATCTTCGGGGTTTCACCAGC
2	AKAP4	ATTCCATCAGCAAGGGGCTC	CTCCTTGGTGTGCCTTAGCA
3	ART4	GGAGGTGGTCACTGAGATTC	GCACGTATTCCGGTAAGG
4	BAGE	TGGCTCGTCTCACTCTGG	TCCTGTTGAGCTGCCGTCT
5	BCAN	GGAGGAGGCGACAAACTTC	GAGCTGTCTCCTTCCAGAACA
6	BSG	CCCTTCCTGGGCATCGT	CGGCGTCGTCATCATCC
7	CA9	GGACATATCTGCACTCCTGC	TGCTTAGCACTCAGCATCAC
8	CCND1	CCTCGGTGTCCTACTTCAAAT	CTCTTTTTCACGGGCTCCAG
9	PROM1	AGTGGCATCGTGCAAACCTG	CTCCGAATCCATTCGACGATAGTA
10	CDC45	GCAGGTGAAGCAGAAGTTCC	GCATGTCCTTCATCCCAAAT
11	NUF2	GAGAAACTGAAGTCCCAGGAAAT	CTGATACTTCCATTCGCTTCAAC
12	CEACAM5	TGTCGGCATCATGATTGG	GCAAATGCTTTAAGGAAGAAGC
13	CSPG4	CCTTTTGGGAGGCCCATGAT	GCAGCCTCAAAAGACACAGC
14	EPCAM	ACTACAAGCTGGCCGTAAAC	AGCCCATCATTGTTCTGGAG
15	EphA2	TCCCTGCTGTGCCATGCT	CCCTCAGCGGAAGTTGCA
16	EZH2	GGCCAGACTGGGAAGAAATC	ACCTCTTGAGCTGTCTCAGT
17	FABP7	AGCCTGGATGGAGACAAACT	TGCCTTCTCATAGTGGCGAA
18	FOSL1	CTGCCGCCCTGTACCTT	TGCTGCTACTCTTGCGATGA
19	GAGE1	TATGCGGCCCGAGCAGTT	CCTGCCCATCAGGACCATC
20	KCNMA1	GACATCACAGATCCCAAAAG	GTGTTGACGGCTGCTCATC
21	SLC1A3	CATCATTGCAGTGGACTGGTTTC	CCCATTTCAACATCTCGGTTCTTC
22	PMEL	ACAGGCCAACTGCAGAGG	CAGTTGGCGCCTGACCAG
23	MGAT5	TCAAAAGGCAGAACCAGTCC	GTGCTGGAGCCATAAACAGT
24	ERBB2	ATACCCTCTCAGCGTACCCTTGT	TCCGGAGAGACCTGCAAAGA
25	HBEGF	TTCTGGCTGCAGTTCTCTCG	AAGTCACGGACTTTCCGGTC
26	HNRNPL	TGGAGCAGAGGCAGCAG	TTTTGTGCGGGTCATCGTAG
27	HMOX1	AGTCTTCGCCCCTGTCTACT	CTTCACATAGCGCTGCATGG
28	TERT	CGTACAGGTTTCACGCATGTG	ATGACGCGCAGGAAAAATG
29	IGF2BP3	AGTTGTTGTCCCTCGTGACC	AGCCTTCTGTTGTTGGTGCT
30	IL13Ra2	GCAATGCACAAATGGATCAGAAG	TGCCAGGTTTCCAAGAACAGAGTA
31	IQGAP1	TGCTGAAGGACTCGTTGCAT	AGATTTCGGCGTTGGTCTGT
32	ITGAV	CGCTTCTTCTCTCGGGACTC	TCACATTTGAGGACCTGCCC
33	KIF1C	ACCGCACCAAGCAAATC	CTCCCTTCTTCCGTCTTCA
34	KIF21B	GTGAACCAGGACAAGACCAG	TGTAGCATGGCATTCTCTCG
35	KIFC3	CTGCGTAAGAAGTGCCACAA	AGGTGGATGATGGAGTCGTC
36	CTAG2	GTGTCCGGCAACCTACTGTT	CACATCAACAGGGAAAGCTG
37	LCK	AGTCAGATGTGTGGTCTTTTGG	CCTCCGGGTTGGTCATC
38	LRRC8A	AGGGAAAGGTGGGCTGCCTTT	ATACTGAAGAGGCAAGCTCCAG
39	MAGEA1	ACTGCAAGCCTGAGGAAGCC	TGGGTTGCCTCTGTCGAGTG
40	MAGEA10	TACTGCACCCCTGAGGAGGTC	TGTGGTGGCAATTCTGTCCTG
41	MAGEA2	ATGCCTCTTGAGCAGAGGAG	GAGCCCTCATCGGATTGTC
42	MAGEA3	GTCGTCGGAAATTGGCAGTAT	GCAGGTGGCAAAGATGTACAA
43	MAGEA4	CCACTACCATCAGCTTCACTTGC	CTTCTCGGAACAAGGACTCTGC
44	MAGEA6	GTCGTCGGAAATTGGCAGT	GCAGGTGGCAAAGATGTACAC
45	MLANA	gctcatcggctgttggtatt	CTGTCCCGATGATCAAACCC
46	MELK	GCCTGCCATATCCTTACTGG	AATCTCCGTTTTGATCCGGG
47	MET	CCATCCAGTGTCTCCAGAAGTG	TTCCCAGTGATAACCAGTGTGTAG
48	MUC1	AATGAATGGCTCAAAACTTGG	CACTAGGTTCTCACTCGCTCAG
49	NLGN4X	AGAATGCCTGCGGAACAAGA	TCCACGAACTTCAGGCCTTC
50	NrCAM	TTGTGCAAAGAGGGAGCATG	GGGCAGTTCCCTGTTGTCCT
51	ANKRD30A	ATCCTAGACTGGCTTCTGCT	ACAAGCATCTCCTGCAATGT
52	CTAG1B	TGTCCGGCAACATACTGACT	ACTGCGTGATCCACATCAAC
53	RPSA	CTGGTCTGAAGGTGTACAGGTGC	CTTAAGAGCCTATGCAAGAACAG
54	PCNA	TCTGAGGGCTTCGACACCTA	CATTGCCGGCGCATTTTAGT
55	PIK3R1	AACGAGTGGTTGGGCAATGA	CCTCGCAACAGGTTTTCAGC
56	PRAME	TCCAGAGCCAGAAGCAG	GGAACAGGTCTACGAGCA
57	PRKDC	ACCTGTTCTGGCAGGATGTC	TCTGAGGACGAATTGCCTTT
58	PTHLH	CCATCCAAGATTTACGGCGA	GGTGGTTCTTTGTGTTGGGA
59	PTPRZ1	ACCCCATCCTCCAGACAACA	GTAGCATGCAAGGCCGAATC
60	RPL19	CTCAGGCTTCAGAAGAGGCT	ATTGGCGATTTCATTGGTCT
61	SART1	AAGCAGCAGCAGGATTTC	TCCAGCAGCCCTTTGTTC
62	SART2	CCCTCTATGAAGGAGTTGCG	GGCCAAAGTGGTTGATGTTG
63	SART3	GAAATGTGCTGCCGTAGA	TGCTGACAAAGACGGTGA
64	SEC61G	GGACTCCATTCGGCTGGTTA	AGCAAATCCTATTGCTGTTGCC
65	SUGT1	CTGACTAAGGCTTTGGAACAGAA	CTGTAAAAGTTTCTAGGGCAGCA
66	SOX10	ATGCCAAAGCCCAGGTGA	TGAGGGAGGTGTAGGCGATC
67	SOX11	ACGGTCAAGTGCGTGTTTCTG	TGCTGGTGCGGTGGTTCCTC
68	SOX2	AAATGGGAGGGGTGCAAAAGAGGAG	CAGCTGTCATTTGCTGTGGGTGATG
69	SOX4	CGTCCTCAGATGACTTTCGG	TCTGGCACTTCCTTCAAACC
70	SPA17	GCTCGGAGAGAAAGGAGGTTC	TACTCCCCCATTCTGCTGGA
71	SPAG9	AGTCATCAGCCCACAAAGTAGCAG	GATTCTCCACCTTCATCACCCATT
72	SPANXB1	TAGTGGTTCGCTACAGGAGGAACGTGA	TTGCCGAAGTTTGAGGGATGTAG
73	STAT3	CCAAGCGAGGACTGAGCATC	CCAGACCCAGAAGGAGAAGC
74	BIRC5	ACTGAGAACGAGCCAGACTT	CGGACGAATGCTTTTTATGTTC
75	TRIOBP	GCCATGACGCCCGATCTG	AGGTGGTGGTGAGCGAGG
76	T/Brachyury	CGCTTCAAGGAGCTCACCA	CGAAGTCCAGCAGGAAGGAG
77	TNC	TGGCATCGGAGAATGCCTTT	CAGCTTCCTCTGGGTTCCTG
78	5T4	GCGGACCCCAGATTAACAAAC	GTGTGGGTACACTTGCTACACC
79	CSAG2	AGTAGACTGTTGAGAGACGCT	TCCACTTCCTCGCCTCTTTG
80	PRSS1	TGCCCCCTTTGATGATGATG	CTGATACCACCCACTGTTCG
81	DCT	CCTGTCTCTCCAGAAGTTTG	CAGAGTCCCATCTGCTTTATC
82	UBE2V1	TCTAATGGAGTGGTGGACCC	CTGTAACACTGTCCTTCGGG
83	NELFA	AACGCCCTGACGACCCT	CGCTCCGCTTCAACTGC
84	WT1	GATAACCACACAACGCCCATC	CACACGTCGCACATCCTGAAT
85	XAGE-1b	TGGATTCTTTCTCCGCTACTG	AAACCAGCTTGCGTTGTTTC
86	CHI3L2	TTGACTGTGGGCGTATC	AGAGGGCTGTTGTGGC
87	CHI3L1	AACGATCACATCGACACCTG	TTGAGACCCAAAGTTCCATC

Based on previous research, we selected 87 TAAs and 8 TME genes that were associated with tumor oncogenicity, proliferation and metastasis or had been used in clinical trials. Specific primers for optimal amplification of these genes were designed using Primer 5 and DNAMAN software and are listed in in Table 2. The specificity of the TAAs and TME genes amplified by RT-PCR using these primers was confirmed by sequencing analysis of the PCR products.

Data Analysis, Signature Training And Statistical Analysis

Heat maps were generated using the MultiExperiment Viewer 4.9.0 software. To code the qRT-PCR data for use with the software, the gene expression levels of the 87 TAAs and 8 TME genes in the GBM tumor tissues and in the normal brain tissues were compared. When the expression level of a given gene in GBM tumor tissue was lower than, equal to, or higher than the mean value of the expression level of the corresponding gene in normal brain tissue, it was defined as < 0, =0, or >0, respectively.

Based on the overall median survival time (12.5 months), 44 GBM patients were divided into two groups depending on whether their survival was <12.5-months or >12.5-months. For both clinical characteristics and gene signatures, we used a linear regression method based on the LASSO algorithm as executed in the R package, and the leave-one-out cross-validation method was used to fit a Cox regression model.18 The clinical characteristics used in the formula included the NLR, the number of eosinophils (EOS) and BAS. GBM patients receiving standard chemotherapy and radiotherapy treatment were assigned a value of 1, otherwise, those patients were assigned a value of 0. The gene expression score was calculated as a linear formula weighted by regression coefficients as follows:

Y1 (Gene expression (GE) score for TAAs) = 0.153*CHI3L1-0.167*EZH2-0.075*PCNA-0.141*PIK3R1-0.046*PRKDC-0.004*SART3-0.121*EPCAM;

Y2 (GE score for TAAs and TME genes) = 0.143*CHI3L1-0.165*EZH2-0.020*PCNA

-0.163*PIK3R1-0.013*PRKDC-0.125*EPCAM-0.099*FOXP3;

Y3 (GE score for TAAs and TME genes, age and treatment) =0.049*CHI3L1-0.133*EZH2-0.066*TRIOBP-0.098*PIK3R1-0.008*PRKDC-0.066*SART3-0.102*EPCAM+0.116*Age-0.131*Treatment;

Y4 (GE score for TAAs and TME genes, NLR and BAS) =0.146*CHI3L1-0.201*EZH2-0.038*TRIOBP-0.092*PIK3R1-0.017*PRKDC-0.110*EPCAM-0.013*FOXP3+0.269*NLR+0.256*BAS; and

Y5 (GE score for TAA and TME genes and all clinical characteristics) = 0.110*CHI3L1-0.150*EZH2-0.066*TRIOBP-0.098*PIK3R1-0.008*PRKDC-0.102*EPCAM+0.237*NLR+0.155*BAS+0.116*Age-0.131*Treatment.

The statistical analysis was performed by using SPSS 19 and GraphPad Prism 7.0. Two-tailed t-tests were used to evaluate the correlation between patients and the other elements. OS was defined as the time from diagnosis of GBM until the patient’s death or the end of follow-up. Univariate analysis and multivariate survival analysis were performed using Kaplan-Meier and Cox regression, respectively.

Results

Expression Of TAAs And TME Genes In GBM Tumor Tissues

We investigated the mRNA expression levels of a panel of 87 TAAs and 8 TME genes in tumor tissues from 44 GBM patients in comparison with 10 normal brain tissues (Figure 1). A two-fold increase in mRNA expression level of a TAA relative to the mean expression level in the 10 normal brain tissues was defined as a positive result. A total of 14 TAAs were identified as positive in GBM tumors at the population level (Figure 1A and B), while the other tested TAAs showed no consistent increase or decrease in GBM patients compared with normal brain tissues. We also evaluated the expression levels of 8 TME genes (Figure 1C and D) and found significantly higher expression of IDO1 and PDL-2 in GBM patients than in normal brain tissues.

Figure 1

Relative expression levels of 87 TAA and 8 TME genes in tumor tissues of 44 GBM patients. Expression levels of each of the 87 TAAs and 8 TME genes in GBM patients were quantified by qRT-PCR and compared with the relevant gene expression level averaged from 10 normal brain tissues. Abbreviations on the y axis of the (A–D) of indicate the individual genes tested. GAPDH was used as the reference gene in this study. Fold change in gene expression level is indicated on the x axis in (A and C). Individual tumor tissues from 44 GBM patients are indicated on the x axis for (B and D). Gene expression levels of 87 TAAs are shown in (A and B), and expression levels of 8 TME genesis shown in (C and D). The red dotted line on the left in (A and C) indicates a two times lower gene expression level compared with the average level in normal tissues. The red dotted line in the middle indicates that the level of gene expression is the same as the average level in normal tissues. The red dotted line on the right indicates a two times higher gene expression level compared with the average level in normal tissues. <0, gene expression level lower than the average level in normal tissue (blue in the heat map). 0, the gene expression level the same as the average level in normal tissue (white in the heat map); >0, gene expression level higher than the average level in normal tissue (red in the heat map).

Clinical Characteristics And Gene Expression Correlated With The OS Of GBM Patients

We next examined the contribution of clinical characteristics and gene expression to OS in 44 GBM patients by univariate analysis. Significant correlations with OS were found for age (p=0.0339, hazard ratio (HR): 1.6161, 95% CI: 1.0372–2.5179), postoperative standard of radiotherapy and chemotherapy (p=0.0082, HR: 0.3347, 95% CI: 0.1488–0.7532), NLR (p=0.0003, HR: 2.8430, 95% CI: 1.6037–5.0399), and BAS (p=0.0152, HR: 2.0712, 95% CI: 1.1509–3.7273) (Table 3). A total of 13 TAA genes with increased expression levels were significantly correlated with OS in GBM patients (Table 3). FOXP3 was the only gene among the 8 TME genes tested for which increased expression was significantly correlated with OS in GBM patients (p=0.0129, HR: 0.7544, 95% CI: 0.6042–0.9420).

Table 3

Correlation Of Clinical Characteristics And Gene Expression Scores With The OS Of 44 GBM Patients

				95% CI
Clinical Characters Or Gene Symbol		P Value	HR	Low	High
Clinical Features	Age	0.0339	1.6161	1.0372	2.5179
	Treatment*	0.0082	0.3347	0.1488	0.7532
	NLR#	0.0003	2.8430	1.6037	5.0399
BAS&	0.0152	2.0712	1.1509	3.7273
TAAs	SURVIVIN	0.0274	0.8114	0.6739	0.9770
	BSG	0.0153	0.6545	0.4648	0.9218
	CDC45	0.0420	0.8152	0.6525	1.0184
	EZH2	0.0209	0.7302	0.5592	0.9535
	MELK	0.0134	0.8274	0.7121	0.9615
	NELFA	0.0220	0.7942	0.6521	0.9673
	PCNA	0.0033	0.5206	0.3368	0.8046
	PIK3R1	0.0054	0.6801	0.5183	0.8925
	PRKDC	0.0022	0.6569	0.5023	0.8592
	SART3	0.0054	0.6219	0.4450	0.8691
	SPAG9	0.0124	0.6288	0.4371	0.9047
	STAT3	0.0153	0.6129	0.4125	0.9105
TRIOBP	0.0435	0.7635	0.5876	0.9922
TME gene	FOXP3	0.0129	0.7544	0.6042	0.9420

Notes: *Standard radiotherapy and chemotherapy after surgery. #Neutrophil to lymphocyte ratio. &Number of blood basophils.

Gene Expression Score Correlated With The OS Of GBM Patients

The median survival time of the GBM patient group with the OS <12.5 months and >12.5 months were plotted using Kaplan-Meier analysis in SPSS 19 (Figure 2A). We evaluated and correlated each of the clinical characteristics and gene signatures with the OS of the two GBM patient groups as described in the Material and Methods section. The overall scores (Table 4) evaluated based on the formulas (Y1-Y5) were significantly different between the two groups, demonstrating the reliability of these formulas (Figure 2B–F). When expression levels of TME genes and clinical features were added one by one to the Y1 formula, there was an increased trend in sensitivity, specificity and accuracy (Table 5).

Figure 2

Correlation of gene expression scores in different models with the OS of GBM patients. The OS rate of the 44 GBM patients was plotted using Kaplan-Meier analysis (A). The median survival (12.5 months) was calculated and used to divide the patients into two patient groups with either <12.5 months or >12.5 months survival time. Five different models (Y1–Y5) using different combination of gene expression scores and clinical characteristics were used to examine the correlations between the two patient groups presented in (B–F). P values were calculated using the student’s t-test. *** indicates P <0.001 and **** indicates P <0.0001.

Table 4

Relative Index Of Different Models Based On The Correlation Between The Elements And The OS Of The 44 GBM Patients

Clinical Characters Or Gene Symbol	Relative Index
Y1	Y2	Y3	Y4	Y5
CHI3L1	0.153	0.143	0.049	0.146	0.110
EZH2	−0.167	−0.165	−0.133	−0.201	−0.150
TRIOBP	–	–	−0.066	−0.038	−0.066
PCNA	−0.075	−0.020	–	–	–
PIK3R1	−0.141	−0.163	−0.098	−0.092	−0.098
PRKDC	−0.046	−0.013	−0.008	−0.017	−0.008
SART3	−0.004	–	−0.066	–	–
EPCAM	−0.121	−0.125	−0.102	−0.110	−0.102
FOXP3	–	−0.099	–	−0.013	–
NLR#	–	–	–	0.269	0.237
BAS&	–	–	–	0.256	0.155
Gender	–	–	0.116	–	0.116
Treatment*	–	–	−0.131	–	−0.131
Intercept	0.050	0.050	0.050	0.050	0.050

Notes: #Neutrophil to lymphocyte ratio. &Number of blood basophils. *Standard radiotherapy and chemotherapy after surgery. Y1: Gene expression score of the indicated TAAs. Y2: Gene expression score of the indicated TAAs and TME gene. Y3: Gene expression score of the indicated TAAs and TME gene and clinical characteristics of gender and treatment. Y4: Gene expression score of the indicated TAAs and TME gene and clinical characteristics of NLR and BAS. Y5: Gene expression score of the indicated TAAs and TME gene and all clinical characteristics.

Table 5

Validation Of The Gene Expression Score Formulas (Y1-Y5)

	Significant coefficients*
Y1	Y2	Y3	Y4	Y5
Sensitivity	0.762	0.762	0.809	0.905	0.905
Specificity	0.684	0.737	0.631	0.737	0.842
Accuracy	0.725	0.750	0.725	0.825	0.875

Notes: *The discrimination coefficients of different models between two patient groups based on the median survival (12.5 months). The larger of the value, and the better of the formula. Y1: Gene expression score of the indicated TAAs. Y2: Gene expression score of the indicated TAAs and TME gene. Y3: Gene expression score of the indicated TAAs and TME gene and clinical characteristics of gender and treatment. Y4: Gene expression score of the indicated TAAs and TME gene and clinical characteristics of NLR and BAS. Y5: Gene expression score of the indicated TAAs and TME gene and all clinical characteristics.

To verify the sensitivity, specificity and accuracy of the gene expression score (Y1-Y5), we calculated gene expression scores for the 44 GBM patients individually, and grouped patients into high and low scoring groups based on the median score. The percentage of surviving GBM patients was significantly different (P<0.05, log rank test) between the high and low scoring groups with all of the five gene expression score formulas (Figure 3).

Figure 3

Correlation of the OS of the 44 GBM patients with high or low gene expression scores. Gene expression scores were calculated based on the level of gene expression of TAAs and TME genes quantified by qPCR as described in Material and Methods. Based on the level of gene expression, the 44 GBM patients were divided into low and high gene expression groups. Correlation of the percentage (on the y axis) of GBM patients with low (green curve) or high (red curve) gene expression scores with survival over time (x axis) was evaluated using 5 gene expression scoring models (Y1–Y5; A–E). P values were calculated using the log rank test and are indicated in the individual plots.

Survival Analysis Of Patients Using TCGA And GEO Databases By Gene Expression Score (Y1-Y3)

Furthermore, to verify the applicability, sensitivity, specificity and accuracy of the formulas (Y1-Y3), gene expression scores were validated against published clinical GBM cohorts from the TCGA (Nature, 2008, n=527, Provisional, n=205) and GEO (GSE4412, n=84).33,34 As no information on NLR, EOS or BAS was available in these databases, we evaluated patients using only the Y1-Y3 formulas. Patients were again divided into high and low scoring groups with respect to gene expression, based on the median scores using the same method as described above (Figure 4). Again, we found significant differences between the two groups for each of the three different databases, as calculated by formulas Y1-Y3, with P values of 0.0033, 0.0018, and 0.0042 for patients in the TCGA (Nature, 2008) data set; 0.0399, 0.0294, and 0.0001 for patients in the TCGA (Provisional) data set; and 0.0139, 0.0095, and 0.0019 for patients in the GSE4412 data set.

Figure 4

Correlation of the OS of GBM cohorts in the TCGA and GEO databases (Nature, 2008, Provisional and GSE4412) with low and high gene expression scores. (A-C), Kaplan-Meier evaluation of OS in the TCGA database Nature, 2008 based on gene expression scores (Y1–Y3); D-F and G-I data from the TCGA (Provisional) and GSE4412 databases, respectively. For all panels, the two groups with scores lower and higher than the median value in (A–C) are indicated by green and red lines, respectively. P values were calculated by using the log rank test, and are indicated in the individual plots.

Discussion

In the present study, we first evaluated the expression levels of 87 TAAs and 8 TME genes in tumor tissues of 44 GBM patients compared with 10 normal tissues. We also established linear risk scores as survival prediction models based on the expression levels of the genes of interest and clinical characteristics for prediction of the prognosis of GBM patients.

Owing to the strong resistance of GBM to conventional therapies such as surgery, chemotherapy and radiotherapy, the median survival time of GBM patients with treatment is approximately only 12.5 months.35 In recent years, an increasing number of immunotherapies targeting human GBM and other solid cancers have been developed. CAR-T cells were generated from patients’ T cells using lentiviral transfection to introduce specific TAAs, leading to cell killing within a short time.36 Various of vaccine based immunotherapies, including DC based vaccines, autologous and allogeneic antigens vaccines, peptides vaccines and viral based vaccines, and the vaccine pulsed with specific TAAs were infused into patients and shown to stimulate autologous anti-tumor immune responses.28,36 The question remained how to predict the prognosis of patients in order to provide better and more effective treatment for GBM patients in such a short time. This study investigated whether prevalent and concomitant patterns of TAAs and TME genes expression in tumor tissues and clinical features of GBM patients could be used not only for prediction of prognosis but also for the design of cocktail immunotherapies (such as a multi-epitope-pulse DC vaccine).37 We determined the gene expression levels of 87 TAAs and 8 TME genes by qRT-PCR. All TAAs and TME genes selected in this study have been reported to be expressed in brain tumors and to induce a series of immune responses in vivo.26,30,38 Moreover, most, if not all, of these genes have already been used in clinical trials of immunotherapy.38,39 We identified 14 TAAs (CHI3L1, CHI3L2, BIRC5, TNC, MELK, CDC45, IGF2BP3, IL13Rα2, NUF2, SOX2, SOX11, HMOX1, EZH2 and FOSL1) with increased gene expression levels of in all GBM tumor tissues; these all have key roles in tumorigenesis, development, invasion and migration of brain tumors. For example, a bifunctional inhibitor of apoptosis protein BIRC5 is highly expressed in many human malignancies including GBM, where it plays important role in the proliferation, drug resistance and anti-apoptosis of cancer cells, and is correlated with the decrease OS.40 BIRC5 is also an ideal target for immunotherapy, and there are more than 60 ongoing or successfully completed clinical trials targeting BIRC5 listed on the Clinical trials website. One clinical study of a peptide vaccine targeting BIRC5 showed partial and even complete remission in participants.19 These TAAs that are highly expressed in almost all GBM patients could be used to design a cocktail DC vaccine. We also found individual differences in expression of some TAAs. For instance, TAAs such as IL13Rα2 was not highly expressed in all GBM patients, although it had been studied in a number of clinical trials involving GBM.41 Brown carried out CAR-T cell immunotherapy targeting IL13Rα2, and some patients showed anti-glioma responses.42 These results indicate that all TAAs that are highly expressed in individuals are of potential importance in clinical immunotherapy, and may have key roles in tumorigenesis, development, invasion and migration.

However, downregulation of immune responses mediated by the TME greatly decreases the effects of immunotherapy.43 In this study, we found that almost of all the GBM tumor tissues had elevated gene expression levels of TME genes (including IDO, TDO, PDL-1, COX2 and FOXP3). These results were consistent with previous observations that these genes have a key role in immune escape, invasion and angiogenesis.15 Thus, if these immunosuppression factors were neglected in immunotherapy, the treatment could be negative affected.

In addition, we examined several clinical features, including as age, gender, treatments, NLR, EOS and BAS. Correlation analysis of these clinical features with the OS of 44 GBM patients suggested that age, chemotherapy, radiotherapy and NLR are important prognostic factors in GBM patients. Similar observations have also been made in different malignant tumors.44–46

The diagnosis of GBM is currently based on the clinician’s experience and judgment, however this is often inaccurate and might change with the patient’s physical condition. The combination of gene expression levels and clinic factors may improve prediction accuracy, and have been used to identify a higher risk of recurrence and death.47 We designed mathematical models using a linear regression method based on the LASSO algorithm, starting with the gene expression levels of TAAs and adding expression levels of TME genes and clinical features one by one to optimize the models. In this way, we developed 5 models (Y1-Y5) with a cut-off value of 0.05 showing improvements with respect to sensitivity, specificity and accuracy (Table 4). These models were further validated using the relevant data in the TCGA and GEO databases suggesting that these formulas could be used objectively and accurately to predict prognosis of patients based on their gene expression scores.

Conclusion

In summary, our study established prognostic prediction models based on a full understanding of gene expression profiles that provides an accurate method for survival prediction and guidance for implementing better treatment strategies. The outcomes of this study will also benefit future personalized prediction and precision immunotherapy for GBM management.