A radiomic signature as a non-invasive predictor of progression-free survival in patients with lower-grade gliomas.

OBJECTIVE: The aim of this study was to develop a radiomics signature for prediction of progression-free survival (PFS) in lower-grade gliomas and to investigate the genetic background behind the radiomics signature.
METHODS: In this retrospective study, training (n = 216) and validation (n = 84) cohorts were collected from the Chinese Glioma Genome Atlas and the Cancer Genome Atlas, respectively. For each patient, a total of 431 radiomics features were extracted from preoperative T2-weighted magnetic resonance images. A radiomics signature was generated in the training cohort, and its prognostic value was evaluated in both the training and validation cohorts. The genetic characteristics of the group with high-risk scores were identified by radiogenomic analysis, and a nomogram was established for prediction of PFS.
RESULTS: There was a significant association between the radiomics signature (including 9 screened radiomics features) and PFS, which was independent of other clinicopathologic factors in both the training (P < 0.001, multivariable Cox regression) and validation (P = 0.045, multivariable Cox regression) cohorts. Radiogenomic analysis revealed that the radiomics signature was associated with the immune response, programmed cell death, cell proliferation, and vasculature development. A nomogram established using the radiomics signature and clinicopathologic risk factors demonstrated high accuracy and good calibration for prediction of PFS in both the training (C-index, 0.684) and validation (C-index, 0.823) cohorts.
CONCLUSIONS: PFS can be predicted non-invasively in patients with LGGs by a group of radiomics features that could reflect the biological processes of these tumors.

Introduction

Gliomas are the most common and fatal primary tumors in the central nervous system (Nuno et al., 2013). Lower-grade gliomas (LGGs), referring to the World Health Organization (WHO) grade II and III gliomas, account for approximately 43.2% of gliomas (Cancer Genome Atlas Research et al., 2015; Jiang et al., 2016). The variable biological behaviors of LGGs result in a wide range of progression-free survival (PFS) times. Accurate prediction of PFS can provide crucial information regarding treatment of gliomas in clinical practice. More specific imaging examinations would be indicated for patients at high risk for tumor progression. Furthermore, identification of a poor PFS helps to determine whether more aggressive treatment should be administered (Zhang et al., 2017).

Magnetic resonance imaging can provide more comprehensive information about tumor heterogeneity than focal tissue samples, and the emerging field of radiomics holds great potential for facilitating better clinical decision-making (Gillies et al., 2016). Radiomics refers to the conversion of digital medical images into mineable high-dimensional data, and its potential application in clinical practice has attracted much attention in recent years. Radiomic analysis has been used to predict the diagnosis, prognosis, response to treatment, and underlying genomic patterns in several types of cancer (Huang et al., 2016; Kickingereder et al., 2016a; Kickingereder et al., 2016b; Li et al., 2016a; Yamamoto et al., 2014). However, whether or not radiomic features have value in prediction of PFS in patients with LGGs is still unclear.

The aim of this study was to develop a novel approach to preoperative, non-invasive, and individualized assessment of PFS in patients with LGGs. We further developed a practical nomogram that incorporated the radiomic signature and other clinicopathologic characteristics for prediction of PFS in patients with LGGs and investigated the biological processes underlying this radiomic signature.

Material and methods

Patients

This study was approved and reviewed by the institutional review board of our Hospital. Two hundred and sixteen patients were enrolled from the CGGA database () as a training set and a further 84 cases were enrolled from TCGA database () as a validation set. All patients in both cohorts met the following inclusion criteria: (a) pathologically confirmed grade II or III glioma according to the typical histological-based WHO classification (Louis et al., 2007); (b) no history of preoperative therapy; (c) availability of preoperative T2-weighted magnetic resonance (MR) images; and (d) availability of data on PFS, clinical characteristics, and genetics. PFS was defined as the time from the date of the initial diagnosis until tumor progression.

Acquisition of MRI data and tumor segmentation

MR images of patients from the CGGA and TCGA databases were obtained from the CGGA imaging database () and the Cancer Imaging Archive (), respectively. The extraction of radiomic features was performed on T2-weighted MR images, because the T2-weighted sequence is well accepted in identifying tumor borders of low-grade gliomas (Kinoshita et al., 2016; Ricard et al., 2007; Wang et al., 2015). We did not use the T1-weighted or contrast enhancement images since it was difficult to identify the tumor borders of LGGs on these sequences. Tumors were segmented on T2-weighted images by two neuroradiologists (XC, JM) using MRIcron software (http://www.mccauslandcenter.sc.edu/mricro). Both neuroradiologists had >15 years of experience in neuroradiology and were blinded to the clinical data. Abnormal hyperintense signals on the T2-weighted MR images were identified as tumor regions and signals from cerebrospinal fluid were excluded. A third senior neuroradiologist (SL) with >20 years of clinical experience in interpretation of brain MRI subsequently re-evaluated the segmented lesions and made the final decision in the event of disagreement.

Extraction of radiomic features

First, normalization (z-score transformation) of image intensity was performed on the whole brain image to transform arbitrary MRI signal intensity values into standardized intensity ranges, thereby avoiding heterogeneity bias. Next, quantitative radiomic features were extracted using the automated approach reported in a previous study which provided a detailed description of each feature in its supplementary material (Aerts et al., 2014). Four hundred and thirty-one radiomic features were extracted for each patient and divided into four groups: (a) group 1, first-order statistics (n = 14) that quantitatively described the distribution of the signal intensity of the images; (b) group 2, shape- and size-based features (n = 8) that quantified the shape and size of the tumor; (c) group 3, textural features (n = 33) that were calculated from the gray-level run-length and co-occurrence matrix in addition to reflecting intratumoral heterogeneity; and (d) group 4, wavelet features (n = 376) that were derived from the features in groups 1 and 3 by wavelet decomposition. The radiomic features were extracted using MATLAB 2014a software (MathWorks, Natick, MA, USA) and are presented in Supplementary Table 1.

Molecular analysis and whole-genome gene profiling

Isocitrate dehydrogenase (IDH) mutations were detected by pyrosequencing in the CGGA training cohort (Zhang et al., 2014) and the molecular profiles (Ceccarelli et al., 2016) of patients in the validation cohort were collected from the TCGA database. Microarray analysis was performed for 47 patients in the training cohort using the Agilent Whole Human Genome Array (Agilent Technologies Inc., Santa Barbara, CA, USA) in accordance with the manufacturer's protocol (Yan et al., 2012). The integrity of total RNA was examined using a 2100 Bioanalyzer (Agilent). Biotinylated cRNA and cDNA were synthesized and hybridized to the array. The data were obtained using the Agilent Feature Extraction Software (version 9.1) and the Agilent G2565BA Microarray Scanner System. Probe intensities were normalized using GeneSpring GX 11.0.

Construction and validation of the radiomic signature

Univariate Cox regression was performed to screen for prognostic radiomic features in the training cohort. P < 0.05 was deemed to be statistically significant. Based on hierarchical clustering, the association between these prognostic features and PFS was represented using a heat map (the radiomic features were normalized by z-score transformation only when the heat map was delineated). Subsequently, the least absolute shrinkage and selection operator (LASSO) Cox regression model was applied to the screened prognostic features to further select the most useful prognostic radiomic features. These features were then integrated into a radiomic signature, and an individualized risk score was calculated from a linear combination of the selected features weighted by their respective coefficients (β):

The patients in the training and validation cohorts were then classified into low-risk and high-risk groups according to a fixed cutoff value. The relationship between the radiomic risk score and PFS was evaluated in the training cohort and then tested in the validation cohort using Kaplan-Meier survival and Cox regression analyses. Finally, the prognostic significance of the radiomic signature was assessed by Kaplan-Meier survival analysis in subgroups of the training cohort. Patients were divided into subgroups according to IDH status (wild-type vs. mutant), age (younger than 40 years vs. 40 years or older), sex (male vs. female), tumor grade (WHO II vs. WHO III), seizure (non-seizure vs. seizure), and presence of oligodendroglioma component (oligocomponent vs. astrocytomas).

Construction of an individualized PFS prediction model

A nomogram was established as an individualized PFS prediction model in the training cohort. The final selection of the model for the nomogram was conducted using a backward step-down selection process based on the Akaike information criterion (Li et al., 2016b). The prognostic performance of the nomogram was estimated in the training cohort and then tested in the validation cohort. The concordance index (C-index), which is often used to assess the discriminative ability of prognostic models in survival analysis, was used as a quantitative measurement of the performance of the nomogram. Calibration curves were constructed to compare the probability values determined by the nomogram with the observed survival fractions.

Radiogenomic analysis

Transcriptome data for 47 patients in the CGGA cohort were used for the radiogenomic analysis. The Pearson correlation algorithm was used to screen genes for their association with the radiomic signature; the association was considered to be statistically significant when the absolute value of Pearson correlation coefficient was >0.4 and the P-value was <0.05. Gene ontology (GO) analysis was performed to investigate the underlying biological processes of the radiomic signature based on DAVID Bioinformatics Resources (http://david.ncifcrf.gov/). The top 200 positive/negative genes that were significantly associated with each feature in the radiomic signature (the absolute value of Pearson correlation coefficient > 0.4 and P < 0.05) were subjected to GO analysis to reveal the underlying biological processes involved in each feature. The radiogenomic analysis was also performed in the validation cohort (84 patients).

Statistical analysis

R version 3.3.2 software (R Foundation for Statistical Computing, Vienna, Austria) was used for all statistical analyses and to generate figures. Univariable and multivariable Cox regression models and Kaplan-Meier survival analysis were performed in the “survival” package. The LASSO Cox regression model and nomogram were constructed using the “glmnet” and “rms” packages, respectively. The processing code for the generation of radiomic risk score and nomogram were provided in the Supplementary material. Differences in clinicopathologic characteristics between the low-risk and high-risk groups were evaluated using the chi-square test. All statistical tests were two-sided, and P < 0.05 was considered to be statistically significant.

Results

Radiomic signature construction and validation

Forty-five prognostic features were screened from 431 radiomic features by univariate Cox regression. The radiomic heat map allowed visualization of the association between these 45 radiomic features and PFS in patients in the training cohort (Fig. 1). Nine features of the 45 radiomic features were further selected using the LASSO Cox regression model (Fig. 2). Detailed descriptions and coefficients of the 9 features are listed in Table 3. A radiomics signature was constructed using these 9 features, and the risk score was calculated by the linear combination of selected features weighted by their respective coefficients (β): risk score = Correlation × (−2.48) + Correlation_HHL × (−2.17 × 10–1) + Kurtosis_HLH × (−3.40 × 10–3) + Median_HLL × (5.95 × 10–2) + Run Length Nonuniformity_HHL × (2.87 × 10–6) + Run Percentage_HLH × (3.33 × 10–2) + Short Run Low Gray Level Emphasis_LLL × 1.13 + Sum Variance_HHL × (−8.04 × 10–4) + Uniformity × (1.92 × 101).

Fig. 1

A heat map of 45 radiomic features screened based on univariate Cox regression. Each column corresponds to one patient in the training cohort, and each row corresponds to one normalized radiomic feature. The PFS bars represents the time of follow-up. The trend line indicates the trend of follow-up, calculated with using the linear regression model. The four radiomic groups are: group 1, first order statistics; group 2, shape- and size- based features; group 3, textural features; group 4, wavelet features.

Fig. 2

Construction of a radiomic signature predictive of progression-free survival using the LASSO regression model. (A) Tuning parameter (lambda) screening in the LASSO regression model. The partial likelihood deviance is generated versus log(lambda), and the lowest partial likelihood deviance corresponds to the optimal number of features. The dotted vertical lines correspond to the optimal lambda value according to the 1 standard error criteria and the minimum criteria. When 9 features remain, the partial likelihood deviance is the lowest. (B) The LASSO coefficient profiles of the 45 radiomic features. A vertical line is drawn at the value chosen by 10-fold cross-validation, which indicates that 9 non-zero coefficients were identified by the optimal lambda value. LASSO, least absolute shrinkage and selection operator.

Table 3

Nine prognostic radiomics features selected for further analysis using the LASSO Cox regression model.

Features	Filter	Descriptions	Coefficients (β)
Correlation	–	Correlation measures the gray level linear dependence between pixels at a specified location relative to each other.	−2.48
Correlation	HHL		−2.17 × 10−1
Kurtosis	HLH	Kurtosis measures the sharpness of the first-order histogram.	−3.40 × 10−3
Median	HLL	Median is the value that divides the upper and lower half of the sorted array of pixel values.	5.95 × 10−2
Run Length Nonuniformity	HHL	Describes the similarity of the lengths of runs throughout the image. This feature is low if the run lengths are similar.	2.87 × 10−6
Run Percentage	HLH	Describes the distribution and the homogeneity of runs of an image in a certain direction. This feature is very high if all gray levels have the run lengths of 1.	3.33 × 10−2
Short Run Low Gray Level Emphasis	LLL	Describes the joint distribution of short runs and low gray level values. This feature is high when the image has many short runs and lower gray level values.	1.13
Sum Variance	HHL	Sum Variance measures the gray-level co-occurrence matrix relationship to distribution of intensity with respect to variance. High Sum Variance corresponds to greater standard deviation of sum average.	−8.04 × 10−4
Uniformity	–	Uniformity describes the uniformity of the image.	1.92 × 101

Abbreviation: LASSO: least absolute shrinkage and selection operator.

With using the third quartile value as a fixed cutoff value, PFS was significantly stratified by the radiomic risk score in both the training (Fig. 3A, P < 0.0001) and the validation (Fig. 3B, P = 0.0233) cohorts. In multivariable Cox analysis, the radiomic risk score was found to be an independent prognostic factor for PFS in both the training (P < 0.001) and the validation (P = 0.045) cohorts (Table 2 and Supplementary Table 3). Finally, the radiomic risk score was found to be associated with PFS in all subgroups in the training cohort (Supplementary Fig. 1).

Fig. 3

Kaplan-Meier survival curves for low- and high- risk groups based on the radiomic risk score. The radiomic signature significantly stratified the progression-free survival in both (A) the training (P < 0.0001) and (B) the validation (P = 0.0233) cohorts. The tick marks represent censored observations.

Table 2

Cox regression analysis for patients in the training cohort.

	Univariate Cox regression			Multivariable Cox regression
	HR	95% CI	P value	HR	95% CI	P value
Age (years)
≤40 vs. >40	1.279	0.902–1.813	0.167
Sex
Male vs. Female	0.791	0.557–1.124	0.191
Histology
A vs. O/OA	0.722	0.656–1.339	0.937
Grade
WHO III vs. WHO II	2.489	1.719–3.604	<0.001	1.945	1.316–2.875	<0.001
Seizure
Yes vs. No	0.573	0.400–0.820	0.002	0.709	0.493–1.021	0.065
IDH Status
MUT vs. WT	0.527	0.368–0.754	<0.001	0.670	0.460–0.974	0.036
Radiomic Risk
High vs. Low (cutoff = −1.755)	2.871	1.980–4.163	<0.001	2.306	1.558–3.411	<0.001

Abbreviations: CI: confidence interval, HR: hazard ratio, MUT: mutant, WHO: World Health Organization, WT: wild-type.

The statistically  significant P value was displayed in bold font (P < 0.05).

Clinicopathologic characteristics

Two hundred and sixteen patients were included in the training cohort. There was no significant difference between the low-risk and high-risk groups with regard to age, sex, histology, whether or not seizure was present, or resection status. However, there were significant differences in WHO grade (P < 0.001) and IDH status (P = 0.001) between the two groups (Table 1). Eighty-four further cases were enrolled in the validation cohort, and no significant difference was found in age, sex, seizure history (not available for only 2 patients), WHO grade, or IDH status between the low- and high- risk groups. There was a significant difference in histology between the two groups (P = 0.046; Supplementary Table 2).

Table 1

Clinicopathologic characteristics of patients with lower-grade glioma in the training cohort (n = 216).

Variable		Low risk (n = 162)	High risk (n = 54)	P value
Age (years)	≤40	100	28	0.201
	>40	62	26
Sex	Male	98	33	0.936
	Female	64	21
Histology	A	59	21	0.745
	OA and O	103	33
Seizure	Yes	86	20	0.066
	No	76	32
	NA	0	2
WHO Grade	WHO II	126	28	< 0.001
	WHO III	36	26
IDH Status	Mutant	120	27	0.001
	Wild type	42	27

Abbreviations: A: astrocytoma, IDH: isocitrate dehydrogenase, NA: not available, O: oligodendroglioma, OA: oligoastrocytoma, WHO: World Health Organization.

Based on the Akaike information criterion, WHO grade, age, seizure, IDH status, and radiomic risk score were selected and integrated into a nomogram (Fig. 4A). The radiomic nomogram yielded a C-index value of 0.684 in the training cohort and 0.823 in the validation cohort. The radiomic nomogram also showed good calibration in both the training and validation cohorts (Fig. 4B and C). As the radiomic risk score was removed from the nomogram, the C-index reduced to 0.668 and 0.815 in the training and in the validation cohort respectively.

Fig. 4

Nomogram for prediction of progression-free survival and the calibration evaluation of the model. (A) A radiomic nomogram that integrates the radiomic risk score with the clinicopathologic characteristics in the training cohort. The C-index values are 0.684 and 0.823 in the training cohort and the validation cohorts, respectively. (B) The calibration curve for the radiomic nomogram in the training cohort. The diagonal gray line represents an ideal evaluation, while the red and blue lines represent the performance of the nomogram. A closer fit to the diagonal gray line indicates a better assessment. (C) The calibration curve for the radiomic nomogram in the validation cohort. WHO, World Health Organization.

The transcriptomic profiles of 47 patients in the CGGA cohort were analyzed to explore the genetic background of the radiomic signature. A heat map was constructed to illustrate the associations between the radiomic risk score, gene expression level, and clinicopathologic characteristics (Fig. 5A). Further GO analysis revealed significant associations of the radiomic risk score with biological processes, including the immune response, programmed cell death, NF-kappaB signaling, and vasculature development (Fig. 5B). Radiogenomic analysis in the TCGA cohort revealed that some of the processes positively associated with high-risk score exist both in the CGGA and TCGA cohorts, such as immune response, NF-kappaB signaling, and angiogenesis (Supplementary Fig. 2). To reveal the underlying malignant biological processes of each feature in the radiomic signature, we analyzed the negatively associated biological processes of Correlation, Correlation_HHL, Kurtosis_HLH and Sum Variance_HHL, given that these four features were negatively correlated with poor PFS. Similarly, the positively associated biological processes of the other five features were analyzed individually because they were positively correlated with poor PFS. Certain individual features were indeed associated with malignant biological processes, including epithelial-to-mesenchymal transition, positive regulation of cell division, the immune response, and positive regulation of cell migration (Supplementary Fig. 3 and Supplementary Fig. 4), all of which could contribute to poor PFS.

Fig. 5

Radiogenomic analysis of the radiomic risk score based on the 47 patients in the training cohort. (A) A heat map of the top 200 genes that were positively associated and the top 200 genes that were negatively associated with the radiomic risk score (P < 0.05; the absolute value of Pearson correlation coefficient > 0.4, Pearson correlation test). Each column corresponds to one patient and each row corresponds to one gene expression level. Associations between the radiomic risk score, gene expression level, and clinicopathologic characteristics are shown. (B) Functional annotation of the radiomic risk score. Gene ontology analysis of the genes that were positively associated with the radiomic risk score revealed that the radiomic risk score was significantly associated with biological processes, including the immune response, programmed cell death, cell proliferation, and vasculature development. WHO, World Health Organization.

Discussion

In this study, many radiomic features were extracted from T2-weighted MR images and a radiomic signature was identified as an approach for stratifying patients into low-risk and high-risk groups. These two groups had significantly different PFS in both the training and validation cohorts. Subsequently, the radiomic signature was integrated into a nomogram along with clinicopathologic characteristics, and an effective tool was constructed for individualized assessment of PFS in patients with LGGs. Moreover, subsequent radiogenomic analysis revealed that the radiomic signature was associated with biological processes, including the immune response, programmed cell death, cell proliferation, and vasculature development.

Previous studies have revealed that PFS can be predicted by several conventional observable MRI features (e.g., tumor contrast enhancement pattern,(Wang et al., 2016a) peritumoral edema (Wang et al., 2016b), and tumor location (Wang et al., 2017)) in certain subgroups of glioma. Additionally, a series of correlations has been found between PFS in patients with gliomas and voxel-based imaging features, such as the maximum relative cerebral blood volume value (Bisdas et al., 2009) and the minimum apparent diffusion coefficient (Romano et al., 2013) of a single voxel in the tumor region. In recent years, radiomic analysis has emerged as a powerful tool for building decision-support models based on high-throughput quantitative features extracted from medical images (Gillies et al., 2016). Radiomics-based PFS-predictive models have been reported for early-stage non-small cell lung cancer (Huang et al., 2016) and advanced nasopharyngeal carcinoma (Zhang et al., 2017), indicating the feasibility of predicting PFS in patients with LGGs by radiomic analysis.

In the current study, 45 of 431 radiomic features were screened by univariate Cox regression analysis, and 9 features that were strongly associated with PFS were further selected from these 45 radiomic features based on the LASSO Cox regression model. The radiomic signature constructed using these 9 features revealed adequate discrimination both in the training (P < 0.0001, log-rank test) and validation (P = 0.0233, log-rank test) cohorts. Intriguingly, most of the 9 prognostic features in the signature were also demonstrated to be able to capture prognostic information in previous studies. For example, it was reported that low Kurtosis and high Uniformity were related to poor outcomes (defined by a shorter time to progression) in patients with primary colorectal cancer (Ng et al., 2013), which is similar to our findings. Specifically, this phenomenon could be biologically interpreted using our radiogenomic analysis because Kurtosis was negatively associated with malignant biological processes, such as positive regulation of cell division and positive regulation of transcription (Supplementary Fig. 3D), while Uniformity was positively associated with mitotic nuclear division and nucleosome assembly (Supplementary Fig. 4D). Another study revealed Sum Variance and Correlation to be two of the most predictive features in differentiating long-term and short-term survival in patients with glioblastoma (Prasanna et al., 2017). Furthermore, the feature Median was found to be associated with PFS in patients with advanced nasopharyngeal carcinoma (Zhang et al., 2017).

In the context of precision medicine, an individualized PFS prediction model that can guide therapeutic strategies is essential. A nomogram is an approach that enables neuro-oncologists to estimate patient survival based on their clinical and biological profiles, which is an advance towards patient-tailored treatment in the changing landscape of neuro-oncology (Bredel, 2008). In the current study, a nomogram that integrates radiomics and clinicopathologic information was established for evaluation of PFS in patients with LGGs for the first time. Favorable results were achieved in two independent datasets (CGGA, C-index 0.684; TCGA, C-index 0.823), indicating that the nomogram is robust and potentially applicable in clinical practice.

Further radiogenomic analysis in our study suggested that the radiomic approach could potentially reflect biological processes and guide treatment for patients with LGGs. For example, the radiomic risk score was found to be positively correlated with certain malignant tumor processes, such as cell proliferation, cell adhesion, and vasculature development. Therefore, more aggressive therapeutic strategies are suggested for patients with high radiomic risk scores. Further, drugs that target the immune system, development of blood vessels, or the nuclear factor kappa B pathway might be helpful for patients with high radiomic risk scores. Although the hypotheses outlined above are preliminary and need to be prospectively evaluated in future studies, our present findings provide an approach for integrating imaging features, clinical characteristics, and genetic information to aid in clinical decision-making with regard to the management of LGGs.

The main limitation of this study is its retrospective design. The imaging protocols used were not fully consistent in that the imaging data were acquired by different MRI scanners in different centers. However, all the imaging data were normalized before extraction of features to reduce bias. Previous studies have already documented the robustness of extraction of radiomic features in terms of repeatability and reproducibility in test/re-test settings (Aerts et al., 2014; Balagurunathan et al., 2014; Fried et al., 2014; Grove et al., 2015; Leijenaar et al., 2013; Parmar et al., 2015). We expect that the performance of the radiomic signature will be further improved by more standardization of imaging data. Additionally, the biological processes underlying the radiomic signature and its components need to be validated in the future by prospectively designed studies.

In conclusion, the present study developed a radiomic signature as a non-invasive approach for preoperative evaluation of PFS in patients with LGGs. A radiomics-based nomogram and subsequent radiogenomic analysis could be useful in precision medicine and improve the therapeutic strategies used in patients with LGGs.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.