Gene expression-based biomarkers designating glioblastomas resistant to multiple treatment strategies.

Despite advances in molecular characterization of glioblastoma multiforme (GBM), only a handful of predictive biomarkers exist with limited clinical relevance. We aimed to identify differentially expressed genes in tumor samples collected at surgery associated with response to subsequent treatment, including temozolomide (TMZ) and nitrosoureas. Gene expression was collected from multiple independent datasets. Patients were categorized as responders/nonresponders based on their survival status at 16 months postsurgery. For each gene, the expression was compared between responders and nonresponders with a Mann-Whitney U-test and receiver operating characteristic. The package 'roc' was used to calculate the area under the curve (AUC). The integrated database comprises 454 GBM patients from 3 independent datasets and 10 103 genes. The highest proportion of responders (68%) were among patients treated with TMZ combined with nitrosoureas, where FCGR2B upregulation provided the strongest predictive value (AUC = 0.72, P < 0.001). Elevated expression of CSTA and MRPS17 was associated with a lack of response to multiple treatment strategies. DLL3 upregulation was present in subsequent responders to any treatment combination containing TMZ. Three genes (PLSCR1, MX1 and MDM2) upregulated both in the younger cohort and in patients expressing low MGMT delineate a subset of patients with worse prognosis within a population generally associated with a favorable outcome. The identified transcriptomic changes provide biomarkers of responsiveness, offer avenues for preclinical studies and may enhance future GBM patient stratifications. The described methodology provides a reliable pipeline for the initial testing of potential biomarker candidates for future validation studies.

Introduction

Glioblastoma multiforme (GBM) is the most common brain tumor among adults making up 47.7% of primary malignant tumors of the brain (1). As part of WHO grade IV gliomas, GBMs are almost uniformly fatal (2). Despite aggressive multimodal treatment, the median survival reaches only 14–17 months in contemporary clinical trials, 12 months in population-based studies (3) and the 5-year survival rates are as low as 5.6% (1).

For newly diagnosed patients, the goal is gross surgical resection with preservation of neurological functions (4). Postoperative therapy depends on the molecular characterization of tumor specimens. However, despite all advances achieved by studying methylation patterns, gene mutations, copy number variations and oncogenic pathway activations (5–7), mechanisms of GBM pathogenesis remain poorly understood. Implementation of molecular genetics into GBM classification is limited to a handful of prognostic or predictive biomarkers with restricted relevance to treatment success (8). A subtype of gliomas with distinct biological and clinical features carries isocitrate dehydrogenase 1 and 2 (IDH1/IDH2) mutations linked to improved prognosis, although without significance in routine clinical practice (9,10), and more than 90% of glioblastomas are IDH-wild type. TERT-promoter mutations are associated with worse survival but are independent of other clinical and molecular factors (11).

The methylation status of the O6-methylguanine-DNA methyltransferase (MGMT)-promoter has implications for therapy response. The therapeutic benefit of an alkylating agent, temozolomide (TMZ), depends on its ability to alkylate/methylate DNA at the N-7 or O-6 positions of guanine residues that triggers cell death. Nevertheless, some tumor cells express a repair enzyme encoded by MGMT. Epigenetic silencing of this gene prevents enzyme synthesis, and such tumors are more sensitive to the effects of TMZ; consequently, MGMT silencing and/or loss of protein expression correlates with improved progression-free and overall survival (OS) in GBM patients after TMZ therapy (12).

In the open-label, randomized, phase III EORTC-NCIC trial with 573 patients, the addition of TMZ to radiation therapy improved median OS compared with radiation therapy alone (14.6 versus 12.1 months, hazard ratio 0.63, 95% confidence interval 0.53–0.75) (13), and survival continued to be superior years after treatment (27.2 versus 10.9% at 2 years and 9.8 versus 1.9% at 5 years) (14). In a subset of patients, MGMT silencing was an independent favorable prognostic factor associated with a median survival of 21.7 months (12), although this has not been confirmed in a prospective manner (12,15). Based on these results, TMZ was approved in conjunction with radiotherapy, followed by adjuvant TMZ in 2005 by the FDA as a new standard of care for newly diagnosed adult glioblastomas (16).

Despite limited benefits for patients with unmethylated tumors, current recommendations for newly diagnosed patients with good performance status up to 70 years with both MGMT-methylated and unmethylated tumors suggest postoperative radiotherapy with concurrent and adjuvant TMZ. For the latter cohort with unknown MGMT-status, better alternatives are being continuously investigated (17); nevertheless, the development of new therapies is time-consuming and challenging (18) as GBMs are wildly heterogeneous (6).

Additionally, not all patients with methylated MGMT-promoter respond initially to TMZ and most of those who respond acquire resistance over time (12). Combination therapies involving bevacizumab in the up-front setting appeared to be promising initially, but benefits were not confirmed in phase III trials: the addition of bevacizumab to radiotherapy—TMZ did not improve OS while increasing the risk of toxicity (19–21). A recent, promising approach compared the lomustine/TMZ plus radiotherapy protocol to the standard of care in newly diagnosed GBM patients with MGMT-silenced tumors in the randomized, phase III CeTeG/NOA-4 trial. The combined treatment improved median OS from 31.4 (95% confidence interval 27.7–47.1) to 48.1 months (22), although the sample size was small and the study had further limitations. Additional clinical trials have been designed to target the most frequent molecular alterations, including the PI3K, EGFR, AKT, mTOR, p53 and RB pathways, usually in combination with TMZ; however, responses in newly diagnosed GBM were either absent or not durable (8); thus, core pathway alterations are challenging for drug design. Novel pathways for molecular targeted therapies and more efficient biomarkers of therapy response are much sought after.

We embraced a novel approach and linked gene expression in GBM patients with response to therapy across the entire transcriptome. For this, we assembled a substantial GBM database and ought to identify differentially expressed genes in tissue samples collected at the surgery in subsequent responders and nonresponders to chemotherapy. Since the majority of GBM patients receive TMZ regardless of their MGMT-status, we focused on upregulated genes associated with TMZ treatment, including patients who received TMZ as a single agent or combined with other substances, such as nitrosoureas. We aimed to identify predictive markers of therapy success, novel targets for personalized treatment strategies, and potential biomarker candidates for the identification of patients who could benefit from participating in clinical trials.

Materials and methods

Database construction

We searched the GEO (http://www.pubmed.com/geo) and the TCGA (http://cancergenome.nih.gov) databases for suitable patient cohorts by using the keywords ‘glioblastoma’, ‘GPL96’, ‘GPL570’ and ‘GPL571’. These three GEO platforms (GPL) refer to three Affymetrix gene array platforms that share identical probe sets to measure gene expression. Non-overlapping probe sets were not utilized in our analysis. Only publications with available raw microarray gene expression data, clinical treatment and response or survival information, and at least 20 patients were included.

The raw CEL files were MAS5 normalized using the Affy Bioconductor library in the R statistical environment (http://www.r-project.org) (23). A second scaling normalization was performed to set the mean expression on each chip to 1000 to reduce batch effects (24).

Clinical data

For each sample corresponding clinical data were retrieved manually, and each dataset was validated by both JTF and BG to ensure the reliable designation of clinical characteristics.

We focused on the association between initial gene expression and subsequent survival in patients treated with various treatment regimes. Because data on RFS were not available for more than a third of patients (35.4%), subsequent analyses were focused on OS. Based on whether the patient was alive at a determined cutoff time (16 months postsurgery, corresponding to the median survival, 95% confidence interval 14.72–17.28), we divided patients into responders or nonresponders for each treatment. In our sample population, the age distribution was relatively wide (ranging from 10 to 89 years). Therefore, we also investigated gene expression differences between responders and nonresponders among younger and older patients (with a cutoff of 55 years of age) separately.

Statistical analysis

Gene expression was accessible for 10 103 unique genes. For each gene, the expression was compared between responders and nonresponders for each individual treatment type with Mann–Whitney U-test and receiver operator characteristics.

Statistical significance was accepted in case of P < 0.05 and fold change (FC) >1.44. FC was calculated by computing the ratio of gene expression between the median expressions in the responder and nonresponder cohorts. In addition, to avoid noise due to genes with very low expression, the mean expression difference between responders and nonresponders was required to reach at least 350. The package ‘roc’ was used to calculate the area under the curve (AUC) and significance (http://www.bioconductor.org) (25). To correct for multiple testing, we calculated the false discovery rate (FDR), and only genes with FDR <10% were considered significant. A chi-square test was performed to examine the relationship between MGMT-quartiles and the number of responders to TMZ therapy and also to compare the number of responders and nonresponders among younger (<55 years of age) and older (≥55) patients.

For the selected genes, survival analysis was performed in the same set of patients by performing Cox proportional hazards regression using all cutoff values between the lower and upper quartiles of expression. Kaplan–Meier plots were drawn to visualize survival differences. For the selection of the strongest biomarker candidates, only genes performing at FDR ≤10% were considered significant in the survival analysis.

Finally, biological processes related to functionally associated gene groups were assessed by gene enrichment analysis performed in the Database for Annotation, Visualization and Integrated Discovery (DAVID) Bioinformatics Resource 6.8 (26).

Results

Transcriptomic database for biomarker selection

Three datasets were eligible for our analysis, and these contain 454 patients diagnosed with GBM. The samples stem from the GEO datasets GSE7696 (27) and GSE108474 (28) and from the TCGA repository (Figure 1A).

Figure 1.

Analysis workflow for the database setup (A). Distribution of patients undergoing various treatment protocols (B) and response to the applied treatment combinations (C). Abbreviations: TMZ-only: patients undergoing single-agent TMZ therapy. TMZ combined: patients treated with TMZ combined with other agents, including nitrosoureas, topoisomerase- or angiogenesis-inhibitors. TMZ + topo-angio: patients receiving TMZ combined only with topoisomerase- or angiogenesis-inhibitors. TMZ + Nitroso-only: patients treated with TMZ combined only with nitrosoureas. TMZ + Nitroso-combined: patients treated with TMZ combined with nitrosoureas (with or without topoisomerase- or angiogenesis-inhibitors). Nitroso-only: patients undergoing treatment with single-agent nitrosoureas. Nitroso without TMZ: patients treated with nitrosoureas with or without topoisomerase- or angiogenesis-inhibitors, excluding TMZ. Topo-angio: patients receiving only topoisomerase- or angiogenesis-inhibitors. Other chemo: patients receiving chemotherapy other than TMZ, nitrosoureas, topoisomerase- or angiogenesis-inhibitors.

Patient distribution was balanced across age cohorts: about half (46.7%) of patients were younger than 55 years of age (<55) while 242 were at least 55 years old at the time of diagnosis (≥55), the entire range spanning between 10 and 89 years (Table 1).

Table 1.

Clinical characterization of the entire GBM dataset and detailed information for patients receiving TMZ in any combination (TMZ-all)

Characteristics	All GBM patients	TMZ-treated patients
Age
Under 55	46.7%	44.5%
55 or older	53.3%	55.5%
Sex
Male	65.6%	67.1%
Female	34.4%	32.9%
OS
(Median OS)	16 months	16.3 months
Dead	92.5%	90.9%
Recurrence-free survival
(Median RFS)	6.9 months	10.7 months
No data	32.2%	28.2%
Treatment
Chemotherapy	100%	100%
Radiation therapy	91.4%	95.3%
TMZ	70.3%	100%
Bevacizumab	12.8%	16.9%
Irinotecan	10.6%	12.2%
Lomustine	11.9%	9.1%
Carmustine	18.1%	13.8%
Dexamethasone	11.0%	8.8%
Any topoisomerase-inhibitor	13.7%	15.4%
Any angiogenesis-inhibitor	25.3%	27.6%
Any nitrosourea	30.0%	23.2%

Nearly two-thirds of patients (65.6%) were males. All patients received some sort of chemotherapy, and the majority of patients (91.4%) also underwent radiation therapy. There were only 15 patients who did not receive radiation therapy but were treated with TMZ.

Altogether 319 patients received TMZ (TMZ-all), out of which 189 patients were treated with TMZ as a single agent (TMZ-only), and 130 patients received TMZ in combination with other therapies, including nitrosoureas, topoisomerase- and angiogenesis-inhibitors (TMZ combined). Out of these 130 patients, 56 were treated with TMZ combined with topoisomerase-inhibitors and/or angiogenesis-inhibitors, but not nitrosoureas. Seventy-four patients received nitrosoureas along TMZ with or without other agents (TMZ + Nitroso-combined), out of which 29 patients were treated with TMZ combined only with nitrosoureas (TMZ + Nitroso-only). One hundred and thirty-five patients received chemotherapy treatment, not including TMZ (Supplementary Table 1, available at Carcinogenesis Online).

Altogether 30% of all patients (n = 136) were treated with one of the following nitrosoureas: lomustine, carmustine, fotemustine, estramustine, laromustine or nimustine, with or without TMZ or combined with other therapies. Out of the 135 patients not treated with TMZ, 48 received nitrosoureas only (Nitroso-only) and 14 patients were treated with nitrosoureas combined with other agents. These two groups were combined as Nitroso without TMZ for subsequent analyses (n = 62) (Figure 1B).

A small portion of patients (13.7%, n = 62) received topoisomerase-inhibitors, including irinotecan, etoposide, topotecan and teniposide, combined with other therapies including TMZ, nitrosoureas or in the case of three patients as a single agent (Supplementary Table 1, available at Carcinogenesis Online). Finally, 25.3% (n = 115) of patients underwent therapy with angiogenesis-inhibitors, including bevacizumab, thalidomide, vandetanib, vatalanib and dexamethasone, with or without TMZ or other therapies, or in 19 cases as single agents (Supplementary Table 1, available at Carcinogenesis Online). Altogether 24 patients received topoisomerase-inhibitors or angiogenesis-inhibitors without TMZ or nitrosoureas (Topo-angio), and 49 patients were treated with chemotherapies other than TMZ nitrosoureas, topoisomerase-inhibitors or angiogenesis-inhibitors (Other chemo). Due to the low number of patients in these treatment combinations, separate analyses were not pursued.

Treatment-dependent response rates

Out of all 454 patients, the 49.1% alive at the 16 months cutoff were assigned to be responders. Among the 189 patients receiving single-agent TMZ, 45.5% were alive within the 16 months cutoff compared with 60.8% of patients treated with TMZ combined. Response rates in the patient cohort treated with TMZ combined only with nitrosoureas (TMZ + Nitroso-only) were 62.1%, while the highest portion of responders were in the treatment category TMZ + Nitroso-combined (68%).

Response rates were lowest among patients receiving single-agent nitrosoureas (Nitroso-only, 43.8%), and similarly low among patients receiving Nitroso without TMZ (45.2%) (Figure 1C). Response rates were in a similar range among patients treated with chemotherapies other than TMZ, nitrosoureas, topoisomerase-inhibitors or angiogenesis-inhibitors (42.9%). Response rates for each treatment category are summarized in Table 2.

Table 2.

Response rates to different chemotherapy regimes

Treatment	Number of patients	% Responders (OS >16 months)
TMZ-all (patients receiving TMZ in any combination)	319	51.7
TMZ-all, younger age cohort (<55)	142	63.4
TMZ-all, older age cohort (≥55)	177	42.4
TMZ-all in a cohort with low MGMT-expression (Q1)	80	65.0
TMZ-all in a cohort with high MGMT-expression (Q2–Q4)	239	47.3
TMZ-only	189	45.5
TMZ combined	130	60.8
TMZ + Nitroso-only	29	62.1
TMZ + Nitroso-combined	74	68.0
Nitroso-only	48	43.8
Nitroso without TMZ	62	45.2
Other chemotherapy	49	42.9

MGMT-expression and age-dependent responses to chemotherapy

MGMT-methylation has repeatedly been described as a prognostic marker in GBM patients treated with alkylating agents, such as TMZ (12,15). MGMT-promoter methylation has been associated with very low MGMT transcript levels (29,30). With no available data on tumor methylation, we used the expression of MGMT (O-6-methylguanine-DNA methyltransferase) as a surrogate marker of gene activity. Based on MGMT-expression, we divided the TMZ-treated population (TMZ-all) into quartiles, and we compared the proportion of responder/nonresponder within each quartile. Almost twice as many responders (52 versus 28) were in the lowest quartile (Q1) compared with nonresponders. Response rates were evenly proportioned in the second quartile (Q2), while in the third (Q3) and fourth quartiles (Q4) nonresponders surpassed the number of responders (Figure 2A). Given the similarity of therapy response across Q2–Q4 we merged data in subsequent analyses. A chi-square test revealed significant difference in the proportion of responders across the four MGMT-quartiles, X2 (1, N = 319) = 8.48, P < 0.037 (Figure 2A).

Figure 2.

Gene expression changes related to TMZ-therapy response. Low MGMT-expression (Q1) is associated with higher response rate to TMZ therapy (TMZ-all) (A). Genes overexpressed in therapy-resistant tumors (B). Higher expression of FCGR2B in tumor specimens is associated with shorter OS and resistance to TMZ-nitrosourea combination therapies (C). Elevated MRPS17 expression is associated with resistance to multiple treatment strategies and is a predictor of worse prognosis in GBMs (D). *indicates significant differences.

Since younger age at diagnosis is a favorable prognostic factor, we investigated the number of responders across the two age groups as a function of MGMT-expression. Among TMZ-treated patients (TMZ-all) significantly more responders were identified among young patients (<55) both below (X2 (1, N = 162) = 5.71, P < 0.017) and above (X2 (1, N = 157) = 7.96, P < 0.005) median MGMT-expression, indicating, that in our dataset age is a more important predictor of therapy responsiveness compared with MGMT-expression.

Genes upregulated in chemotherapy-resistant tumors

Initial resistance to TMZ therapy poses a significant problem; therefore, we focused on gene expression patterns in tumor specimens of subsequent nonresponders. Identification of genes was conducted separately for each treatment combination and patient cohort listed in Table 2. The identified genes were subjected to survival analysis to assess the relation to OS, and only genes performing at FDR ≤10% in the Cox regression were considered significant. Significantly upregulated genes in nonresponders associated with worse clinical outcome corresponding to each treatment subcategory are summarized in Supplementary Table 2, available at Carcinogenesis Online.

The initial overexpression of the CSTA gene in patient samples was indicative of subsequent resistance to all treatment combinations containing TMZ, with the strongest predictive value in the cohort treated with TMZ + Nitroso-only (AUC = 0.714, P = 0.008). CSTA was also the single upregulated gene associated with subsequent resistance to TMZ therapy combined with other agents (TMZ combined). In a subset of patients treated with single-agent TMZ, we identified four genes (MRPS17, CSTA, LOX and PDPN) significantly overexpressed in initial samples of subsequent nonresponders also associated with worse outcome, although with AUC <0.7. A similar pattern of gene upregulation was identified among nonresponsive patients with highly expressed MGMT (Q2–Q4), where MRPS17, LOX, CSTA, IGFBP3 and ADM were significantly overexpressed (Supplementary Table 2, available at Carcinogenesis Online). Upregulation of the MRPS17 gene was associated with a lack of response to multiple treatment combinations, including treatment with single-agent TMZ (TMZ-only) and single-agent nitrosoureas (Nitroso-only) (Figure 2B and D), with a stronger effect in the latter (AUC = 0.72, P = 0.0011).

The highest proportion of responders (68%) were among patients treated with TMZ combined with nitrosoureas (TMZ + Nitroso-combined), where FCGR2B upregulation provided the strongest predictive value (AUC = 0.717, P < 0.00038) (Figure 2C). Ten genes were upregulated in patients resistant to TMZ + Nitroso-only, where DYNLT3 upregulation showed the strongest association with treatment resistance (AUC = 0.766, P = 0.00023). ALOX5, ALOX5AP, CSTA and SNX10 genes were overexpressed both among nonresponders treated with TMZ + Nitroso-only and in the treatment group TMZ + Nitroso-combined.

In patients treated with single-agent nitrosoureas, initial overexpression of 10 genes (ELOVL2, MEOX2, MRPS17, PIPOX, RPH3A, TMEM158, TMEM22, ARAP3, CNGA3 and FGFR3) was associated with worse survival, with AUC values above 0.7 for all, but one gene. MEOX2 and ELOVL2 provided a particularly strong predictive value of treatment resistance (AUC = 0.801, P < 0.00001 and AUC = 0.761, P = 0.000054, respectively).

Altogether 19 significantly upregulated genes were identified in samples of patients younger than 55 years of age and nonresponsive to TMZ, out of which CSTA, LOX, MRPS17 and PDPN were already depicted as overexpressed in nonresponders treated with TMZ-only. In the subset of patients with the lowest MGMT-expression (Q1) we discovered three upregulated genes (MDM2, MX1 and PLSCR1) associated with worse survival in therapy-resistant tumors, PLSCR1 being the strongest candidate (AUC = 0.712, P = 0.00036). These three genes were also upregulated in the younger cohort (<55) of TMZ-treated patients (Figure 2B). The list of significantly overexpressed genes associated with subsequent treatment resistance is summarized in Supplementary Table 2, available at Carcinogenesis Online.

In the ≥55 population, none of the identified upregulated genes fulfilled the statistical filtering criteria because gene expression values of the significant genes were too low. Thus, in the older population, we were unable to identify any age-specific potential biomarkers of TMZ resistance.

We subjected the significantly overexpressed genes to gene enrichment analysis performed using DAVID. Leukotriene biosynthesis and metabolism were significantly enriched in resistant tumor samples derived from patients treated with TMZ-containing therapies, represented by the expression of ALOX5, ALOX5AP and PLA2G5 genes (P = 0.0007).

Genes upregulated in responders to chemotherapy

Our second aim was to uncover predictive biomarkers of therapy effectiveness. We compared the expression of 10 103 genes between specimens of responders and nonresponders and focused on genes significantly (P < 0.05 and FC >1.44) overexpressed in responders with the mean expression difference ≥350, and FDR <10%. We selected genes significantly associated with OS at FDR ≤10%.

Among patients treated with TMZ-only, higher expression of three genes (DLL3, FERMT1 and PCSK1N) was linked to better survival outcome (Figure 3A), although none of the genes reached an AUC of 0.7. Upregulated DLL3 expression in tumor samples was associated with sensitivity to both single-agent TMZ therapy (TMZ-only) and TMZ combined with nitrosoureas (TMZ + Nitroso-only) (Figure 3B and C), with a higher predictive value in the latter (AUC = 0.697, P = 0.0092). Overexpression of DLL3 and SOX4 was associated with improved response to subsequent TMZ + Nitroso-combined. DLL3 overexpression was also associated with improved outcomes in patients with high MGMT-expression (Q2–Q4). In patients treated with single-agent nitrosoureas (Nitroso-only) or Nitroso without TMZ, higher expression of the CD24 gene was associated with better survival.

Figure 3.

Genes overexpressed in chemotherapy sensitive GBM tumor specimens. Upregulation of DLL3 is associated with responsiveness to all TMZ-based treatment strategies (A). High DLL3 expression is associated with improved survival outcome (B). Consistently upregulated DLL3 expression in tumor samples is associated with sensitivity to both TMZ-nitrosourea combinations and single-agent TMZ therapy (C).

A separate analysis of the two age cohorts (<55 and ≥55) in the TMZ-all treatment group revealed two upregulated genes (DLL3 and PCSK1N) linked to better survival outcome in younger patients who responded to therapy (Supplementary Table 3, available at Carcinogenesis Online). Curiously, no upregulated genes could be identified in specimens of responders in the older (>55) population receiving TMZ-containing regimes.

Independent validation of the results

Finally, we have extended our online accessible biomarker validation tool with the dataset utilized in the present study. The registration-free online interface enables independent validation of results presented in the current manuscript as well as the investigation of new future biomarker candidates. The analysis homepage can be accessed at www.rocplot.org/gbm/.

Discussion

Favorable prognostic factors of GBMs are limited to maximum safe resection, good performance status, completion of radiation and chemotherapies, and young age at diagnosis (4,31), and MGMT silencing predicts greater benefits of TMZ treatment (12). To improve prediction models, we assembled a dataset of GBM specimens and analyzed gene expression patterns before therapy initiation to elucidate potential mechanisms of subsequent treatment response. In the absence of data on MGMT-methylation, patients were classified based on MGMT-expression status (29,30). Age outweighed the effects of MGMT-expression, with enhanced responsiveness to TMZ in younger patients. Moreover, our results suggest a striking age-dependent heterogeneity: no consistently upregulated genes could be identified neither in nonresponders nor in responders within the older age group, contrary to expression patterns observed in younger patients. Older age has already been associated with more aggressive clinical behavior in GBMs, and tumorigenic pathway activations vary with the age of the patient (32). The heterogeneity of gene expression in older age groups suggests the possibility of age-dependent GBM subtypes, demanding further observations.

Most responders (68%) were in the treatment group receiving TMZ combined with nitrosoureas (with or without other agents), where FCGR2B (Fc Fragment of IgG Receptor IIb) upregulation provided a strong negative predictive value. Increased FCGR2B expression as part of a local immune signature has been previously associated with high-risk GBMs (33). Our results suggest that this immune phenotype and the enhanced local immune response forces the combined TMZ and nitrosourea treatment ineffective. Thus, patients with a high initial expression of FCGR2B may be good candidates for alternative treatment strategies and participation in clinical trials.

Lack of response to multiple treatment strategies was associated with elevated expression of MRPS17, proposing the existence of common resistance mechanisms against different treatment regimens. Strikingly, MRPS17 (mitochondrial ribosomal protein 17) is located on the short arm of chromosome 7 (7p11.2), the most frequently amplified chromosomal region in GBM, also including the EGFR gene at position 7p11.2 (34). It is well established that EGFR amplifications upregulate the RAS/RAF/MAPK and PI3K signaling pathway translating into increased proliferation and tumor cell survival. The amplification of MRPS17 has also been described in GBM (35). The involvement of MRPS17 in oxidative phosphorylation (OxPhos) may have consequences for tumor cell survival (36). Elevated MRPS17 expression along with upregulation of other MRPs is characteristic of epithelial breast cancer cells suggesting that cancer cells amplify OxPhos and fuel cancer metabolism (37). Studies on intrinsic and acquired resistance indicate that multiple rather than single molecular events and altered cell signaling are accountable for the lack of therapy effectiveness (38). The amplification/upregulation of MRPS17 at the same chromosomal position as EGFR likely affects cancer cell metabolic activity along functions related to proliferation and tumor cell survival, offering a potential to target multiple vulnerabilities with drug combinations. Our results are strikingly similar to findings formulated in another study integrating genomic aberrations and gene expression, and well-known drivers in GBM, such as EGFR may be synergistically acting on metabolic processes via coaltered, nearby genes (39).

In our dataset CSTA (Cystatin A) overexpression was associated with subsequent resistance to any therapy incorporating TMZ, including the combination of TMZ with nitrosoureas. CSTA functions as a regulator of proteolytic enzymes, cysteine cathepsins, playing a significant role in the invasion and growth of brain tumors (40). CSTA upregulation has previously been described in human malignant gliomas: CSTA positive cells in GBM tumor samples were located close to tumor blood vessels, particularly in leukocytes and inflammatory host cells, possibly reflecting the level of inflammatory cells in the tumor tissue (41). CSTA expression displayed a significant correlation with markers of invasive/migratory GBMs, CD68 and CXCR4, supporting its role in cysteine cathepsin-mediated proteolysis of the extracellular matrix, promoting tumor cell migration and invasion (41). High expression of CSTA has been previously identified as a significant prognostic factor of shorter survival in gliomas (41). Elevated initial CSTA expression in tumor samples may therefore signal an inflammatory tumor environment and an already ongoing process of cell invasion into the surrounding brain tissue, making treatments likely ineffective. In tumor samples with subsequent resistance to TMZ, gene ontology revealed enrichment of pathways related to leukotriene biosynthesis and metabolism. Leukotrienes are proinflammatory lipid mediators, functioning as critical immune modulators of leukocyte migration, and are implicated in a number of inflammatory disorders, including cancer (42). Activated leukotriene metabolism in TMZ-resistant tumors also supports the likelihood of the role of inflammation in treatment resistance.

GBMs encompass biologically distinct tumor groups with markedly different locations, age of onset and clinical outcome (6,43,44). The majority (~90%) are classified as primary GBM that occur without evidence of a less malignant precursor lesion in older patients, whereas secondary GBMs progress from low-grade gliomas and occur in younger patients, with a significantly better prognosis (45). In our dataset, overexpression of three genes delineated a subset of younger patients with low MGMT-expression resistant to TMZ within a population generally associated with a favorable outcome. These genes have already been implicated in glioma etiology: MX1 is suggested to be part of the interferon-related gene signature for DNA damage resistance that predicts poor survival, particularly within the proneural subtype consisting of secondary GBMs (46). MDM2 negatively regulates and destabilizes the TP53 protein (47), and amplification and overexpression of the oncogene is observed in 8–10% of GBMs (48). MDM2 also provides a potentially druggable target: the inhibitor AMG232 exhibited in vitro activity against GBM cell lines and stem cells with higher selectivity against p53 wild type over p53 mutant cells (49). Upregulated PLSCR1 among low-grade glioma patients is related to a worse prognosis and a higher risk of cancer recurrence (50). Our results, thus, reconfirm the importance of the identified genes as promising potential treatment targets.

Upregulation of DLL3 was repeatedly identified in tumor specimens of subsequent responders to TMZ single-agent or combination therapies, and elevated DLL3 expression was associated with improved clinical outcome. DLL3 (Delta-like canonical Notch ligand 3) is a member of the DSL family of Notch ligands and inhibits Notch-pathway activation (51). Based on gene expression GBMs have been grouped into proneural, neural, classical and mesenchymal transcriptomic subtypes that correlates with responsiveness to therapy, with greatest benefits in classical (52), and the best prognosis among the most differentiated proneural subtypes (53). DLL3 is a signature marker of proneural GBMs, characterizing mostly younger patients associated with improved prognosis (53), suggesting that DLL3 activity on Notch-signaling may limit tumor growth. Comparison across more than 20 tumor types revealed the highest DLL3 expression in low-grade gliomas (54), with homogeneous expression in IDH-mutant low-grade gliomas in contrast to patchy or nonexistent expression in IDH-wild type GBMs (55). Our entire dataset consisted of GBMs with no low-grade gliomas present among samples; moreover, information about IDH mutation status was not available. According to Phillips et al. almost a third of GBM samples could be classified as proneural (53); thus, DLL3 overexpression may be a predictor of good prognosis not only in low-grade gliomas but in GBMs as well. Our data further expand the role of DLL3 overexpression as a predictor of TMZ effectiveness, with implications for potential combination therapies: glioma cell lines with high DLL3 expression could successfully be targeted in vitro by the anti-DLL3 antibody–drug conjugate rovalpituzumab tesirine (Rova-T) (55).

Analyses of longitudinal genetic changes in matched de novo versus recurrent GBM patient tumor, together with studies on recurrence-initiating stem-like cancer cells have provided a wealth of information regarding molecular pathways and individual genes involved in tumor recurrence and TMZ resistance (56,57). The discrepancy between described previously genes involved in resistance and recurrence and our results stems from two inherently different approaches: we analyzed tumor specimens at the time of surgery, before any therapy initiation, moreover we focused at initial gene expression instead of the nature of genomic alterations.

There are a number of limitations of our study. First of all, the majority (91.4%) of patients received radiation therapy along with chemotherapy; therefore, the observed effects might be associated with radiation instead of treatment with TMZ or nitrosoureas. Only 15 patients were treated with TMZ who did not receive radiation, not permitting further analyses between responders and nonresponders. Nevertheless, we were still able to depict group differences related to various treatment regimens.

Another major limitation is the lack of validation of the findings with experimental results such as qPCR or IHC from clinical specimens with known treatment history. Nevertheless, our final motivation was to develop a broadly accessible transcriptome-level online validation tool for the initial testing of predictive biomarker candidates in glioblastomas. Our system provides highly reliable data for the selection of the most robust biomarker candidates for subsequent validation studies (58). The findings of the present analysis illustrate our system’s applicability in future validation studies.

In summary, our unveiled gene expression patterns reconfirmed previously identified biomarkers, revealed novel markers of responsiveness and uncovered potential ‘resistant phenotypes’ contributing to the inefficiency of multiple treatment strategies. Our results narrow the list of genes that deserve attention in future validation studies in GBM treatment response, help to identify candidates for more intense interventions while also provide novel drug targets for personalized treatment strategies.

Click here for additional data file.