Cluster Analysis of DSC MRI, Dynamic Contrast-Enhanced MRI, and DWI Parameters Associated with Prognosis in Patients with Glioblastoma after Removal of the Contrast-Enhancing Component: A Preliminary Study

BACKGROUND AND PURPOSE: No report has been published on the use of DSC MR imaging, DCE MR imaging, and DWI parameters in combination to create a prognostic prediction model in glioblastoma patients. The aim of this study was to develop a machine learning – based model to ﬁ nd preoperative multiparametric MR imaging parameters associated with prognosis in patients with glioblastoma. Normalized CBV, volume transfer constant, and ADC of the nonenhancing T2 high-signal-intensity lesions were evaluated using K-means clustering. MATERIALS AND METHODS: A total of 142 patients with glioblastoma who underwent preoperative MR imaging and total resection were included in this retrospective study. From the normalized CBV, volume transfer constant, and ADC maps, the parametric data were sorted using the K-means clustering method. Patients were divided into training and test sets (ratio, 1:1), and the optimal number of clusters was determined using receiver operating characteristic analysis. Kaplan-Meier survival analysis and log-rank tests were performed to identify potential parametric predictors. A multivariate Cox proportional hazard model was conducted to adjust for clinical predictors. RESULTS: The nonenhancing T2 high-signal-intensity lesions were divided into 6 clusters. The cluster (class 4) with the relatively low normalized CBV and volume transfer constant value and the lowest ADC values was most associated with predicting glioblastoma prognosis. The optimal cutoff of the class 4 volume fraction of nonenhancing T2 high-signal-intensity lesions predicting 1-year progression-free survival was 9.70%, below which the cutoff was associated with longer progression-free survival. Two Kaplan-Meier curves based on the cutoff value showed a statistically signi ﬁ cant difference ( P ¼ .037). When we adjusted for all clinical predictors, the cluster with the relatively low normalized CBV and volume transfer constant values and the lowest ADC value was an independent prognostic marker (hazard ratio, 3.04; P ¼ .048). The multivariate Cox proportional hazard model showed a concordance index of 0.699 for progression-free survival. CONCLUSIONS: Our model showed that nonenhancing T2 high-signal-intensity lesions with the relatively low normalized CBV, low volume transfer constant values, and the lowest ADC values could serve as useful prognostic imaging markers for predicting survival outcomes in patients with glioblastoma.

therapy. 5 However, GBM frequently recurs as tumor cells infiltrate beyond the CEL components, leading to a dismal prognosis with an average overall survival rate of approximately 16 months. 6 Tumor cells that infiltrate past the CEL margin are visualized as hyperintense lesions on T2WI and are known as nonenhancing T2 high-signal-intensity lesions (NE-T2HSILs). 7 It is known that the preoperative T2-hyperintensity lesions that surround the CEL components of GBM are mainly composed of NE-T2HSILs and vasogenic edema. 7 The importance of the residual NE-T2HSIL has recently been recognized, and studies have found that NE-T2HSIL affects the prognosis of patients with GBM. [7][8][9] The preoperative prognosis of GBM in terms of the NE-T2HSIL, therefore, enables more aggressive management of the tumor if necessary, which could lead to survival benefits. 7 Advanced MR imaging detects various aspects of tumor pathophysiology and enables noninvasive visualization of the tumor. Because the survival and proliferation of tumor cells are highly related to angiogenesis and an increase in vascular permeability, 10 perfusion and diffusion imaging is expected to provide additional useful information about the NE-T2HSILs of GBM. Dynamic contrast-enhanced (DCE) MR imaging is a standard technique used to assess the integrity of the BBB on the basis of T1 enhancement. An important perfusion-related parameter measured in DCE MR imaging studies is the volume transfer constant (K trans ). K trans is defined as the rate at which the contrast agent leaks into the extravascular extracellular space per volume of tissue. 11 DSC MR imaging, another perfusion-weighted MR imaging technique, measures capillary perfusion on the basis of the susceptibility effect in T2*-weighted images. Normalized CBV (nCBV), the most common perfusion parameter used in DSC MR imaging studies, reflects the presence of blood vessels in each individual voxel. 12 Last, DWI evaluates the random motion of water molecules. Quantitative analysis of the DWI enables the calculation of the ADC, which is known to have an inverse correlation with tissue cellularity. 13 Extracting and combining MR imaging features from multiple modalities is labor-intensive, and achieving satisfactory diagnostic accuracy remains a major challenge. With recent developments in artificial intelligence, machine learning (ML) techniques have been steadily applied in glioma imaging studies and are revealing ways to solve these problems. ML techniques can efficiently process complex imaging data, identify meaningful disease patterns, and thus help radiologists make precise predictions about the progression of the disease. 14 Therefore, the application of ML techniques in imaging analysis will enable the effective integration of complementary imaging information obtained from using multiple MR imaging modalities.
Although previous studies have been conducted to combine several MR imaging modalities for the improvement of the prognosis of GBM, no report has been published on the use of DSC MR imaging, DCE MR imaging, and DWI parameters in combination to create a prognostic model. 15,16 Our study specifically focused on the NE-T2HSILs of GBM because of their importance as a site of recurrence and their relationship to prognosis, especially in patients with GBM who underwent GTR of the contrastenhancing (CE) components. Our aim was, therefore, to find nCBV, K trans , and ADC parameters associated with prognosis in NE-T2HSIL on T2-FLAIR after GTR of the contrast-enhancing components in patients with GBM using ML-based cluster analysis.

Patients
This retrospective study was approved by the institutional review board of the Seoul National University Hospital (IRB No. 1811-164-992), and the requirement for informed consent was waived. In this study, a total of 273 patients who were diagnosed with GBM from April 2010 to December 2019 were initially reviewed. Among the 273 patients, patients were selected according to the selection criteria outlined below.
The inclusion criteria were as follows: 1) patients older than 18 years of age, 2) diagnosed with GBM according to the 2016 World Health Organization Classification of Tumors of the Central Nervous System, 3 3) who had preoperative conventional MR imaging including 3D CE-T1WI and T2-weighted FLAIR imaging, 4) who had preoperative advanced MR imaging, including DCE MR imaging, DSC MR imaging, and DWI, and 5) who underwent standard treatment, which includes GTR followed by concurrent chemoradiation therapy with 6 cycles of adjuvant temozolomide. The exclusion criteria were as follows: 1) patients with loss of raw data (n ¼ 22), 2) lost to follow-up (n ¼ 41), 3) who had undergone partial resection (n ¼ 12), 4) who had undergone biopsy only (n ¼ 14), 5) who did not complete standard treatment (n ¼ 12), 6) who had unreadable data (n ¼ 28), and 7) a loss of clinical information (n ¼ 2). Under these inclusion and exclusion criteria, a total of 142 patients (85 men and 57 women; age range, 22-84 years) were ultimately enrolled in this study (Online Supplemental Data).
All included patients were divided into the progression group (n ¼ 113) or the nonprogression group (n ¼ 29) 1 year after the operation. The patients periodically underwent clinicoradiologic follow-up after completion of the standard treatment and were diagnosed with disease progression if they met at least 1 criterion of the Response Assessment in Neuro-Oncology (RANO) criteria. The RANO criteria are as follows: 1) $ 25% increase in the sum of the products of perpendicular diameters of enhancing lesions with the smallest tumor measurement, 2) any new lesion, 3) clear clinical deterioration not attributable to cause other than the tumor, and 4) clear progression of nonmeasurable disease. 17

Mask Segmentation and Advanced Image Processing
The 3D CE-T1WIs of all patients were registered and resliced into isometric voxels of 1 Â 1 Â 1 mm 3 to achieve spatial alignment and correct motion artifacts across consecutive images. The T2 FLAIR, nCBV, K trans , and ADC maps were then coregistered and resliced to the isovoxel CE-T1WI using rigid transformations with 6 df in the SPM package (Version 12; www.fil.ion.ucl.ac.uk/ spm/).
A 3D U-Net-based deep learning model was created using the open data set and images from several organizations, and the tumor segmentation masks were generated on the basis of 3D CE-T1WI and FLAIR images. Using the tumor masks, we analyzed nCBV, K trans , and ADC values in the NE-T2HSILs, and the segmentation of the NE-T2HSIL was validated by an experienced neuroradiologist (with 19 years of experience in neuro-oncologic imaging). The process of mask segmentation and image processing is shown in the Online Supplemental Data.

K-Means Clustering Analysis and Optimization of the Number of Clusters
Using the tumor-segmentation masks, we performed voxelwise Kmeans clustering in the NE-T2HSILs based on nCBV, K trans , and ADC maps. The K-means clustering module in the scikit-learn Python package (https://scikit-learn.org/stable/index.html) was used. K-means clustering is an iterative analysis that partitions a data set into K clusters, optimizing the similarities among data in a group while maintaining the greatest possible separation between different groups. K initial means are randomly selected, and each time new data are added to the program, they are assigned to the cluster with the closest mean. In the next step, the mean of the selected cluster is recalculated, and this process is repeated until all the data have been added to the program. 18 All voxels from all segmented masks in the training set were first combined and then divided into multiple clusters based on the similarities across data points in the same cluster and the differences across data points in different clusters. To determine the optimal number of clusters, we ran the program several times from 3 to 6 clusters. Finally, the best cluster number was chosen, which best discriminates the 1-year progression-free survival (PFS) in the training set on the basis of a univariate analysis.

Patient Population: Training and Test Sets
We randomly assigned all patients to a training or test set to prevent overfitting in the K-means clustering analysis. We divided the training set and test set at a ratio of 1:1. Each training set consisted of 71 patients, and the test set consisted of 71 patients. The patients in the training and the test sets were equally balanced with respect to the 2 prognosis groups. Next, the K-means clustering performance was evaluated using the parametric data of all patients in the test set.

Clinical Predictors and Outcome Definition
Clinical predictors were obtained from all patients' medical records, including sex, age, isocitrate dehydrogenase isozyme 1 (IDH1) mutation status, and O 6 -methylguanine-DNA methyltransferase (MGMT) promoter methylation status. The primary end point of this study was PFS. For those patients who were diagnosed with tumor progression according to the RANO criteria, PFS was calculated from the day of the operation until the day of progression. For those patients who showed no progression during the follow-up period, PFS was monitored at the time of the last follow-up MR imaging examination, and PFS was estimated by survival analysis.

Variable Selection
The selection of variables that are significantly relevant to the patient's survival outcome is extremely important in developing a prognostic model. The variables we analyzed in this study included the volume fractions of each parametric cluster in the NE-T2HSIL and several clinical predictors, including sex, age, IDH1 mutation status, and MGMT promoter methylation status. We first applied these variables in the univariate analysis and then included the variables that were found to be statistically significant in the multivariate analysis.
MR imaging protocol, imaging processing, and statistical analysis are summarized in the Online Supplemental Data.

Demographic Data of the Study Population
A total of 142 patients who were diagnosed with GBM were enrolled in our study. The demographic data of all enrolled patients are summarized in Table 1. During the follow-up period, 113 (79.6%) patients experienced progression, and 29 (20.4%) patients did not experience progression 1 year after the operation. Patient demographics in the training and test sets are summarized in the Online Supplemental Data.

Optimization of the Number of Clusters Using the Training and Test Sets
To determine the optimal number of clusters, we randomly divided all enrolled patients into a training set and a test set at a ratio of 1:1. The initial K-means clustering was performed on the parametric data of all patients in the training set. All voxels from the segmented NE-T2HSIL masks were divided into clusters on the basis of the nCBV, K trans , and ADC maps. A range of cluster numbers was entered, and finally, the cluster number of 6 (k ¼ 6) was chosen because it best discriminated the 1-year PFS in the training set on the basis of univariate analysis (Online Supplemental Data). Receiver operating characteristic analysis was performed for both the training and test sets, and the areas under the curve were 0.69 and 0.67, respectively. The area under the curve values suggest that the division into 6 clusters is an acceptable distinction (Online Supplemental Data).

Class 4r Correlated with PFS in All Patients
After we validated the number of clusters, the training and test sets were combined for survival analysis. The 2D clustering plots based on the nCBV, K trans , and ADC values are demonstrated in the Online Supplemental Data. The average percentages of voxels in each of the 6 clusters were 4.99%, 18.42%, 8.41%, 27.34%, 29.24%, and 19.60% in numeric order. The parametric information for each cluster is shown in the Online Supplemental Data. Among the 6 clusters, "class 4" was the cluster with the relatively low nCBV and K trans values and the lowest ADC values in the segmented NE-T2HSIL mask. We found that a higher class 4 volume fraction was associated with a shorter PFS in patients with GBM. We defined the term "class 4r" as the volume fraction of cluster 4 in the total NE-T2HSIL of each patient. A video clip with a 3D clustering plot can be accessed in the Online Supplemental Data.

Stratification Based on Class 4r and Survival Analysis in All Patients
Our study suggested that class 4r in the NE-T2HSIL was associated with 1-year PFS prediction in GBM. The optimal cutoff of class 4r for stratifying the shorter and longer PFS groups was 9.70% of the NE-T2HSIL volume fraction. The cutoff value distinguished the 2 PFS groups with a significant difference in the log-rank test (P = .037). Thus, our results showed that patients with GBM with a class 4r of ,9.7% had a significantly longer PFS time than patients with a class 4r of $9.7%. The Kaplan-Meier survival curves and a table with the "numbers at risk" are shown in Fig 1.

Univariate and Multivariate Analysis in All Patients
The univariate Cox proportional hazard analysis showed that clinical variables such as age, IDH1 mutation status, and MGMT promoter methylation status were independent predictors of 1-year PFS in patients with GBM. In addition, among the MR imaging parametric variables, only an increase in class 4r was identified as an adverse predictor of PFS (Table 2). In the multivariate Cox proportional hazard analysis, a statistically significant difference in PFS between patients with low and high class 4r values (hazard ratio, 3.04; 95% CI, 1.00-9.10; P value ¼ .048) was observed, which was independent of prognostic genetic factors, including IDH1 mutation status and MGMT promoter methylation status ( Table 3). The concordance index value for PFS was 0.699 (standard error, 0.025), indicating that the volume fraction of regions with the relatively low nCBV and K trans values and the lowest ADC value is a moderately good

DISCUSSION
In this study, we developed a multiparametric prognostic model for patients with GBM who received GTR of the CEL mass. We found that a higher volume fraction of voxels with relatively low nCBV and K trans values and the lowest ADC value in NE-T2HSIL showed a significant association with worse GBM prognosis: A low K trans value indicates preserved vascular permeability due to subtle BBB damage, 19 a low nCBV value implies a lack of tumor angiogenesis, 20 and a low ADC value reflects an increase in tumor cellularity. 13 Thus, regions of relatively low nCBV and K trans values and the lowest ADC value represent the hypoxic hypercellular regions within the NE-T2HSIL. We believe that a high content of these regions in the NE-T2HSIL can be a significant prognostic factor in predicting the survival outcomes of GBM. Few studies have used DSC MR imaging, DCE MR imaging, or DWI to evaluate MR imaging-derived parameters in the NE-T2HSIL to predict the prognosis of GBM. Derived from DWI, a low ADC value can be used as a marker indicating regions with a high content of microscopic tumor cells. 13 On T2WI, however, these microscopic tumor cells are in a mixture of vasogenic edema, where high ADC values are present. 21 Hence, the low ADC level helps locate the regions with high microscopic tumor cell content from the vasogenic edema mixture. 21,22 Consistent with our finding, Lee et al 21 also suggested that the minimum ADC value of NE-T2HSIL may indicate the infiltration of neoplastic cells in peritumoral edema.
Derived from DCE and DSC MR imaging, K trans and nCBV values reflect the extent of perfusion and tissue vascularity. Our finding of low K trans and low nCBV levels can be interpreted as hypoxic tumor regions, which are key to the highly infiltrative  nature of GBM. [23][24][25] Local hypoxic tumor conditions prompt the activation of hypoxia-inducible factor 1, leading to vascular endothelial growth factor production. The vascular endothelial growth factor then promotes angiogenesis and invasion of the normal brain parenchyma, which ultimately lead to tumor progression. 25,26 However, our results contradict those of previous studies. For example, Kim et al 27 recently suggested the potential of a high K trans value within the NE-T2HSIL in predicting the prognosis of GBM. According to that study, patients with GBM with higher K trans values in the 99th percentile had a worse prognosis. With respect to rCBV, Jain et al 9 reported that higher mean rCBV values in NE-T2HSIL are associated with a worse prognosis in GBM. The reason for the discrepancies between our current data and these 2 previous studies is that these 2 reports did not include the ADC parameter, a tissue cellularity imaging marker. We believe that in the NE-T2HSIL, a low ADC value resulting from hypercellularity is more important than hypervascularity or hyperpermeability for the prognosis; additionally, hypoxic hypercellular areas seem to have an aggressive potential.
We believe that this discrepancy is caused by the following 2 factors: First, while the previous studies used only 1 or 2 MR imaging parameters, our study included 3 MR imaging parameters (ADC, K trans , and nCBV) at the same time to create a prognostic model of GBM. Therefore, we detected the most prognostic components reflecting tumor cellularity, angiogenesis, and vascular permeability in patients with GBM. Second, in previous studies, the NE-T2HSILs were analyzed as a whole, with a representative parametric value being determined independent of their composition. Our study, however, divided the NE-T2HSILs into edematous parts and microscopic tumor cells on the basis of DSC MR imaging, DCE MR imaging, and DWI, taking into account the different components that make up the NE-T2HSILs in GBM.
Despite the discrepancy between our results and previous studies, we suggest that regions with low parametric values within the NE-T2HSILs may relate to the early stages of pseudopalisade development. Pseudopalisades are unique pathologic features that differentiate GBM from low-grade gliomas. 28 Although the detailed mechanism for the formation of pseudopalisades remains unclear, a new model proposed by Rong et al 26   insight into our results. According to this model, the high incidence of intravascular thrombosis in patients with GBM can directly trigger or spread hypoxia within the tumor. After vascular insult, however, the BBB and vascular structure remain intact to some extent. This feature keeps permeability relatively low until pseudopalisading cells migrate away from the vascular occlusion site, secrete proangiogenic factors, and develop central necrosis. 26 Therefore, tumor voxels with relatively low nCBV and K trans values and the lowest ADC value could represent NE-T2HSILs with large numbers of pseudopalisading cells in the early stages of development.
Our study has several limitations. First, it is based on a retrospective design that is inevitably prone to selection bias. However, we tried to minimize the selection bias by collecting most of the data from April 2010 to December 2019 in the Seoul National University Hospital. Second, all patient data were collected from a single medical center. Future studies with multicenter data could help validate and generalize our results. Third, the size of our study population was relatively small, especially the number of samples in the nonprogression groups, so a future study is warranted in a large population for the validation of our results. Last, although we have suggested a plausible correlation between pathologic features and imaging features, this study lacks pathologic validation because the NE-T2HSILs are not the usual target of GBM.

CONCLUSIONS
The multiparametric model showed that the NE-T2HSILs with the relatively low nCBV and K trans values and the lowest ADC value could serve as useful imaging markers in predicting survival outcomes in patients with GBM. Our investigation could help radiologists locate hypoxic, hypercellular regions of the NE-T2HSIL, consider more aggressive treatment to prevent early tumor recurrence, and ultimately improve the overall survival of patients with GBM.
Disclosure forms provided by the authors are available with the full text and PDF of this article at www.ajnr.org.