A gene signature predicting prognosis of patients with low- grade gliomas receiving temozolomide therapy

doi:10.21203/rs.3.rs-2822634/v1

Download PDF

Research Article

A gene signature predicting prognosis of patients with low- grade gliomas receiving temozolomide therapy

https://doi.org/10.21203/rs.3.rs-2822634/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

Temozolomide (TMZ) has been used as a first-line therapy against low-grade gliomas (LGGs) combined with other chemotherapy drugs. However, there has been no reliable index predicting TMZ response of patients with LGGs. In this study, we aim to investigate the relationship between gene expressions and the prognosis of TMZ therapy in LGGs. We integrated transcriptome and clinical data of 171 LGGs from the Chinese Glioma Genome Atlas (CGGA). Consensus LASSO Cox regression was used to identify 14 key genes related to different clinical outcomes under TMZ chemotherapy. We constructed and evaluated a risk score based on the 14 genes. Patients with LGGs of lower risk scores (low-risk group) generally had better survival than those LGGs of higher risk scores (high-risk group), which is independent of clinicopathological factors. High-risk patients showed activation of innate and humoral-type immunity. The prognostic contribution of the risk score was validated in an independent validation cohort of 65 patients. Besides, combined with three independent predictors (grade, IDH1 mutation status, and chr1p19q co-deletion status), we further developed a nomogram to predict the benefit from TMZ treatment in LGGs. Our results indicate that transcriptome-based index can optimize the treatment strategy for patients with LGGs under TMZ therapy.

low-grade gliomas (LGGs)

Temozolomide (TMZ) therapy

genes

risk scores

prognosis

Gliomas, the most common primary tumors of the central nervous system, are grouped by the World Health Organization (WHO) with histopathology from grade I (lowest grade) through grade IV (highest grade)[1,2]. The low-grade gliomas (LGGs), a diverse group of grade I and II, generally have a lower degree of malignancy, a higher treatment effectiveness and a better prognosis in comparation with high-grade gliomas (grade III and IV). However, LGGs, incurable disease, have high risk of recurrence because of the characteristics of diffuse growth and the potential for progressing to a higher tumor grade[3,4]. The median overall survival for LGGs patients ranges from 2 to 15 years based on many factors, such as tumor histopathological, molecular phenotype, growth rate and therapeutic strategy[5].

Nowadays, Temozolomide (TMZ) is used as the first-line therapy against LGGs combined with other chemotherapy drugs. It shows good tolerability with limited side effects[6], however, treatment effectiveness differs among patients receiving the same therapeutic strategy. Some data already suggested a close connection between gene mutations, malignant transformation and TMZ response[7,8]. Some previous studies managed to portray the clinical behavior of LGGs related with genomic mutation. Other studies constructed a prognostic model using biological features (DNA repair ability as well as cell cycle processes) that are crucial to tumor recurrence[9,10]. Yet the relationship between genes and the prognosis of patients receiving Temozolomide (TMZ) adjuvant chemotherapy has been left undiscovered. It remains a challenge to effectively assess the prognosis of low-grade gliomas .

The Chinese Glioma Genome Atlas (CGGA) project includes the long-term follow-up data of LGGs patients with TMZ treatment, providing the possibility for establishment of treatment prognosis scores. In this study, we gathered a large number of samples to identify the connection between genes and the prognosis of patients receiving TMZ, therefore confirming that prognosis scores, consist of some key genes, are indeed related to the prognosis of TMZ therapy. Our results provide further evidence of assessing the Prognosis of TMZ treatment through key genes for the better benefit of patients.

1. Data sources

We curated the clinical and gene expression data of LGG patients from the Chinese Glioma Genome Atlas. The first batch of data from the CGGA database was used as a training set to develop TMZ treatment prognosis scores, and the second batch of data from the CGGA database was used to evaluate the prognostic performance of risk scores. Only patients with complete overall survival time and transcriptomics data were retained. The gene expression data were log2- transformed for subsequent analysis.

2. Prognostic gene analysis

The Cox regression was used to screen genes associated with the prognosis of TMZ treatment after correct significant clinical variables. Genes with p values less than 0.01 were considered to be potentially associated with treatment prognosis. Then these genes were then classified as protective or risky according to whether HR was less than 1. Gene set enrichment pathway analysis was conducted based on Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG), and Reactome database by using online tools Metascape. The protein-protein interaction network analysis was carried out by the STRING database, and only the high confidence edges were retained. The Molecular Complex Detection (MCODE) algorithm has been applied to identify densely connected network components. Pathway enrichment analysis was also applied to each MCODE component independently, and the three best-scoring terms by p-value have been retained as the functional description of the corresponding components.

3. Establishment of treatment prognosis scores

These genes which were potentially associated with treatment outcomes were subjected to penalized multivariate Cox proportional hazards survival modeling using an algorithm for variable selection based on L1-penalized (Lasso) estimation. A penalty parameter, λ1, reflecting the predictive potential and calculated by cross-validation, was inflicted upon the gene expression signals during survival modeling. Using this model, some key genes were selected via 10-fold cross-validation in the training sets. Due to the randomness of cross-validation, the selected key genes were variable, so we used consensus Lasso regression to identify key genes. That is, genes retained at high frequency after multiple run Lasso regression can be considered to be the most influential for treatment prognosis.

The order of frequency represented the degree of influence of these genes, and then these genes were incorporated into the Cox model in turn. We then used the area below the time-dependent ROC to evaluate the prognostic efficacy of the model. When the ROC reaches a certain point, the model is optimal and contains the least genes, which will be identified as the key genes for prognosis. Finally, we used 10-fold cross-validation in the training set to determine the weight coefficients of these genes, and the cumulative sum of the product of gene weight and corresponding expression was used as the TMZ prognosis score. And the clinical variables were pooled and a visualized nomogram was drawn using R package rms.

4. Evaluation of treatment prognosis scores

We systematically assessed the prognostic efficacy of risk scores in training sets and independent validation sets. Patients were firstly divided into high- and low-risk groups using themedian risk score as the cutoff. Then the Kaplan Meier method was used to draw the survival curves of the two groups, and the log-rank test was used to judge the difference. The prediction accuracy was evaluated by time-dependent AUC obtained by R package timeROC. We used calibration curves to assess the agreement between nomogram predicted risks and observed risks, and performed decision curve analyses to measure the net benefit of identifying truly high-risk patients who should intervene because of their use of nomogram as a screening tool.

5. Statistical analysis

R (V 4.1.3) software was used for the above statistics and analysis. The images were stitched together by Adobe Illustrator software. Wilcoxon rank-sum test was used for comparison of box plots. Spearman coefficient was used for correlation analysis. The Chi-square test (Fisher's exact test if necessary) was used to compare the clinical features of the two groups. Multivariate logistic regression analysis was used to evaluate the clinical features that affected clustering. Kaplan Meier method was used to plot the survival curve. Cox analysis was used to assess characteristics associated with survival. All hypothesis tests were two-sided tests. P-values of multiple tests were corrected by the FDR method.

1. Analysis of genes related to the prognosis of LGG patients after TMZ treatment

We obtained 171 low-grade gliomas (LGG) patients receiving temozolomide (TMZ) chemotherapy from the CGGA database as a training set to screen for genes associated with TMZ treatment prognosis (Figure 1A). Univariate Cox regression showed that tumor stage, radiation therapy status, IDH1 mutation status, and chr1p19q co-deletion status had an impact on overall survival (Figure 1B, Supplementary Figure 1). Therefore, we used Cox regression to correct these clinical variables to obtain genes associated with the prognosis of TMZ treatment. A total of 488 risk genes and 128 protective genes passed the screening threshold (p < 0.01; Figure 1C, Supplementary Table 1). Functional enrichment of these genes revealed that protective genes were associated with nervous system development, while risk genes were associated with multiple pathways associated with tumor progression including cell cycle, immune response, EMT, and extracellular matrix organization (Figure 1D, Supplementary Table 2).

We further analyzed the internal correlations of risk genes, protective genes, and non- significant genes, and found that the overall correlations of protective genes and risk genes were higher than those of non-significant genes, suggesting that these genes formed specific functional modules (Figure 1E). At the same time, the overall correlation of protective genes was significantly higher than that of risk genes, suggesting that risk genes were involved in multiple functions. In addition, we analyzed the protein network formed by risk and protective genes respectively. The degree distribution of the two networks showed that the network of risk genes reflected multi-modularity (Figure 1F). Moreover, the centrality measures of the two networks are different in distribution. The mean betweenness centrality of the network of risk genes was low (Figure 1G), and the mean closeness centrality was high (Figure 1H). Network module analysis revealed that risk genes formed multiple functional modules, including lipid metabolism, cell cycle, interferon response, and Wnt signaling pathway (Supplementary Figure 2). These results suggest that TMZ prognostic risk genes are involved in various carcinogenic processes.

2. Development of risk score associated with the benefit of TMZ treatment in LGG

We used consensus LASSO Cox regression to screen key genes for prognostic genes and repeated 1000 iterations to obtain the occurrence frequency of each gene (Figure 2A), of which 14 genes appeared more than 950 times. We then included high-frequency genes in the model in turn and used a 10-fold cross-validation method to evaluate the prognostic model performance by the area under time-dependent ROC at varying follow-up times. The results showed that the prognostic model had the highest predictive power when predicting 3-year survival, while its predictive power was relatively poor for 1-year survival (Figure 2B). Based on the overall prediction of different survival conditions, we selected a model consisting of 14 genes. For these 14 genes, we used Cox regression to determine the contribution of these genes based on a 10-fold cross-validation method. These features, ranked by median weight, included PAX3 (1.06), ADM2 (0.76), GRB14 (0.44), CYP4F12 (0.37), HELZ2 (0.24), RGS16 (0.19), IL13RA2 (0.15), CHRNA1 (0.06), IGF2BP3 (-0.08), HOXA4 (-0.15), IL34 (-0.28), KCNIP3 (-0.32), CACNA1H (-0.33), and GNAL (-0.52) (Figure 2C). Finally, the cumulative sum of the product of gene weight and corresponding expression was used as the TMZ prognosis score.

We assessed the prognostic efficacy of risk scores in training sets. The distribution of risk scores and survival status indicated that patients with lower risk scores generally had better survival than did those with higher risk scores (Figure 2D). When patients were divided into the high-risk groups based on the median risk score, significant differences were shown in survival time between the two groups, which median survival was 1680 days (95% CI 855-2696) for the high-risk group and not reach for the low-risk group (HR 8.58, 95% CI 4.32- 17.05; p<0.0001, log-rank test; Figure 2E). And we found the excellent prognostic accuracy of the risk scores with time-dependent ROC analysis at varying follow-up times, especially for 3- and 5-year survival predictions (Figure 2F).

To confirm that the risk scores had similar prognostic values in different populations, we applied it to the independent validation set of 65 patients from CGGA batch 2 data. There were also differences in the distribution of risk scores and survival status in this data set (Figure 2G). We observed a significant difference in survival time between the high and low-risk groups, which median survival was 1048 days (95% CI 827-1654) for the high-risk group and not reach for the low-risk group (HR 4.94, 95% CI 2.41-10.12; p<0.0001, log-rank test; Figure 2H). Risk scores also showed excellent predictive accuracy in independently validated cohorts, with the AUC of survival years at a different time above 0.8 (Figure 2I).

3. The prognostic power of risk scores was independent of clinicopathological factors

After multivariable adjustment by clinicopathological variables including grade, radio status, IDH mutation status, and chr1p19q co-deletion status, the risk score remained a powerful and independent factor in the training cohort (HR 2.73, 95% CI 2.09-3.58, p < 0.001; Figure 3A). We also noted similar results in the independent validation cohort (HR 1.28, 1.05-1.55, p = 0.013; Figure 3B). The IDH1 mutation status and chr1p19q co-deletion status are key genomic markers related to the prognosis of glioma treatment. When stratified by IDH1 mutation status, the risk score was still a clinically and statistically significant prognostic model in the train and independent validation cohort (Figure 3CD). We got a similar conclusion, in which patients with chr1p19q wild-type and the high-risk score had the worst survival when patients were stratified according to the co-deletion status of chr1p19q (Figure 3EF). We also focused on the prognostic effect of tumor stage on risk score, and the risk score can distinguish patients with different stages from two groups with significant differences in survival (Figure 3GH). Further, we explored the prognostic performance of risk scores in glioblastoma (GBM). Although there are certain differences between LGG and GBM in the genome and other aspects, we found that the risk score still showed certain prognostic ability, and patients with lower risk scores tended to have better survival benefits (p = 0.062 in train data sets and p = 0.015 in independent validation sets; Figure 3IJ).

4. The risk score was related to the immune cell infiltration of LGG

We then analyzed the relationship between risk score and immune microenvironment. We use a variety of algorithms including ESTIMATE, XCELL, and CIBERSORT to measure the composition of the tumor's immune microenvironment. We first analyzed the association of key genes in the model with tumor purity, immune score, stromal score, and different immune cell infiltration proportion. The HOXA4, CYP4F12, and ADM2 were positively correlated with the immune and stromal score but negatively correlated with tumor purity (Figure 4A). The immune cells that are positively associated with these genes include macrophages, naïve CD8+ T cells, and activated dendritic cells (Supplementary Figure 3, Figure 4B).

The GNAL, on the other hand, was negatively correlated with the immune score, stromal score, and some myeloid lymphocyte infiltration (Figure 4AB). We then assessed differences in immune microenvironments between high and low-risk groups. Patients with high-risk scores had higher immune and stromal scores and lower tumor purity (Figure 4C), indicating that unfavorable prognosis in the high-risk group may be associated with the variation in tumor immune microenvironment. Further immune cell infiltration analysis showed that the activated dendritic cells, M1 type macrophages, B cells, and Th2 subset of CD4+ T cells were higher in patients with high-risk scores, while the B cell plasma and naïve CD8+ T cell were lower in patients with high-risk score (Supplementary Figure 4, Figure 4D). These results indicated that high-risk patients showed activation of innate and humoral-type immunity.

5. Development of a nomogram for predicting the benefit of TMZ treatment in LGG

To better predict the prognosis of LGG patients after receiving TMZ therapy in the clinic, a prognostic nomogram was developed by integrating risk score and three independent predictors of mortality including grade, IDH1 mutation status, and chr1p19q co-deletion status from the above analyses into a multivariate Cox regression model (Figure 5A). The calibration plot showed that the nomogram performed well in predicting patient survival according to an ideal model in training and independent validation sets (Figure 5BC). To evaluate the clinical benefits of the nomogram model, we performed decision curve analysis in both the training and independent validation sets. With 3-year survival as the endpoint, the curve showed that nomogram presented more clinical net benefits than several competing intervention strategies, namely, intervention for all, intervention for none, and intervention based on different clinical indicators (Figure 5DE).

Temozolomide chemotherapy is the main treatment for LGGs, and its therapeutic effect varies greatly among different patients. Hence, we systematically assess the influence of transcriptomeson the prognosis after Temozolomide chemotherapy, attempting to find the relationship between some key genes and prognosis.

Tumor stage, radiation therapy status, IDH1 mutation status, and chr1p19q co-deletion status, these elements all have been identified as important factors related to patients’ overall survival. Therefore, we managed to find other associated factors that independent of clinicopathological features. Genes associated with TMZ treatment prognosis is emerging into our sights. In our study, we classified the genes into risk genes, protective genes and non-significant genes. The protective genes were associated with nervous system development, while risk genes were associated with tumor progression into higher grades. Finally, we found out 14 key genes from 488 risk genes and established the TMZ-prognosis score, which was the cumulative sum of the product of gene weight and corresponding expression.

We identified risk score that was associated with the prognosis of Temozolomide therapy in LGGs. The numerical of risk scores were related to the prognosis of patients--patients with lower risk scores generally had better survival and significant differences were shown in survival time between the low and high gropes. Two large randomized EORTC studies already attempted to predict the prognosis of LGGs though risk scores, which include 5 factors. A patients with low-risk prognostic score (<two factors) generally had a longer-survival than with a high-risk score(3-5 factors )[11].

However, no research has figured out the relationship between key genes and TMZ related prognosis. Our result fills the gap. We found that for the high-risk group, median survival was 1680 days; while not reach for the low-risk group. What’s more, we suggested that this module had better performance in 3- and 5- year survival predictions than 1- or 7-years survival predictions, although the AUCs of survival years at a different time are all above 0.8. However, for some LGGs patients with low risk scores, we were unable to accurately evaluate their survival years. Whether can we incorporate more genes to solve this problem, is a challenge need us to face.

As we know, earlier clinical stage means better prognosis. Data suggests that LGGs patients without 1p and 19q co-deletion is a favorable news, because they may have a better response to chemotherapy and longer overall survival[12]. Radiation therapy status is also an important factor, which has great influence on overall and progression-free survival times[13]. Although it has no predictive value for chemotherapy response, IDH1 mutation is an independent prognostic factor in LGGs[5,14,16]. The prognostic power of risk scores was independent of the above related clinicopathological factors. However, data from other messages suggests that besides those elements, prognosis also depends on other elements like age at diagnosis, Karnofsky performance status (KPS), first presenting symptom, extent of resection and certain molecular markers[17]. There is one limitation of this study that should be noted here. In our results, because of limited data from the CGGA database, we only make more influencial factors as independent factors, while some little influencial factors are not taken into account. Our risk scores may perform better in mostly LGGs patients, but it may has accurate prediction for some patients with early age or basic diseases.

Interestingly, we found the association between risk score and immune microenvironment, suggesting that the indirect effect of the immune microenvironment on prognosis. High-risk patients are more likely to show activation of innate and humoral-type immunity,which may cause worse prognosis. Besides, to better predict the prognosis of LGG patients after receiving TMZ therapy in the clinic, a prognostic nomogram was developed, which integrate risk score and three independent predictors of mortality including grade, IDH1 mutation status, and chr1p19q co-deletion status into a multivariate Cox regression model. It provides comprehensive and individualized analysis of different patients to assess the prognosis. We would like to provide free online software for patients and health care providers.

In summary, we not only found an association between key genes and prognosis but also established a valuable risk scores tool for LGGs after receiving Temozolomide chemotherapy. This tool provides an individualized estimate of survival based on specific genes of different patients. By predicting the prognosis of patients, clinicians can better design follow-up mode and optimize treatment to improve survival outcomes. To facilitate clinical use of this nomogram, free online software for its implementation will be provided.

Ethical Approval

This study did not involve ethical approval.

Competing interests

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.

Authors' contributions

MXL had contributions to the conception and design of the study; WYZ performed the data analysis, interpretation, and manuscript drafting; LGQ and DJY contributed to the data acquisition and assessment. All authors contributed to manuscript reviewing, approved the final version of the manuscript, and agreed to be accountable for all aspects of the study.

Funding

The authors declare no conflict of interest.

Availability of data and materials

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

D. N. Louis, A. Perry, G. Reifenberger. The 2016 World Health Organization Classification of Tumors of the Central Nervous System: a summary. Acta Neuropathol 2016 Vol. 131 Issue 6 Pages 803-20
D. A. Forst, B. V. Nahed, J. S. Loeffler. Low-Grade Gliomas. Oncologist 2014 Vol. 19 Issue 4 Pages 403-13
A. F. Haddad, J. S. Young, J. Y. Oh. The immunology of low-grade gliomas. Neurosurg Focus 2022 Vol. 52 Issue 2 Pages E2
I. K. Mellinghoff, B. M. Ellingson, M. Touat. Ivosidenib in Isocitrate Dehydrogenase 1-Mutated Advanced Glioma. J Clin Oncol 2020 Vol. 38 Issue 29 Pages 3398-3406.
A. F. Hottinger, M. E. Hegi and B. G. Baumert. Current management of low-grade gliomas. Curr Opin Neurol 2016 Vol. 29 Issue 6 Pages 782-788.
S. G. R. McDuff, J. Dietrich, K. M. Atkins. Radiation and chemotherapy for high-risk lower grade gliomas: Choosing between temozolomide and PCV. Cancer Med 2020 Vol. 9 Issue 1 Pages 3-11.
B. E. Johnson, T. Mazor, C. Hong. Mutational analysis reveals the origin and therapy-driven evolution of recurrent glioma. Science 2014 Vol. 343 Issue 6167 Pages 189-193.
Y. Yu, J. Villanueva-Meyer, M. R. Grimmer. Temozolomide-induced hypermutation is associated with distant recurrence and reduced survival after high-grade transformation of low-grade IDH-mutant gliomas. Neuro Oncol 2021 Vol. 23 Issue 11 Pages 1872-1884.
G. Li, F. Wu, F. Zeng. A novel DNA repair-related nomogram predicts survival in low-grade gliomas. CNS Neurosci Ther 2021 Vol. 27 Issue 2 Pages 186-195.
E. H. Bell, P. Zhang, E. G. Shaw. Comprehensive Genomic Analysis in NRG Oncology/RTOG 9802: A Phase III Trial of Radiation Versus Radiation Plus Procarbazine, Lomustine (CCNU), and Vincristine in High-Risk Low-Grade Glioma. J Clin Oncol 2020 Vol. 38 Issue 29 Pages 3407-3417
M. Deacu, S. Popescu, A. Docu Axelerad. Prognostic Factors of Low-Grade Gliomas in Adults. Curr Oncol 2022 Vol. 29 Issue 10 Pages 7327-7342.
J. S. Smith, A. Perry, T. J. Borell. Alterations of chromosome arms 1p and 19q as predictors of survival in oligodendrogliomas, astrocytomas, and mixed oligoastrocytomas. J Clin Oncol 2000 Vol. 18 Issue 3 Pages 636-45.
C. Leighton, B. Fisher, G. Bauman. Supratentorial low-grade glioma in adults: an analysis of prognostic factors and timing of radiation. J Clin Oncol 1997 Vol. 15 Issue 4 Pages 1294-301.
M. Y. Li, Y. Y. Wang, J. Q. Cai. Isocitrate dehydrogenase 1 Gene Mutation Is Associated with Prognosis in Clinical Low-Grade Gliomas. PLoS One 2015 Vol. 10 Issue 6 Pages e0130872.
I. Duran and J. J. Raizer. Low-grade gliomas: management issues. Expert Rev Anticancer Ther 2007 Vol. 7 Issue 12 Suppl Pages S15-21.
H. Gittleman, A. E. Sloan, J. S. Barnholtz-Sloan. An independently validated survival nomogram for lower-grade glioma. Neuro Oncol 2020 Vol. 22 Issue 5 Pages 665-674.

No competing interests reported.

Download PDF

Version 1

posted

You are reading this latest preprint version

A gene signature predicting prognosis of patients with low- grade gliomas receiving temozolomide therapy

Status:

Version 1

Abstract

Figures

Introduction

Methods

Results

Discussion

Declarations

References

Additional Declarations

Supplementary Files

Status:

Version 1