Radiomics-based patient inclusion model improves clinical trial performance

Patient populations

This study was approved by the University of South Florida Institutional Review Board. The analysis includes patients who participated in the TH CR-406/SARC021 multicenter clinical trial of Doxorubicin plus Evofosfamide (Dox+Evo) versus Doxorubicin alone (Dox) in locally advanced, unresectable or metastatic soft-tissue sarcoma. Full trial protocol and results were published by Tap et al. (6). A total of 640 patients were enrolled. The primary endpoint of the trial was OS. CT images obtained prior to treatment were available for analysis in 588 patients.

Patient data and CT images

Patient covariates and CT image were obtained from the Sarcoma Alliance for Research through Collaboration (SARC). Patient data and CT images, obtained prior to treatment, was available for 588 patients. The CT images were uploaded into HealthMyne Quantitative Imaging Decision Support (QIDS) software (https://www.healthmyne.com, Madison, WI), where a radiologist with 10 years (S.F.) identified and segmented all visible lesions. 346 patients were found to have at least one lesion in the lung, the most common metastatic site in the considered cohort (followed by liver lesions, identified in 106 patients), as anticipated (7). Only lung patients were included in the study to enable comparison of image features between individuals, and hence the use of radiomics. Of these patients, 303 had contrast enhanced CT scans of the lung which could be analyzed, and these were used for quantification. This total cohort of 303 patients used in this study was randomly divided 70:30 into training and test sets using the sample function in R version 4.0.2.

Radiomic feature extraction

The images were pre-processed as described in the Supplementary Methods. For each patient, a total of 166 features were calculated for each patient using standardized algorithms from the Image Biomarker Standardization Initiative (IBSI) v5 (12). The radiomic features included statistical, histogram, shape & size, Grey Level Cooccurrence Matrix (GLCM), Grey Level Run Length Matrix (GLRLM), Grey Level Size Zone Matrix (GLSZM) and Neighboring Grey Tone Difference Matrix (NGTDM) features, as well as 16 peritumoral features as described before (13). Laws and Wavelet features were not extracted due to their poor reproducibility (14). As standard in radiomic studies (15), to ensure the radiomic signatures provide additional information compared with tumor volume, the features strongly correlated to volume (Pearson | r|>0.8) were excluded from further analysis, while volume itself was included. Spatial stability of the features was assessed, and unstable features excluded. 100 features remained following the exclusion, as detailed in Supplementary Table 1.

Feature Selection

The goal of this analysis was to identify the radiomic features and patient covariates differentially associated with OS in the two treatment groups, which was the primary endpoint in the original trial (6). First, univariable Cox proportional hazards regression analysis was used to assess the degree and direction of statistical association of each feature and covariate with post-treatment OS, separately in Dox and Dox+Evo arms. For each arm, features and covariates were considered promising in either of the two scenarios: (i) They showed association (p<0.05) with survival in one group AND no association (p>0.20) in the other group, or (ii) they showed significant association (p<0.05) with survival in both groups in opposite directions (HR >1 in one group and <1 in the other). Multiple comparison correction was not used, as the statistical significance assessment was carried out in the next step.

Correlation between the remaining features was calculated (Pearson’s correlation coefficient for continuous and Chi Square independence test statistics for categorical variables). For significantly correlated (p<0.05) feature groups, feature with lowest univariable Cox regression p value in the corresponding treatment group was retained as a representative of the group, and others excluded to avoid redundant information. Of the remaining features and covariates, the one with lowest p-value ratio in the two treatment groups (low divided by high) was used in model training.

Final model construction

The two final sets of features and covariates predicted to be most informative of the differential response to Dox or Dox+Evo were used to build the corresponding separate multivariable Cox proportional hazards regression models. Risk scores that are log-transformed relative risks of death were calculated using the predict.coxph function in R for all patients in the model training cohort and used to determine threshold for patient virtual inclusion and exclusion from the trial. The process of determining the optimum risk score threshold is described in the results section.

Risk score values were predicted for all patients in the test cohort from the final multivariable Cox models constructed as above. The threshold values found to result in best separation of the treatment arms as found in the training set were applied to enrich the test cohort, and survival was compared between the treatment arms in the included subset of the test cohort using log-rank test.

For the Dox model, where patients with high risk score values were expected to perform poorly under Dox, and thus be more likely to favor Dox+Evo treatment, a search for the optimal threshold was performed iteratively including sub-cohorts of patients with risk score above 1^st, 2^nd,3^rd etc to 97^th percentile of the total training cohort, evaluating survival difference between the treatment arms in terms of Cox regression p-value and hazard ratio each time. Thus, the entire range of possible thresholds was interrogated, to check if such selection can lead to significant treatment arm separation and identify an optimal threshold value.

Patients

Clinical covariates included in the analysis are listed in Table 1. Description of the clinical features is included in Supplementary Table 2. Presence of lung metastasis was associated with significantly poorer overall survival in the entire cohort of 588 patients (p=0.048, HR=1.24 (1.00-1.54)). This was not the case for patients with liver metastases, which was the second most common metastatic site (p=0.41). Among patients with lung metastases, no significant survival difference was observed between treatment groups (p=0.8), similarly to the entire cohort (p=0.45).

No significant difference in OS was observed between the full training and test cohorts (p=0.40, median OS: 17.4 (15.2-20.6) vs. 20.4 (14.0-26.9) training vs test). No significant differences were also seen between Dox and Dox+Evo treatment groups in the training (p=0.67, HR=1.07 (0.78-1.49) median OS: 18.3 (12.6-21.2) vs. 17.2 (15.2-22.1) months Dox vs. Dox+Evo) or test cohorts (p=0.30, HR=1.31 (0.46-1.23) median OS: 23.3 (16.5-31.8) vs. 14.9 (11.1-27.2) months Dox+Evo vs. Dox).

Feature selection

Results of univariable Cox proportional hazards regression, performed separately in the Dox and Dox+Evo treatment groups to identify those radiomic features and clinical covariates differently associated with OS, are shown in Table 2. Among clinical covariates, the histological classification of the primary tumor, tumor grade and prior radiotherapy were significantly associated with survival in the Dox and not in the Dox+Evo groups. Following elimination of correlated features, histology and prior radiotherapy remained, and histology was chosen to be included in the model due to lower p value ratio (0.015 vs. 0.029). There were no clinical features significantly related to survival in the Dox+Evo group. A history of smoking was significantly associated with shortened survival (p=0.04, HR=1.62 (1.01-2.61)) in the Dox group yet was insignificantly associated with longer survival in the Dox+Evo group (p=0.17, HR=0.73 (0.46-1.15).

No features were found to be significant in the Dox+Evo and not in the Dox group. Seven uncorrelated radiomic features were found significantly associated with survival in the Dox and not in the Dox+W.D.T was involved in collectioEnvo group: Volume Density, Number of Connected 3D Components, Inverse Variance, Correlation, Short Run Emphasis, and Small Zone Emphasis. Of these features, Short Run Emphasis, a measure of heterogeneity, showed the lowest ratio of p value in the Dox to the p value in Dox+Evo groups, and was chosen for training a prediction model of post-treatment survival.

Multivariable model

The three features identified above (histology, non-smoking history, and radiomic Short Run Emphasis) were combined in a multivariable Cox model trained on the Dox cohort of the training data set, producing a highly significant signature (p=0.0001) of survival. Details of the model are shown in Supplementary Table 3. No corresponding model was developed in the Dox+Evo group, as no clinical or radiomic features specific to this treatment arm were identified. The Dox model was used to predict risk scores for the entire training set cohort, including the Dox+Evo group, providing a predicted measure of risk of death if the Dox treatment was applied to all patients. Because patients with highest risk scores for Dox monotherapy are expected to benefit the most from the alternative (Dox+Evo) treatment, they should be included in the trial.

On the other hand, patients with a low risk score should be excluded and undergo Dox monotherapy instead. Such patient enrichment strategy for the trial should result in an improved treatment benefit of Dox+Evo in the included patients. To assess this, we performed the log-rank test for difference in survival between Dox vs. Dox+Evo as a function of the risk score threshold for the remaining patients whose score was above that threshold. The schematic of the process is shown in Figure 1. As described in Methods, the threshold separating high-from low-risk groups was incrementally increased to identify an optimum that reached a significant difference in OS, while including the largest fraction of patients. The results of this analysis are shown in Figure 2, demonstrating that increasingly different OS can be observed for the two treatment groups when patients with low-risk scores are excluded from the analysis (Figure 2A).The smaller p-values encountered with increasing thresholds were consistent with decreasing HR <1 (Figure 2B), showing increasingly significant treatment benefit of Dox+Evo vs. Dox with more stringent inclusion criteria. A threshold of 1.00 allowed inclusion of 42% of the initial training cohort in the trial and showed a significant advantage of Dox+Evo over Dox [p=0.017, HR=0.57 (0.36-0.90)]. This result was visualized in divergent Kaplan Meier curves for the treatment groups in the included patients and longer survival in the Dox+Evo group (Figure 2C, median survival 15.8 (14.2-21.4) vs 9.1 (7.6-13.7) months Dox+Evo vs Dox) and a reverse trend for the excluded patients (Figure 2D, median survival 20.9 (14.7-27.2) vs 27.0 (20.1-Not estimable(N.E)) Dox+Evo vs Dox, p=0.036), with Dox treated group showing significantly longer survival, compared to no difference in the whole training cohort (Figure 2E). These figures make apparent that the most significant difference (p<10⁻⁵) between included and excluded groups is their response to Dox monotherapy. Indeed, the difference in survival between included and excluded groups in the Dox+Evo trial arms was insignificant (p=0.48). A model using only clinical features as input was also trained and its performance evaluated as above, showing significant (p=0.0007) advantage in survival of the Dox+Evo vs. Dox, at and above a risk score threshold of 1.23 corresponding to a lower inclusion rate of 32% of initial training cohort (Supplementary Figure 1). Although the Short Run Emphasis feature by itself did not significantly discriminate (minimum p-value 0.097 including only patients with SRE<0.84, Supplementary Figure 2), it added to the significance of clinical features and thus increased the number of potential patients on trial from 32% to 42%.

• Download figure
• Open in new tab

Figure 1. Patient inclusion model.

Patient selection into the trial based on Dox group survival was executed according to the following method: firstly (1) radiomic and clinical features associated in training cohort with survival in Dox but not Dox+Evo treatment group were included in a multivariable Cox regression model (2), trained on Dox treated patients. The risk score assigned by the model to each training set patient was then used as a biomarker for inclusion into analysis, iteratively calculating the p value and hazard ratio for survival comparison between treatment arms depending on minimum risk score threshold (3). If available, threshold corresponding to significant (p value<0.05) treatment benefit of Dox+Evo at highest percentage of patients included was chosen (4), and subsequently tested in the test cohort (5), with risk scores assigned by the multivariable Cox model developed in step (2). A corresponding model can also be developed based on Dox+Evo group survival, using a maximum risk score threshold.

• Download figure
• Open in new tab

Figure 2. Multivariable Cox model enables selection of patients who benefit from Evofosfamide+Doxorubicin in training cohort.

Quantification of the p value of overall survival difference in the training cohort between the Evofosfamide+Doxorubicin (Dox+Evo) vs. Doxorubicin alone (Dox) treatment arms depending on the minimum risk score for patient inclusion, as predicted by the model (A), shows a risk score threshold of 1.00 at which Doxorubicin +Evofosfamide (Dox+Evo) group shows significantly longer OS (p<0.05). Exclusion of patients with high risk scores leads to monotonic decrease in the hazard ratio (B), and the 1.00 risk score threshold corresponds to inclusion of 42% of patients in the trial (indicated by red dotted line). The Kaplan-Meier plots by treatment arms show significantly better OS in the included (C) and significantly worse OS in the (D) excluded patients for the Dox+Evo treatment compared to Dox only. In all training set patients (E) no difference between the arms was observed.

Model testing

The multivariable Cox model trained in the above section was used to predict risk scores for all patients in the test cohort. Similar to the training cohort, an increase in minimum risk score threshold for inclusion lead to a monotonic decrease in p value (Figure 3A) and HR (Figure 3B) for the overall survival comparison between treatment groups. Applying the threshold of 1.00 determined a priori in the training set as the optimum threshold for inclusion showed a significantly better survival in the Dox+Evo vs the Dox treated group [p=0.002, HR=0.29 (0.13-0.63) Figure 3C]. As in the training cohort, this was associated with an increased median survival of 22.8 (12.3-N.E) for Dox+Evo vs. 6.3 (3.1-13.7) for Dox. As shown in Figure 3D the differences in the two treatment arms for the remaining excluded patients was insignificant for both OS (p=0.72) and median survival 26.0 (18.4-N.E) vs. 27.2 (18.6-N.E), similarly to the starting whole test cohort (Figure 3E). The threshold corresponded to inclusion of 38% of the initial test cohort. As in the training set (Figure 2C,D), the selection by risk score threshold separated patients who did and did not respond to Dox (p<10⁻⁴), whereas it did not discriminate (p=0.54) responses of the Dox+Evo group. The plot of p value vs inclusion threshold (Figure 3A) shows a matching profile of improving treatment benefit of the Dox+Evo treatment (because of decreasing effectiveness of Dox) with increasing risk score, further supporting the model and the use of radiomics in patient selection.

• Download figure
• Open in new tab

Figure 3. Results in the test cohort confirm the validity of the model.

Risk scores predictions in the test cohort based on a multivariable Cox model trained on Dox treated training cohort patients can be used to identify patients who would benefit from Dox+Evo treatment. Graph in (A) shows that increasing the minimum risk score of patients included in the analysis leads to a stronger difference in survival between the treatment groups, as described by the p value of the comparison. For the risk score threshold of 1.00, a highly significant difference is observed (red point and dotted line), which corresponds to a decreased hazard ratio of the combination vs standard therapy (B). These differences are apparent from the Kaplan-Meier curve in the included patients (C) showing significantly longer survival in the Dox+Evo group, while the excluded patients (D), or all test set patients (E) show no difference in survival between treatment groups.

Radiomic model interpretation

Given the complexity of the question, directly interpreting the imaging information in the model may be challenging. The analysis of p value vs. feature value threshold (as shown in Figures 3 and 4 for risk scores) for the Short Run Emphasis (SRE) feature, shows that the treatment benefit of Dox+Evo is maximized if only patients with tumors of low SRE are included (minimum p= 0.097, HR=0.62 (0.35-1.09) at SRE<0.845, Supplementary Figure 2). Taking the entire dataset into account, the p value at this threshold reaches p=0.01 signifying the relevance of this radiomic feature in patient selection. The multivariable model developed above also favors low SRE values, as shown in Figure 4A for both training and test cohorts. Comparing representative tumors with extreme SRE values reveals visual differences. In line with the model, a patient censored after over 2.5 years showed low SRE in the lung lesion (Figure 4B) at baseline; visually associated with regularity and homogeneity of the mass. Conversely, another patient deceased on Dox therapy less than 5 months after enrollment presented a lung lesion with high SRE and highly heterogeneous appearance (Figure 4C). While these show extremes, the value of using a quantitative SRE threshold is to identify those patients whose scans may be more equivocal.

• Download figure
• Open in new tab

Figure 4. Differences in radiomic features can be apparent visually.

The model for selection of patients likely to benefit from Evofosfamide treatment favored low Short Run Emphasis (SRE) radiomic feature for proposed inclusion into the trial. As shown in the violin plot (A), significantly lower SRE is observed in the included vs. excluded patient groups both in training and test cohorts. Qualitatively, a representative tumor with low Short Run Emphasis (SRE, B) appears more regular and homogeneous in a contrast enhanced CT scan than a corresponding tumor with similar volume (15.0 vs 16.5ml respectively), and relatively high SRE (C), which shows higher intratumor heterogeneity. In the violin plot a solid line indicates median while dotted lines indicate 25^th and 75^th percentile. **** p<0.0001