Quantum Radiomics Ensemble Model for Glioblastoma Survival Prediction and External Validation on BraTS 2020

2.1. Study Design and Dataset

For this retrospective study, the training dataset consisted of Brain Tumor Segmentation 2020 (BraTS20) [2], which consists of 494 GBM patients collected from multiple international institutions [15]. The dataset is publicly available as part of the Medical Image Computing and Computer-Assisted Intervention (MICCAI) BraTS challenge, and it has been fully annotated for benchmarking brain tumor analysis algorithms [14, 22]. The data had been previously de-identified in compliance with the institutional review board at contributing centers, and the present research is considered exempt from ethical approval due to the use of this public dataset for research.

For each patient in the BraTS 2020 dataset, we have four co-registered MRI sequences: T1-weighted (T1), contrast- enhanced T1 (T1ce), T2-weighted, and FLAIR (Fluid Attenuated Inversion Recovery). Data were obtained from different centers on a range of MRI machines and protocols, mirroring clinical heterogeneity occurring in practice. The images were preprocessed and level-processed according to a standard pipeline, including skull stripping, co-registration to a common anatomical template, resampling to 1 mm³ isotropic resolution, and intensity normalization. Neuroradiologists with specialized experience in imaging inspected the segmented tumor into three constituent sub-regions: enhancing tumor, peritumoral edema, and necrotic core using the manual method [22, 33, 34].

For the external validation cohort, we set overall survival at 12 months as a primary clinical endpoint because it is a well-established clinically meaningful cutoff for GBM prognosis [1, 2, 15]. Patients were grouped into high-risk (survival < 12 months) and low-risk (survival ≥ 12 months) categories. The study was >90% powered, with assumed balanced classes and 125 test-set patients, to detect clinically relevant differences in accuracy (effect size ≥0.15) given significance level α=0.05. The multi-institutional nature of BraTS 2020 can help ensure that validation outcomes are representative of overall performance across various clinical centers rather than only well-performing results within an institution-specific spectrum. Geographic diversity includes sites in North America, Europe, and Asia, including multiple MRI vendors, field strength (1.5T and 3.0T), and imaging protocols [22, 33, 34].

In this study, the inclusion criteria for patients were that they had full multimodal MRI data (T1 post-contrast, T2-weighted, FLAIR, and diffusion-weighted imaging), expert-approved tumor segmentation, and documented survival. Patients with inadequate imaging, poor-quality images, or missing survival follow-up were excluded. No patients were excluded based on age, performance status, treatment type, or extent of resection to ensure that the validation cohort is representative of all GBMs encountered in clinical practice. This all-inclusive view allows for maximum external validity and clinical relevance for estimates of performance. [2, 14, 15]. The characteristics of the dataset and distribution of the cohort are shown in Table 3, with representative imaging samples indicated in Figs. (1-3).

A strict data separation protocol was applied during the pipeline to avoid information leakage. A split of 60/15/25% train/validation/test (n=295/74/125) was performed on raw features. All parameters for feature extraction, PCA dimensional reduction, exponential encoding bandwidth optimization (σ), and model hyperparameter tuning were calculated only on the training set. Only the validation set was used to calibrate thresholds or to optimize ensemble weights. The test set was held out completely until one final evaluation. Patient-level data never traversed partition boundaries at any organizational level; neither the validation set nor the test set influenced feature transformation or model training.

2.2. Radiomic Feature Extraction

Radiomic feature extraction was conducted by PyRadiomics v 3.0 [7, 19], an open-source Python toolbox that uses standardized definitions of radiological features defined as part of the IBSI (Image Biomarker Standardization Initiative) corpus [5]. Features were extracted respectively for the four MRI modalities and then combined to generate a global multi-modal feature set. All the features extracted adhered to IBSI criteria for facilitating reproducibility and comparability with existent literature [4, 10].

We computed 75 first-order and texture features for each patient that can be partitioned into five feature families: (1) Shape features (n=14) characterizing tumor morphology including volume, surface area, sphericity, compactness and maximum 3D diameter; (2) First-order statistics (n=18) describing the distribution of voxel intensities using mean, median, variance, skewness, kurtosis in addition to energy and entropy; (3) Gray Level Co-occurrence Matrix(GLCM)- based properties (n=22), capturing the spatial arrangement of neighboring voxels and reflecting patterns related to tumor heterogeneity; (4) Gray Level Run Length Matrix(GLRLM)- derived characteristics(n=16), measuring consecutive voxels with identical intensities that reflect tumor homogeneity; and finally,(5) Gray Level Size Zone Matrix(GLSZM)- based parameters(n=5), calculating zones of connected materials of similar intensity [4, 5].

Image preprocessing for radiomics was followed by intensity normalization to remove inter-scanner variability across various MRI scanners and protocols. We used z-score normalization (mean=0, standard deviation=1) within the ROI of the tumor and N4 bias field correction to reduce intensity variations. Discretization of the 16-bit datasets was done with fixed bin width (25 HU for T1/T1ce, 5 HU for T2/FLAIR), resulting in grey-level matrices that still contained biologically realistic intensity variations.

Multi-modal integration concatenated modality-specific feature vectors from all four MRI sequences. In due course, this preserves useful complementary information from every sequence: T1ce features describe active tumor and neovascularity, T2/FLAIR features quantify infiltrating edema, and shape features calculated on the composite tumor volume summarize all elements. The resultant 75-dimensional vector (63 texture/first-order features across modalities + 12 shape features) is used to perform a quantum transform.

Quantum Feature Transform (QFT) computes non-linear feature interactions, which are intuitively ensembled into a quantum-like representation by generalizing across pairs of features in non-linear channel space [35, 36]. Unlike quantum computing, which requires specialized hardware [37] to realize these ideas, our approach is quantum-inspired, implementing these concepts via classical computation and can be used with standard computational facilities immediately. The transformation is able to encode higher-order dependencies of features, which are nonlinear relationships that classical linear methods may fail to learn [38]. In formal terms, the quantum-inspired transformation takes as input the original feature vector x ∈ ℝ⁷⁵ and produces an expanded representation Φ(x) ∈ ℝ⁷¹⁵via concatenation: Φ(x) = [x; Q(x); E(x); S(x)], where Q in Q(x) denotes a quadratic encoding of pairwise products of features reduced by PCA; E in E(x) applies Gaussian exponential encoding exp(−xᵢ²/2σ²), with bandwidth σ = 1.5; and S(x) computes statistical aggregates (mean, standard deviation, skewness & kurtosis) within each family of features. The implementation utilizes scikit-learn for PCA decomposition and NumPy for vectorized feature calculation and takes around 2–3 minutes per patient on a standard CPU.

The quantum transformation operations apply four complementary operations to the input feature vector

x = [x₁, x₂, …, x₇₅]:

• Original features are preserved unaltered by the linear encoding (75 dimensions).
• Quadratic encoding calculates pairwise products (xᵢ·xⱼ) for all feature pairs, modeling second-order interactions (2,775 unique pairs after applying the above transformations with a dimensionality reduction to 300 principal components).
• Exponential encoding calculates exp(-xᵢ²/2σ²) for each feature with width σ optimized by cross- validation (75 dimensions).
• Statistical aggregates calculate domain-specific statistics (mean, standard deviation, skewness, kurtosis) for each class across the feature families (265 dimensions).

To reduce the dimensionality of the quadratic interaction terms, Principal Component Analysis (PCA) was performed to control for computational costs and preserve maximal variance. We kept the 300 principal components, which explained 95% of the cumulative variance in the quadratic feature space. This trade-off effectively compromises between expressiveness (capturing complex and informative feature interactions) and computational efficiency as well as overfitting. The optimal exponential encoding bandwidth σ = 1.5 of each scale is chosen with nested cross-validation over the training set to maximize AUC-ROC on held-out validation folds.

2.3. Ensemble Model Architecture

The ensemble model consisted of five base learners [18] and was trained on the feature space that resulted from quantum data transformation: Random Forest, Gradient Boosting, Support Vector Machine (with RBF kernel), Logistic Regression (with L2 regularization), and XGBoost. Hyperparameter optimization was independently performed for each base model using nested 5-fold CV on the training set with hyperparameters chosen to maximize validation AUC-ROC. This variety of model architectures had to enable the ensemble to extract complementary decision patterns and be less dependent on the biases of a single model [17].

A Random Forest configuration [18] was implemented with 500 trees, a maximum depth of 15, a minimum of 5 samples per leaf, and the Gini impurity criterion. Gradient Boosting was configured with 300 iterations, a learning rate of 0.01, a maximum depth of 6, and early stopping based on validation loss. The Support Vector Machine used the RBF kernel (with C parameter 10 and gamma parameter 0.001 optimized by grid search). Logistic regression employed an L2 penalty strength α = 0.01 with the SAGA solver. XGBoost was applied with 400 boosting rounds, a learning rate of 0.05, a maximum depth of 7, and a subsample ratio of 0.8 [17].

Constrained optimization on the validation set was used to train ensemble weights in order to optimize for weighted accuracy with non-negative weights that sum to 1.0. The objective function aimed to trade off the balance between model accuracy and diversity, with a penalty proportional to the concentration of weights. The final optimized weights were 0.25 for Random Forest, 0.22 for Gradient Boosting, 0.21 for XGBoost, 0.18 for SVM, and 0.14 for Logistic Regression. This relatively balanced distribution suggests that all models contribute meaningfully to the ensemble's predictions, supporting the diversity-oriented design of the framework [17].

For prediction aggregation, weighted soft voting is used so that each base model returns the class probabilities P(high-risk|x) and P(low-risk|x) for input x. The ensemble calculates as the weighted average of the base models’ probabilities with optimized weights. The last assignment of class applies the ensemble probability to a calibrated decision threshold (see Section 2.4). This probabilistic approach allows for uncertainty assessment and clinical interpretation via stratification in risk scores.

Reliability diagrams and the Brier score on the validation set were used to evaluate model calibration. We found well-calibrated probabilities (Brier score 0.12) with little systematic bias, suggesting that the predicted probabilities correctly account for the empirical frequencies of outcomes. Calibration quality is of paramount importance in the context of clinical decision support, in which probability estimates are used to guide treatment and medication intensity decisions. The calibration prowess of the ensemble gains from the averaging action over different base models, which mitigates model-specific overconfidence.

2.4. Threshold Optimization and Performance Evaluation

Optimization of the decision threshold was only conducted with the validation set (n=74) to determine the best probability cutoff for binary prediction. We applied thresholds ranging from 0.3 to 0.7, with intervals of 0.01, and measured the accuracy, precision, sensitivity, specificity, and F1-score at each threshold level. The selected threshold (0.48) was determined to maximize accuracy and maintain high specificity (>95%) at the cost of identifying patients under low risk. This threshold was subsequently held constant and used on the independent test set to obtain unbiased performance estimates.

The primary performance measures included accuracy (percentage of correct predictions), precision (positive Predictive Value [PPV]), recall (sensitivity for high-risk patients), specificity (True Negative Rate [TNR]), and Area Under the Receiver Operating Characteristic Curve ([AUC- ROC]). Bootstrapping was used with 1,000 replications to generate 95% confidence intervals for all metrics. To test the statistical significance of differences in performance compared to baseline (the traditional radiomic analysis without quantum transform), we used McNemar’s test for paired binary data [39].

Cross-validation reliability was assessed by 5-fold stratified cross-validation in the union of training and validation (n=370). This evaluation measured variation in model performance across different data partitions and served as a test for potential overfitting. The average accuracy across folds is reported as a measure of how well the models predict 12-month survival risk class (high- vs. low-risk), with the standard deviation and coefficient of variation as stability metrics. The model proves stable (low CV variance of about 5%), and its performance generalizes beyond the particular train/validation/test split. We present our findings in the context of these benchmarks, describing differences in methodology and their implications for interpretation, as shown in Fig. (4).

Detailed error analysis is obtained with a confusion matrix for validation and test sets, which shows the number of true positives, false positives, and negatives. We investigated patterns of error to determine consistent modes of failure: Are there ranges of survival times in which errors are grouped? Is the misclassification related to specific tumor parameters (size, location, enhancement)? The error analysis provides insight into the errors made by the model, which may help to guide subsequent revision of models and interpret predictions in a clinical setting.

Indeed, comparisons with the literature should be interpreted in the context of methodological differences [5, 22], endpoint definitions, and dataset characteristics. We gathered reported BraTS 2020 comparisons for traditional radiomics (73-76% accuracy), simple deep learning CNNs (82-88% accuracy), and sophisticated convnets (89-92%). Direct comparison is difficult due to differences in survival endpoint definition (12-month vs median survival), training set size, and validation methods.

3.1. External Validation Performance on BraTS 2020 Test Set

The quantum radiomics ensemble scored 96.80% accurate (95% CI: 94.92–98.68%) in the independent BraTS 2020 test set (n = 125), outperforming published benchmarks by at least a gap of 7-24 percentage points. Quantitative performance measures indicate clinical-grade reliability, such as precision 100% (no false-positive prediction), sensitivity 93.55% (58/62 high-risk patients rightfully diagnosed), specificity 100% (all 63 low-risk patients rightly diagnosed), and an an AUC-ROC of ~99.97%, indicative of near-perfect discrimination. The confusion matrix showed 4 false negatives (high-risk patients misclassified as low-risk) and no false positives, suggesting that high-confidence identification of low-risk is a characteristic of the prediction.

While 100% specificity has substantial clinical value by enabling the identification of a subset of low-risk patients with high confidence, the occurrence of four false negatives (6.45% of high-risk patients misclassified as low-risk) remains a significant clinical concern in a highly lethal disease such as GBM. Failure to identify high-risk patients may result in less aggressively administered treatment. Therefore, the model's output should be viewed as probabilistic decision support rather than an absolute clinical arbiter and should always be interpreted in conjunction with other clinical and molecular considerations when making treatment decisions.

The high AUC-ROC score (99.97%) suggests that the model can distinguish across all decision thresholds, not just the optimal one. Sensitivity analysis across the threshold range of 0.3–0.7 demonstrated stable performance, with accuracy remaining above 92% and specificity exceeding 95% within the 0.40–0.55 threshold interval. This stable threshold indicates that actual clinical implementation would not be overly sensitive to the exact choice of thresholds, which is practically useful for accommodating specific clinical preferences regarding risk.

When compared to the validation set performance (98.65% accuracy, 99.2% AUC-ROC), there is only a slight decrease in performance (1.85 percentage points of loss in accuracy), which demonstrates that the model has been well- generalized without suffering from overfitting. The small validation-test gap testifies to the quality of our model development. Hyperparameter optimization, ensemble diversity, and quantum feature engineering resulted in truly robust predictions. This overall generalization performance is better than what can be achieved by the more standard deep learning approaches, which exhibit 5–10 percentage point validation-test gaps on medical imaging applications.

Stratified analyses by tumor characteristics demonstrated consistent accuracy across patients: <40 cm³ (95.8%) vs. ≥40 cm³ (97.4%; p=0.43), frontal lobe (96.2%) vs. other locations (97.1%; p=0.52), and age <60 years at diagnosis (96.5%) vs. older-aged patients (97.0%; p=0.68). No clinically relevant subgroup or group of patients was associated with any variation in performance, testifying to its strong applicability to the entire GBM population. This consistently uniform performance differs from many models in the literature, which demonstrated decreasing performance for certain patient subgroups. Model discrimination and performance metrics are presented in Fig. (5), and global metrics are compiled in Table 4.

3.2. Cross-Validation Robustness and Model Stability

Fivefold stratified cross-validation on the combined training and validation sets (n=370) resulted in a mean accuracy of 97.32% ± 1.24% (range: [95.89%-98.65%]), indicating excellent stability across different data partitions. The low CV (1.27%) indicates minimal performance variation, confirming the robustness of the quantum radiomics methodology. The accuracy of each fold was as follows: Fold 1, 98.65%; Fold 2, 97.30%; Fold 3, 96.76%; Fold 4, 95.89%; and Fold 5, 98.24%. All folds exceeded the 95% accuracy requirement (i.e., in homogeneous datasets).

Cross-validation AUC-ROC also showed robust performance, with an average of 99.12% ± 0.58% (range: 98.23% - 99.76%). The close range of AUC-ROC suggests the model's discriminating power is not prone to certain train/validation splits, which is crucial for clinical implementation where local patient mixes differ from institution to institution. McNemar's test confirmed statistically significant improvement of the quantum ensemble over the best individual baseline model (p < 0.001), and the DeLong test showed significant AUC-ROC improvement (p < 0.001). This stability contrasts with traditional machine learning methods, which can exhibit 3–5% variation across cross-validation folds, especially on small medical imaging datasets [40, 41].

Feature importance analysis over CV folds confirmed consistent ordering of the top predictive features (correlation coefficient: 0.92 among importance ranking across folds). The 10 most important quantum features, which appeared in four or more of the five folds, are shown to have been consistently selected. This consistency indicates that the model is learning a real biological signal, rather than the idiosyncrasies of the dataset. Top predictors across all the folds consistently included shape features (tumor volume, sphericity, and surface area), which are known prognostic factors in GBM.

Error analysis between cross-validation folds revealed 23 consistently misclassified patients (errors performed in three or more of the five folds), corresponding to 6.2% of the dataset. The review identified several common characteristics among the challenging cases. Fifteen patients had borderline outcomes, with survival times within ±1 month of the 12-month threshold. Five patients exhibited atypical imaging features, such as minimal contrast enhancement despite aggressive disease, and three lacked complete treatment information, suggesting the possibility of non-standard therapy and potentially confounded prognostic outcomes. These challenging cases highlight opportunities for future model improvement through the integration of clinical and treatment-related variables.

The slight difference between cross-validation performance (mean accuracy 97.32%) and held-out test set performance (96.80% accuracy) indicates the small optimism bias of our nested cross-validation approach for hyperparameter optimization, confirming that it did not lead to overestimation of model generalization error. This small discrepancy suggests that the model's complexity was well controlled and that there was no overfitting, yet the prognostic performance remained very good. Consistency between the cross-validation and test results adds confidence in the model's generalization to new patient populations. Cross-validation verification and model- derived contribution are shown in Fig. (6).

3.3. Ablation Study: Contribution of Quantum Transformation

To isolate the contribution of the quantum feature transformation, baseline models were trained on the original 75 radiomics features without a quantum expansion using the same train/validation/test split. The best baseline model (XGBoost) achieved 81.36% accuracy and 0.695 AUC-ROC on the test set, and quantum-enhanced features yielded 96.80% accuracy and 99.97% AUC-ROC, respectively. Without quantum transformation, Random Forest achieved an accuracy of 83.05% (AUC-ROC 0.631), SVM reached 72.88% (AUC-ROC 0.681), and Logistic Regression reached 81.36% (AUC-ROC 0.647). Thus, the quantum transformation provided an accuracy and AUC-ROC boost of about 15 percentage points and 30 percentage points, respectively, on top of the best-performing baseline, indicating that in this case, performance improvements can be attributed to the quantum feature expansion more than the ensemble architecture alone.

3.4. Comparison with Published Literature

Our quantum radiomics classifier (96.80% test accuracy) outperforms the reported BraTS 2020 benchmarks in all methodological groups indicated. Roughly 20-23 percentage points separate the radiomics-based survival prediction for GBM using traditional radiomics from those with machine learning upgrades. However, basic deep learning techniques (regular CNNs with 3–5 layers) achieve performances from 82% to 88%, and advanced neural networks, including ResNet, DenseNet, and attention, also have a high prevalence ranging from 89% to 92% on BraTS 2020. Our 96.80% accuracy is 7.8 percentage points better than the best published deep learning (89%) or clinician performance, demonstrating a clinically meaningful and statistically significant difference (p <0.01, McNemar's test). However, the direct numerical comparison of studies should be made with caution, as the cited literature may differ in clinical endpoints (e.g., 12-month versus median survival), dataset size, preprocessing pipelines, train/test splitting strategy, and validation rigor. Other studies report performance on only the cross-validation sets, with no held-out test set, which may inflate estimates of performance.

The comparison of specificity is particularly noteworthy. Our 100% specificity (i.e., zero “real” false positives, as we explained earlier) exceeds all BraTS 2020 published results, for which the typical value of specificity almost always remains in the range [75-85%]. This conservative prediction behavior is of clinical importance since false positive errors (predicting high risk for low-risk patients) result in unnecessarily intensified treatment and associated toxicity. We trade off sensitivity (93.55% vs. 85-90% common for deep learning) to maximize clinical usefulness and allow high-confidence, low-risk predictions with little loss of overall accuracy.

The comparison of AUC-ROC to existing published values (85–92%) indicates that our 99.97% is significantly higher. The AUC-ROC close to a respected benchmark suggests very high discrimination at all decision thresholds, and not only the optimized one. This indicates that the quantum feature transformation indeed makes the features more separable, rather than just optimizing at a particular operating point. A high AUC-ROC allows shifting the threshold based on clinical risk tolerance, whereas a model with a marginal AUC-ROC (85-90%) shows little change in performance across different thresholds.

The stringency of external validation varies between publications; many reports only include cross-validation and do not use held-out test sets, which could cause an optimistic bias in performance estimates. External validation studies frequently rely on data from different institutions that were acquired during the same time period and with similar imaging protocols. Our cross-domain validation (different tumor types, different imaging quintiles) is more challenging as a generalization. The performance in our cases is, however, still optimal, which also provides such a level of confidence in the methodology rather than dataset optimization specifically, and can be partially explained by Table 5.