Targeting Treatment Resistance in Head and Neck Squamous Cell Carcinoma – Proof of Concept for CT Radiomics-Based Identification of Resistant Sub-Volumes

Purpose Radiomics has already been proposed as a prognostic biomarker in head and neck cancer (HNSCC). However, its predictive power in radiotherapy has not yet been studied. Here, we investigated a local radiomics approach to distinguish between tumor sub-volumes with different levels of radiosensitivity as a possible target for radiation dose intensification. Materials and Methods Of 40 patients (n=28 training and n=12 validation) with biopsy confirmed locally recurrent HNSCC, pretreatment contrast-enhanced CT images were registered with follow-up PET/CT imaging allowing identification of controlled (GTVcontrol) vs non-controlled (GTVrec) tumor sub-volumes on pretreatment imaging. A bi-regional model was built using radiomic features extracted from pretreatment CT in the GTVrec and GTVcontrol to differentiate between those regions. Additionally, concept of local radiomics was implemented to perform detection task. The original tumor volume was divided into sub-volumes with no prior information on the location of recurrence. Radiomic features from those sub-volumes were then used to detect recurrent sub-volumes using multivariable logistic regression. Results Radiomic features extracted from non-controlled regions differed significantly from those in controlled regions (training AUC = 0.79 CI 95% 0.66 - 0.91 and validation AUC = 0.88 CI 95% 0.72 – 1.00). Local radiomics analysis allowed efficient detection of non-controlled sub-volumes both in the training AUC = 0.66 (CI 95% 0.56 – 0.75) and validation cohort 0.70 (CI 95% 0.53 – 0.86), however performance of this model was inferior to bi-regional model. Both models indicated that sub-volumes characterized by higher heterogeneity were linked to tumor recurrence. Conclusion Local radiomics is able to detect sub-volumes with decreased radiosensitivity, associated with location of tumor recurrence in HNSCC in the pre-treatment CT imaging. This proof of concept study, indicates that local CT radiomics can be used as predictive biomarker in radiotherapy and potential target for dose intensification.


INTRODUCTION
Head and neck squamous cell carcinoma (HNSCC) accounts for approximately 4% of all malignancies across Europe and the USA (1,2). For locally advanced HNSCC, the standard of care is definitive radiotherapy, whenever possible combined with concurrent chemotherapy. Despite advances in treatment using modern radiotherapy delivery techniques, local recurrences still occur in up to 50% of patients and pose the predominant pattern of failure in HNSCC (3). Therapeutic options for recurrent HNSCC are mainly palliative and comprise salvage surgery, re-irradiation and systemic therapyhowever, outcome is poor, and treatment often is associated with significant morbidity (4). Therefore, improvement of primary radiochemotherapy is essential.
Several early clinical trials have studied the role of radiation dose intensification to the primary tumor and metastatic lymph nodes in order to improve local control for HNSCC (5)(6)(7). However, concern still exists about excessive acute and late toxicities of this approach and therefore no large randomized trial has been conducted so far. In addition, intensification of radiotherapy by using altered fractionation schemes only lead to very modest improvement of outcome (8). Thus, a dose of 70 Gy delivered over 7 weeks to the entire tumor is widely considered a standard in patients undergoing chemoradation for HNSCC (9). Heterogeneity within a tumor has been recognized as a characteristic of head and neck cancer (10), with some sub-regions being more resistant to radiochemotherapy. Consequently, the strategy of treatment intensification to these subvolumes could lead to better outcomes in terms of local control and subsequently overall survival, without a significant increase in treatment-related toxicities. Recent advancement in radiation technology in principle allows such an escalation of radiation dose to tumor sub-volumes. However, identification of these subvolumes is a crucial step within this therapeutic concept.
In contrast to selective biopsy specimens obtained from a small area of the tumor, the use of imaging as biomarker has the advantage to analyze the entire three dimensional tumor volume. The feasibility of delivering a dose boost, so-called dose painting, to tumor sub-volumes has been previously demonstrated, mostly based on functional, metabolism-related (fluorodeoxyglucose) or hypoxia-related (fluoromisonidazole) positron emission tomography (PET) imaging (11)(12)(13). In recent years, highthroughput, multidimensional and quantitative images analysis (radiomics) revealed that relevant biological information can be extracted also from routinely acquired, easy to perform morphological imaging, such as regular computed tomography (CT) (14)(15)(16). Many investigations have studied the potential of radiomics to predict the risk for metastatic spread, progression-free survival, overall survival or biological phenotypes (17)(18)(19) showing encouraging results. In the context of radiotherapy, several studies indicated an association between CT heterogeneity and local tumor control (20)(21)(22)(23)(24). However, there is so far scarce data on radiomics as a method to individually tailor dose distribution to specific subvolumes within a tumor (25)(26)(27).
Here, we investigated CT radiomics for identification of radioresistant sub-volumes of the primary tumor leading to persistence or recurrence after curative radiochemotherapy. Reliable detection of these resistant sub-volumes may allow for a tailored treatment by increasing radiation dose to these parts.

Study Population
We retrospectively analyzed patients with primary locally advanced HNSCC (cT3/4 or cN+) treated at our institution between June 2004 and October 2015, who experienced a local in-field tumor recurrence. Local in-field tumor recurrence was defined as a recurrence occurring within the high-dose target volume (excluding lymph nodes), and had to be confirmed by biopsy. Only patients that received a definitive high dose radiation treatment with an equivalent total dose of 68 -70 Gy in 2-2.11 Gy fractions and a concomitant systemic therapy with either cisplatin and/or cetuximab during the course of radiotherapy were included. A further requirement was a contrast-enhanced CT (CE-CT) imaging before treatment and a FDG-PET/CT at the time of local recurrence. Contrast protocol was not standardized, however patients with contrast visible only in thyroid were not included. Follow up of patients was done according to the institutional routine. The first FDG-PET/CT for treatment response was usually done three months after completion of treatment. A subset of patients included in this study and now reviewed retrospectively has been treated within a prospective trial back then (28).
Data analysis was approved by the Swissethics and was carried out in accordance with Swissethics guidelines and regulations. Patients gave informed general consent.
We identified in total 66 patients fulfilling the above mentioned inclusion criteria. However, twenty-one cases had to be excluded, as a reliable image registration between follow-up 18 Ffluorodeoxyglucose ( 18 F-FDG) PET and CE pretreatment planning CT was not possible. Five additional patients were excluded due to the small size of recurrence (volume < 27 voxels), which precludes a reasonable radiomics analysis. Thus, in total 40 patients were available and were split into a training cohort, consisting of retrospectively collected data (n = 28), and a validation cohort, comprising patients from a prospective phase II clinical study (n = 12) (28). In the majority of cases, the primary tumor was located in the oropharynx, for further details see Table 1.

Contrast-Enhanced Planning CT Image Acquisition
Iodine contrast was injected intravenously prior to CT imaging. The contrast protocol varied over the years of data acquisition. Images were acquired on five different scanners (Siemens Somatom Definition AS n = 20, Siemens Volume Zoom n = 13, Siemens Somatom Plus n = 5, GE Healthcare, Discovery STE n = 1 and GE Healthcare, Discovery LS n = 1) with 120 or 140 kV tube voltage, 2 -3.27mm slice thickness and <1mm in-plane voxel size. The tube current varied between 120-450 mAs, and angular exposure adaption was applied in 7 scans (17%). Images were reconstructed using filtered back projection and soft kernel (B30).

Definition of Primary Tumor and Recurrence Volume
The primary tumor was contoured based on the CE planning CT assisted by the co-registered pretreatment diagnostic FDG-PET/ CT and MRI scans if available. However, in contrary to the FDG-PET, which was available in all the patients, a pretreatment MRI was only present in 5 cases. The recurrence region was contoured as high FDG uptake in the follow-up FDG-PET/CT at the time of the detection of the recurrence (Supplementary Figure S1). High FDG uptake region was considered as 40% SUVmax threshold-based sub-volume in the recurrence PET-CT. Details of the PET scanning protocol can be found in the supplement (Supplementary Table S1). The contours were then transferred to the initial planning CT by rigid registration with focus on the tumor region (Eclipse v.15, Varian Medical Systems, USA). Two structures were created for further analysis: the overlap of tumor recurrence and primary tumor (GTV rec ) and the primary tumor contour minus the recurrence (GTV control ). Prior to radiomics analysis both GTV rec and GTV control were postprocessed, by removal of contours on slices affected by metal artifacts.

Radiomics Analysis for Differentiation of Controlled to Non-Controlled Tumor Sub-Volumes and Detection of Radioresistant Volumes
Radiomics analysis was performed using Z-Rad software implementation (Python v 2.7). This implementation was benchmarked in the Image Biomarker Standardization Initiative (29,30). Images were resampled to cubic voxels (3.3 mm, largest voxel dimension in the dataset) using trilinear interpolation. Hounsfield unit range of -20 to 180 HU was set to limit the analysis to soft tissue. Due to the small size of analyzed volumes, the feature set was limited to intensity (n = 25) and texture features (n = 136). The full list of analyzed features is presented in the Supplementary Table S2. Volumes with less than 27 voxels were excluded from analysis.
In total, three different radiomics models were built ( Figure 1): The "bi-regional radiomics" analysis aims to differentiate tumor regions, which stay controlled after treatment (GTV control ) and the non-controlled region, which resulted in a recurrence detected on the FDG-PET CT (GTV rec ). For detection of radioresistant volumes within the primary tumor without prior knowledge on their location and prior to any treatment two "local radiomics models" were built. Here, the predefined region of interest (ROI), corresponding to the primary GTV and from now on named GTV, is divided into smaller sub-regions, which are later used as a mask for feature extraction. As a consequence, instead of one vector of features describing the GTV, local radiomics returns a set of radiomics-based parametric maps. These maps visualize for example changes in heterogeneity or contrast across the GTV. The local features are calculated based only on the voxels within the GTV, whereas the neighboring voxels (e.g. healthy tissue, manually excluded artifacts) are set to 'not a number' (discarded). We have implemented two approaches for the definition of sub-regions ( Figure 1): • Division of the GTV into a fixed number of sub-regions with the center of the grid attached to the center of the GTV (in this study 8 sub-regions, the number was chosen to provide sufficient information in 3D and to ensure large enough neighborhood to define texture). The same number of subregions was analyzed per patient and their volume depended on the GTV size. • Division of the GTV into sub-regions with fixed size. In the second approach, the placement of the grid is optimized to cover a GTV volume as large as possible. The full coverage is rarely feasible, because in order to avoid meaningless features calculated only on a few voxels, we have set a threshold T = 25 of minimum number of ROI voxels (voxels with value different than 'not a number') in the sub-region. In this study, we decided for 5x5x5 voxels grid with 5 voxel shift. There was no overlap between sub-regions to ensure independent description of the sub-regions. The 5 voxels grid size was a trade-off between ensuring a large enough neighborhood to define texture and a high spatial resolution of the model. Other sizes of the grid were not tested. In contrary to the first approach different numbers of subregions were analyzed per patient but their size was constant.

Statistical Analysis
The classification and detection tasks were performed using logistic regression. In the bi-regional model, the two labels were assigned to GTV control and GTV rec volumes. For local radiomics models, the recurrent label was assigned to regions with more than 50% of voxels overlapping with GTV rec . The 50% threshold was chosen arbitrarily, other thresholds were not tested. The remaining sub-regions were labeled as control.
The following feature reduction and binary classification procedure was used. Features with high and moderate correlation to number of voxels in the local region were excluded (r > 0.5). Principal component analysis combined with univariate logistic regression was used for dimensionality reduction (20). The number of retained principal components was computed using the Horn method. Radiomic features were grouped based on their Spearman correlation to principal components. Per group, only the feature with the largest area under operator receiver characteristic curve (AUC) was chosen, given that p-value < 0.05. The final model was built in the multivariable logistic regression with backward selection of variables based on the Akaike Information Criterion. The final model was tested in the separate set as described in the sections below. The confidence intervals were computed with 2000 stratified bootstrap replicates.
For the bi-regional radiomics model and the first local radiomics model based on 8 sub-volumes, training was performed on the training cohort and the results validated in the validation cohort as previously specified in the`study population`section. The 75th percentile threshold of predictions in the training cohort was set as classification threshold to optimize the sensitivity and specificity of the model.
The second local radiomics model with a fixed 5x5x5 grid size, was trained and tested in leave-one-out cross-validation (LOOCV) on the patient level. The full split into training and validation cohort was not feasible due to the presence of metal artifacts in some images. These artifacts, especially when present in the middle of the GTV, strongly influence optimal placement of the grid and thus for this experiment only patients with no visible artifacts in the GTV were selected. Data from 23 patients with no visible CT artifacts were analyzed. The final model was trained on 12 patients (= 48%), who had both classes of subvolumes (at least one sub-region with 50% contribution of recurrent voxels), which was a requirement of LOOCV.

Differentiation of Controlled to Non-Controlled Tumor Regions
On average, the volume of GTV rec was 20% (range: 2%-71%) of the initial primary tumor volume. A logistic regression model with backward selection of variables showed a good discrimination between GTV rec and GTV control , in both the training and the validation cohort, see Figure 2A  In the bi-regional radiomics, features were extracted from GTV rec and GTV control and only a differentiation between recurrent and controlled sub-volumes was performed; (B) in two local radiomics models, a detection task was performed and thus sub-volumes were defined without any prior information on the location of recurrence. In the classification task, a sub-volume was classified as recurrence (X) if more than 50% of the voxels belonged to the original contour of the recurrence (red). CI: 0.72 -1.00). This bi-regional radiomics model consisted of two features: "GLRLM gray level non-uniformity" and "merged GLCM sum entropy" -indicating that sub-volumes showing higher density heterogeneity have a higher propensity of recurrence. In the validation cohort, the model achieved high sensitivity (= 0.75) and specificity (= 0.83).

Detection of Radioresistant Sub-Volumes -Local Radiomics Model
The first model with division of the GTV into 8 equi-volume regions comprised three features: "GLCM cluster shade", "GLDZM gray level variance" and "histogram median". This local radiomics model showed slightly inferior performance (AUC = 0.70; 95%CI: 0.53-0.86) in the validation cohort compared to the bi-regional radiomics model ( Figure 2B). In the validation cohort, in all patients in whom recurrence subvolumes were detected, at least one of the sub-volumes was correctly identified. The median size of detected recurrence was 43% of the entire recurrence volume. This result is linked to the threshold for recurrence sub-volume definition (50%) and thus can be further improved in the future studies.
For the second definition with division of the GTV into subvolumes of an equal size (5x5x5 voxels) the rigid grid placement was chosen on an individual patient basis in order to maximize the coverage of studied tumor. 12 patients (= 48%) had both classes of sub-volumes (recurrence and tumor control, based on a 50% criterion). The average AUC in leave-one-out crossvalidation was 0.68. In three cases the AUC was below 0.5, indicating worse than random prediction. This was tracked back to a small number of analyzed sub-volumes (e.g. only two subvolumes) or a high recurrence involvement (more than 25%) in most of the sub-volumes, meaning that in this tumor no true controls with no contribution of radioresistant clones were present. Similar features were selected in different cross validation folds: GLSZM zone size entropy (n=12), GLCM joint maximum (n=9) and histogram range (n=8).
Detailed summary of the all models together with numbers of sub-volumes used for model training and validation is presented in Table 2. Number of sub-volumes in the local radiomics analysis per patient is summarized in the Supplementary Figures S2 and S3.

DISCUSSION
HNSCCs, exhibiting a high heterogeneity in CT images, were previously shown to respond worse to radiochemotherapy (20)(21)(22)(23)(24). Several prognostic radiomic signatures have been proposed recently, but studies on predictive signatures with potential impact on treatment optimization are scarce (31). In this study, we were able to differentiate intratumoral levels of radiosensitivity by means of CT radiomics. Further, we proposed an algorithm for pretreatment detection of radioresistant regions.
Identification of treatment resistant tumor sub-volumes by means of medical imaging has been previously investigated, mainly by detection on hypoxic imaging but also different imaging modalities. Tumor hypoxia is a known adverse prognostic factor for local control after radiotherapy of HNSCC (32). However, only two studies tried to correlate spatial location of tumor recurrence and initial hypoxic region, and only one of them showed overlap of those regions (33,34). A correlation of high FDG uptake parts of a tumor to regions of more resistant tumor subvolumes leading to recurrences has been shown (35,36). Accordingly, FDG-based dose escalation strategies have already been exploited and proven to be feasible (37). However, also other imaging modalities like diffusion-weighted MRI or dynamic contrast-enhanced MRI, were shown to be able to identify risk factors for worse outcomes in head and neck cancer (38). Comparison studies have shown that the volumes defined by However, in contrast to our study, their results were only tested in leave-one-out cross-validation, owing to a small number of patients. In addition, the study did not clearly define if CT information was used together with PET information for definition of the primary tumor contour. Delineation based solely on PET images may lead to neglecting tumor regions with low FDG uptake (41). Our local radiomics results indicate that sub-volumes with higher CT heterogeneity are more radioresistant, which is in agreement with previous studies showing that higher tumor-wide CT heterogeneity is linked with reduced local control rate (20)(21)(22)(23)(24). Local radiomics in PET/CT has also been investigated in the context of outcome prediction in nasopharyngeal cancer, showing a higher prognostic value than with models based on entire tumor volume analysis (26). However, this study did not use local radiomics for sub-regional detection of treatment resistance sub-volumes.
With the introduction of intensity modulated radiation therapy (IMRT) the concept of delivery higher doses to head and neck tumors has regained interest. IMRT allows not only for better sparing of OARs, but also enables delivery of simultaneously higher doses to selected areas (42), so called dose painting. The identification of radioresistant tumor sub-volumes that require higher radiation doses, is an unmet clinical need. The ESCALOX study currently investigates a dose escalation up to 80.5 Gy using a simultaneous integrated boost (SIB) to the whole primary tumor and large involved lymph nodes (43), whereas another group performed a planning study with dose escalation only to hypoxic areas within the tumor as defined by 18 F-Fluoroazomycin arabinosid (FAZA) PET/CT and found it to be feasible while respecting the maximum OAR constraints (44). Resistance of head and neck cancer cells to radiotherapy is not conditioned by one single biological feature but rather driven by several different mechanisms (39). Consequently, a method capturing all these distinct features is desirable. The underlying biology for the distinct radiomics signature is not simple nor exactly known for the individual case. But, a correlation of radiomic signature with underlying genomic alterations and biological phenotype has been shown (45). In our study, the selection of the radioresistant subvolumes is meant to be done by radiomics analysis, which does not account for any particular biological background/histopathological difference. However, a strong relationship between medical images, or more precisely, the extracted, quantitative imaging features, and the underlying tumor phenotype and biology was shown previously (46,47). Thus, in contrast to the above mentioned study, which selected the radioresistant parts solely based on hypoxic regions, our radiomics approach covers a broader spectrum of underlying biological phenotypes/alterations. Still, the histopathological factors behind the selected sub-volumes in our study are unknown, and thus definition of an adequate dose boost is not straightforward and should be a subject of further research.
This study was performed on relatively small patient cohort. However, for two out of three proposed models (bi-regional and local radiomics with fixed number of sub-volumes), a successful model validation in a priori defined cohort was performed. Remarkably, the performance of the models was higher in the validation cohort, in comparison with the training cohort. The validation cohort is a prospective cohort with standardized imaging and treatment protocolin comparison to the training cohort, which was retrospectively collected. Therefore, higher model performance in the validation might be a result of premature training or better data quality. In the head and neck region, HPV status is known to influence CT values distribution within primary tumor (48). Unfortunately, HPV status was unknown for a large proportion of patients (35%). Consequently, an analysis on HPV status as an effector is not possible due to the high number of missing values in the individual patient cohorts. The analyzed cohort is also heterogeneous in terms of tumor subsites, and theirs impact on the performance of the model should be further investigated. The presence of artifacts did not allow for a full Area under receiver operating characteristic curve (AUC) and 95% confidence intervals. The second local radiomics model was tested in the leave-one-out cross-validation, thus no results for the training cohort are shown and the validation AUC is the average over the folds. *Most frequently chosen features over different folds. # Median and range of the number of subvolumes analyzed over different folds.
validation of the fixed grid size local model. With the introduction of iterative metal artifact reduction reconstruction algorithms in CT few years ago, this limitation of the model is addressed for future studies (49). Availability of larger data collections may also permit for testing different settings in the local radiomics analysis, such as variable size of the grid or smaller value of the translation vector, allowing for an overlap between sub-volumes. Today, it is not clear which methodology provides optimal results. In our study, we assumed that all of the analyzed tumor sub-volumes are independent and thus, no overlap was allowed. Additionally, classification was performed using logistic regression. In the future, segmentation algorithms may be tested to improve the predictive power of local radiomics.
Additional validation of the proposed approach may be derived from surgical cohorts, where pretreatment local radiomics maps may be correlated with full-mount tumor histopathology. Alternatively, if contrast-enhanced CT is present at the time of recurrence similarity of the radiomics features in the corresponding areas in the two, sequential investigations could be measured to evaluate if the signature remains stable over time.
This study has an inherent selection bias, since only patients with observed recurrences were analyzed. In the real pretreatment classification, a standard radiomics analysis can be used prior to local radiomics to preselect patients with high risk of tumor recurrence, as shown by other studies (20)(21)(22)(23)(24).
In conclusion, this is the first study indicating that tumor radioresistance can be localized on pretreatment CT images with validation of the radiomics model in an independent cohort. This is a proof of concept study, indicating that local CT radiomics can be used as predictive biomarker in radiotherapy and potential target for dose intensification.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by Ethics Committee canton Zurich. The patients/ participants provided general consent to participate in this study.

AUTHOR CONTRIBUTIONS
MBperformed data analysis, supervised implementation of local radiomics, wrote part of the manuscript and reviewed the manuscript. MP -contributed to study design, performed image registration and contouring, wrote part of the manuscript and reviewed the manuscript. OR -provided expertise, contributed to study design, reviewed contours and registration, and reviewed the manuscript. TF, HS, EH-S, LR and LB -performed image registration and contouring, reviewed the manuscript. MSperformed standard radiomics analysis, reviewed the manuscript, AAimplemented local radiomics, reviewed the manuscript. MH, MG and ST-L -provided expertise, contributed to study design and reviewed the manuscript. All authors contributed to the article and approved the submitted version.

FUNDING
This project was supported by the Swiss National Science Foundation Sinergia grant (310030_173303). The clinical study used as one of the cohorts was supported by a research grant from Merck (Schweiz) AG. The funder bodies were not involved in the study design, collection, analysis, interpretation of data, the writing of this article or the decision to submit it for publication.