Radiomic Features Associated With HPV Status on Pretreatment Computed Tomography in Oropharyngeal Squamous Cell Carcinoma Inform Clinical Prognosis

Purpose There is a lack of biomarkers for accurately prognosticating outcome in both human papillomavirus-related (HPV+) and tobacco- and alcohol-related (HPV−) oropharyngeal squamous cell carcinoma (OPSCC). The aims of this study were to i) develop and evaluate radiomic features within (intratumoral) and around tumor (peritumoral) on CT scans to predict HPV status; ii) investigate the prognostic value of the radiomic features for both HPV− and HPV+ patients, including within individual AJCC eighth edition-defined stage groups; and iii) develop and evaluate a clinicopathologic imaging nomogram involving radiomic, clinical, and pathologic factors for disease-free survival (DFS) prediction for HPV+ patients. Experimental Design This retrospective study included 582 OPSCC patients, of which 462 were obtained from The Cancer Imaging Archive (TCIA) with available tumor segmentation and 120 were from Cleveland Clinic Foundation (CCF, denoted as SCCF) with HPV+ OPSCC. We subdivided the TCIA cohort into training (ST, 180 patients) and validation (SV, 282 patients) based on an approximately 3:5 ratio for HPV status prediction. The top 15 radiomic features that were associated with HPV status were selected by the minimum redundancy–maximum relevance (MRMR) using ST and evaluated on SV. Using 3 of these 15 top HPV status-associated features, we created radiomic risk scores for both HPV+ (RRSHPV+) and HPV− patients (RRSHPV−) through a Cox regression model to predict DFS. RRSHPV+ was further externally validated on SCCF. Nomograms for the HPV+ population (Mp+RRS) were constructed. Both RRSHPV+ and Mp+RRS were used to prognosticate DFS for the AJCC eighth edition-defined stage I, stage II, and stage III patients separately. Results RRSHPV+ was prognostic for DFS for i) the whole HPV+ population [hazard ratio (HR) = 1.97, 95% confidence interval (CI): 1.35–2.88, p < 0.001], ii) the AJCC eighth stage I population (HR = 1.99, 95% CI: 1.04–3.83, p = 0.039), and iii) the AJCC eighth stage II population (HR = 3.61, 95% CI: 1.71–7.62, p < 0.001). HPV+ nomogram Mp+RRS (C-index, 0.59; 95% CI: 0.54–0.65) was also prognostic of DFS (HR = 1.86, 95% CI: 1.27–2.71, p = 0.001). Conclusion CT-based radiomic signatures are associated with both HPV status and DFS in OPSCC patients. With additional validation, the radiomic signature and its corresponding nomogram could potentially be used for identifying HPV+ OPSCC patients who might be candidates for therapy deintensification.

Purpose: There is a lack of biomarkers for accurately prognosticating outcome in both human papillomavirus-related (HPV+) and tobacco-and alcohol-related (HPV−) oropharyngeal squamous cell carcinoma (OPSCC). The aims of this study were to i) develop and evaluate radiomic features within (intratumoral) and around tumor (peritumoral) on CT scans to predict HPV status; ii) investigate the prognostic value of the radiomic features for both HPV− and HPV+ patients, including within individual AJCC eighth edition-defined stage groups; and iii) develop and evaluate a clinicopathologic imaging nomogram involving radiomic, clinical, and pathologic factors for disease-free survival (DFS) prediction for HPV+ patients.
Experimental Design: This retrospective study included 582 OPSCC patients, of which 462 were obtained from The Cancer Imaging Archive (TCIA) with available tumor segmentation and 120 were from Cleveland Clinic Foundation (CCF, denoted as S CCF ) with HPV+ OPSCC. We subdivided the TCIA cohort into training (S T , 180 patients) and validation (S V , 282 patients) based on an approximately 3:5 ratio for HPV status prediction. The top 15 radiomic features that were associated with HPV status were selected by the minimum redundancy-maximum relevance (MRMR) using S T and evaluated on S V . Using 3 of these 15 top HPV status-associated features, we created radiomic risk scores for both HPV+ (RRS HPV+ ) and HPV− patients (RRS HPV− ) through a Cox regression model to predict DFS. RRS HPV+ was further externally validated on S CCF . Nomograms for the HPV+ population (M p+RRS ) were constructed. Both RRS HPV+ and M p+RRS were used to prognosticate DFS for the AJCC eighth edition-defined stage I, stage II, and stage III patients separately.

INTRODUCTION
The rise in the incidence of human papillomavirus (HPV)related cancers has caused a significant epidemiological shift (1) in oropharyngeal squamous cell carcinoma (OPSCC). It is estimated that HPV causes more than 70% of OPSCC cases in the United States (2). HPV+ OPSCC differs from its HPV− counterpart in response to treatment and disease aggressiveness (3). In order to account for this, the most recent American Joint Committee on Cancer (AJCC) eighth edition tumor staging system was modified to incorporate HPV status, with different staging systems for HPV (p16) positive and negative tumors (4). HPV+ patients tend to respond better to definitive radiotherapy or combined chemoradiotherapy protocols and are less likely than HPV− patients to develop disease recurrence and metastases. Thus, it has become critical to develop biomarkers within the HPV+ and HPV− populations for risk stratification.
Treatment of OPSCC patients is on the cusp of a paradigm shift; current clinical trials are geared toward reducing treatment toxicity for HPV+ patients without compromising survival outcomes since low-risk patients typically could benefit from a lower dose of radiotherapy or less invasive surgical operations (treatment deintensification) (5). However, it is challenging to tailor the most optimal treatment strategy for each patient. Although high T (T3/ T4) and N (N3) stages as well as tobacco use are clinically accepted risk factors for HPV+ OPSCC patients, these categorical predictors neglect the oncogenic differences between individual patients (6). A recent phase II randomized controlled trial NRG HN002 reported that patients who are grouped as candidates for treatment deescalation using a combination of clinicopathologic factors did not meet the goal of 2-year disease-free survival (DFS) in the radiotherapy only arm (7). This exhibits a clear unmet clinical need for the development of objective biomarkers to identify patients who could truly benefit from treatment de-escalation.
On the other hand, the unmet need for HPV− patients is to precisely identify patients at high risk of developing local or regional failure after treatment, patients who might be candidates for targeted treatment escalation (8). While the eighth edition of the AJCC staging modifications that separate HPV+ and HPV− patients is a major advance (9), it still may not be sufficient for accurate risk stratification as it de-emphasizes the importance of nodal metastasis by diagnosing most of the new p16+ OPSCC patients into stages I and II (10).
Radiomics is the process of computational extraction of large numbers of quantitative imaging features, such as texture features, from routine radiologic scans (e.g., MRI and CT) for characterization of the disease (11). These features are able to detect subtle changes in imaging intensity patterns within a local region which in turn may help better describe the cancer phenotype as well as the tumor microenvironment. While radiomic features in the immediate vicinity outside the tumor have shown significant value in differentiating disease subtypes for lung (12) and breast (13) cancers, we are not aware of any work that has attempted to collectively evaluate the role of textural patterns from both within (intratumoral) and outside the tumor (peritumoral) to predict HPV status or to identify their association with disease-specific survival in OPSCC.
In this study, we sought to explore the prognostic value of both intratumoral and peritumoral HPV status-associated radiomic features on CT scans and compared and combined them with clinical and pathologic factors on over 500 OPSCC patients. The prognostic radiomic biomarker for HPV+ OPSCC was validated both internally (237 patients) and externally (120 patients) on two different cohorts. We also evaluated the utility of the radiomic signature to prognosticate DFS within each individual AJCC eighth edition-defined stage group. Finally, this study also involved creation and validation of a clinicopathologic nomogram for estimating DFS for HPV+ OPSCC patients. Figure 1 shows the overall methodology comprising radiomic feature extraction and selection, prognostic signatures, and radiomic nomogram construction and validation.

Patients
Two OPSCC cohorts were included in this study: The Cancer Imaging Archive (TCIA, n = 462) OPC-Radiomics cohort (10)  and the Cleveland Clinic Foundation (CCF) cohort (n = 120). All patients had undergone pretreatment radiotherapy planning CT.
TCIA is an open archive of deidentified cancer-specific medical images and associated clinical metadata accessible for public download (14). Clinicopathologic and outcome information for patients in the CCF cohort were collected after obtaining approval from the Institutional Review Board of Cleveland Clinic. Demographic data are shown in Table 1. For the TCIA cohort, a total of 473 patients with OPSCC treated with curative intent at the Princess Margaret Cancer Center between 2005 and 2010 were reviewed. Histopathologic confirmation was used for the diagnosis of OPSCC and p16 immunohistochemistry was used to assess HPV status. Patients were triaged using inclusion criteria that involved the availability of i) radiotherapy planning CT scans with matched clinical information (HPV status by p16 immunohistochemistry, survival information) and ii) binary mask for gross tumor volume (GTV). Fourteen patients with the following criteria were excluded: i) CT images containing artifact (n = 6); ii) number of voxels within tumor is less than 200, which was deemed to be insufficient for feature extraction (n = 5); and iii) tumor mask contains normal brain tissue (n = 3). Following the patient exclusion criteria, 462 patients from the TCIA cohort and 120 patients from the CCF cohort were included for subsequent radiomic analysis. The flowchart for patient enrollment is illustrated in Supplementary Figure S1.

CT Imaging
The CT images for the TCIA cohort were acquired (10) from one of the following CT scanners: General Electric Discovery ST, General Electric Lightspeed Plus, or Toshiba Medical Systems Aquillion ONE. CT scans were acquired in helical mode with a slice thickness of 2.5 mm (General Electric) or 2 mm (Toshiba), at 120 kVp and 300 mAs tube current. Image resolution was 1 mm for all the scans. The CT images for the CCF cohort were acquired from either The General Electric Medical System or The Siemens Medical System. CT scans were acquired in helical mode with a slice thickness of 3 mm, at 120 kVp and 235 mAs tube current. Image resolution is between 0.4 and 0.5 mm for most of the patients, with an image matrix of 512 × 512.

Intratumoral and Peritumoral Compartment Definitions
The binary intratumoral masks which outlined the primary GTV were obtained using the Radiation Therapy Structures  (RTSTRUCT) for the TCIA cohort (10). Primary tumors on the CCF cohort were manually segmented by two board-certified head and neck radiologists JL (with 5 years of clinical expertise) and SS (with 6 years of clinical expertise) across all of the twodimensional CT sections using a hand-annotation tool in axial view. Morphologic dilation operations were then performed for all patients on intratumoral masks to define the annular ring region outside the tumor up to a radial distance of 15 mm based on previous studies in lung (15) and breast cancer (13), where peritumoral margins >15 mm were not associated with disease recurrence. The intratumoral masks were then subtracted from the dilated masks to obtain the peritumoral regions, which were then subdivided into three peritumoral rings of 5-mm-radius increments. Implementation details on peritumoral masks are provided in Appendix E1, section 1.

Statistical Analysis
Within the TCIA cohort, 462 patients were randomly allocated to S T (180 patients: 100 HPV+ and 80 HPV−) and S V (282 patients: 237 HPV+ and 45 HPV−) in an approximately 3:5 ratio. For both HPV+ and HPV− patients, the ratio of non-censored patients in S T and S V was kept balanced. The clinical end points of interest for this study were HPV status and DFS. DFS was defined as the time interval from the radiotherapy end date to the date of either last follow-up (censored) or local, regional, distant failure and death (event), whichever happened first. The difference of continuous variables (i.e., age) between cohorts (i.e., S T , S V , and S CCF ) was determined using the analysis of variance (ANOVA) and the association between categorical factors was estimated using the chi-square test.
Since there were no HPV− patients in the CCF cohort, we used only the TCIA cohort for HPV status prediction. A machine learning classifier was first constructed using a combination of intratumoral and peritumoral features for the prediction of HPV status using S T . To remove redundant features, all possible pairs of features in S T were tested for correlation by calculating the Spearman correlation coefficient (SCC). For any pair of features with SCC greater than 0.80, the feature with the higher Wilcoxon rank sum p-value was removed. A linear discriminant analysis (LDA) machine-learning classifier was subsequently trained in conjunction with the minimum redundancy-maximum relevance (MRMR) (20) feature selection approach using a 100-run, 3-fold cross-validation setting. The top 15 most frequently selected radiomic features (F t ) that best discriminated between HPV+ vs. HPV− across all iterations were identified from S T . An unadjusted p-value <0.05 (using two-sided Wilcoxon rank sum tests) was employed to indicate statistical significance. These features were then evaluated via the LDA classifier in terms of HPV status prediction on S V using the area under the receiver operating characteristic curve (AUC) metric.
We also developed dedicated radiomic risk score classifiers for patients within the individual HPV+ (RRS HPV+ ) and HPV− (RRS HPV− ) categories. The least absolute shrinkage and selection operator (LASSO) method was applied to F t using the cases in S T for both HPV+ and HPV− patients. After identifying the top ranked features, the corresponding LASSO coefficients were used for constructing risk classifiers for the HPV+ (RRS HPV+ ) and HPV− (RRS HPV− ), respectively. Both RRS HPV+ and RRS HPV− were calculated for each patient via a linear combination of selected features that were weighted by corresponding coefficients: where n (ranging from 0 to 15) is the number of features selected by LASSO for HPV+ or HPV− patients, x i refers to the HPV status-associated feature value, and b i is the corresponding weighted coefficient. The potential association of RRS HPV+ and RRS HPV− with DFS was first assessed in S T and then evaluated in S V . The prognostic ability of RRS HPV+ was further externally validated in S CCF . Patients were classified into high or low risk based on the median value of the RRS HPV+ or RRS HPV− in S T , which was then applied to S V and S CCF . Kaplan-Meier survival curves were used to visualize the survival rate for the high-and low-risk groups. At any given point on the survival curve, the probability that a patient in either the high-risk or low-risk group remains alive is presented (21). The log-rank test and hazard ratio were used to compare the survival differences between the two groups. The same Kaplan-Meier survival analyses were further performed for each cancer stage group defined by the AJCC eighth edition on the dataset combining S V and S CCF . The final values for n were determined when the hazard ratio for high risk over low risk reached the highest in S V for HPV+ patients (Supplementary Figure S4). The value of the tuning parameter in the LASSO-Cox model (l) was averaged out by 10 crossvalidation to minimize the error within S T . Constructions of RRS HPV+ and RRS HPV− were performed using in-house software implemented in the MATLAB R2019b platform (MathWorks).
Univariate Cox proportional hazards analysis on the effect of RRS HPV+ , RRS HPV− , and the individual clinicopathologic variables (gender, smoking status, drinking status, T stage, N stage, and AJCC eighth edition of overall stage) on DFS was evaluated. Variables significant in univariate analysis were included for multivariable Cox proportional hazards analysis to investigate the relationships between the various covariates (including the RRS HPV+ , RRS HPV− , and the clinicopathologic variables).
To further investigate the independent prognostic value of the RRS HPV+ with existing clinical factors (gender, smoking, and drinking status) and pathological staging factors (T stage, N stage, and the AJCC stage eighth edition), we constructed nomogram models for HPV+ patients comprising i) only the clinical factors (M c ) of gender, smoking, and drinking history; ii) only the pathologic staging factors (M p ) of AJCC eighth edition overall stage and N stage; and iii) both the pathologic T stage and the RRS HPV+ (M p+RRS ) using S T and validated them in S V and S CCF . The prognostic ability of M p+RRS was compared against M c and M p for HPV+ populations in terms of the concordance index (C-index). The nomograms were validated using 1,000 bootstrap resampling to calculate C-index with confidence intervals. Calibration curve analysis was performed to compare the nomogram-predicted DFS with the actual DFS. Decision curve analysis was adopted to calculate the net benefit for M p+RRS in comparison with M c and M p , for verification of the clinical usefulness of the nomogram (22). Nomogram construction, calibration plot, and decision curve were implemented using the "rms" and "SvyNom" packages under R statistical software (version 4.0.3; R Foundation for Statistical Computing, Vienna, Austria).

Clinicopathologic Characteristics
The clinical and pathologic characteristics of patients in the TCIA training set (S T ), TCIA internal validation set (S V ), and CCF external validation set (S CCF ) are summarized in Table 1. No significant differences were found in most features among S T , S V , and S CCF . All patients included in this study underwent radiotherapy. Eighty-seven of the 180 patients in S T (48.3%) and 145 of the 282 patients in S V (51.4%) were also treated with chemotherapy. The median DFS for HPV+ patients was 6.19 years in S T and 6.24 years in S V . The median DFS for HPV− patients was 1.27 years in S T and 1.55 years in S V . Thirty-three

Experiment 1: Prediction of HPV Status
Within S T , 454 uncorrelated intratumoral and peritumoral features were obtained after feature pruning. From these uncorrelated features, the top 15 were identified for predicting HPV status by MRMR feature selection, of which 11 were peritumoral features and 4 were intratumoral features ( Table 2). Changes in classification performance on account of differently selected features in S T are provided in Supplementary Figure S2. Standard deviation of CoLlAGe sum of variance from the 0-5-mm rings ( Figure 2C) was identified as the most discriminating peritumoral feature (p < 0.001), while the median of Haralick correlation-info1 ( Figure 2B) was identified as the most discriminating feature (p < 0.001) within the tumor. Both features were differentially expressed on CT scans between HPV+ and HPV− patients, and representative examples are illustrated via colormaps of the feature representations overlaid on the tumor areas or the annular ring regions around the tumor (Figure 2A). Noticeably, three of the four selected Haralick features were measures of correlation applied with different statistics on both intratumoral and peritumoral regions, suggesting that there might be significant pixel correlation-related pattern differences between the HPV+ and HPV− patients. We also observed consistently higher Laws feature values for HPV− patients compared with HPV+ across S T and S V . These features quantify intensity smoothness, abrupt edge changes, and ripple patterns on CT scans. Interestingly, all the selected Laws and  Figure 2D illustrates that using combined peritumoral and intratumoral features improved AUC compared with only using the intratumoral features and using only the peritumoral features for predicting HPV status.

Experiment 2: Prognosticate DFS for Both HPV− and HPV+ Patients, Including Within Individual AJCC Eighth Edition-Defined Stage Groups
A three-feature radiomic signature was identified as having the best prediction of DFS for both HPV+ and HPV− patients (Supplementary Figure S4). The radiomic risk scores for HPV+ patients (RRS HPV+ ) and HPV− patients (RRS HPV− ) across S T , S V , and S CCF are illustrated in Supplementary Figure S5. Details of the features selected are provided in Table 2, with their coefficients in the Cox model provided in Supplementary  Figure S6. RRS HPV+ was constructed using two peritumoral and one intratumoral features, while RRS HPV− was constructed using three peritumoral features. The median RRS HPV+ value (−0.0809) in S T was used as the cutoff threshold for defining high-and low-risk groups, resulting in statistically significant The prognostic ability of the radiomic risk score (RRS HPV+ ) was further evaluated for the patients within the AJCC eighth edition-defined different stage groups in S V +S CCF . The HRs of predicting DFS using RRS HPV+ for stage I ( Figure 4A Results of the univariable analysis are shown in Table 3. T3 stage, a moderate or heavy drinking history, and a higher RRS HPV+ were significantly associated with worse DFS for the HPV+ population in S T . N2 stage and a higher RRS HPV− were significantly associated with worse DFS for the HPV− population in S T . In multivariable analysis, RRS HPV+ (DFS hazard ratio,   Table 5). In the decision curve analysis, M p+RRS yielded a better net benefit compared with M c or M p individually when the threshold probability <0.35 ( Figure 5D).
When evaluating the prognostic ability of M p+RRS for HPV+ patients in S V +S CCF within each AJCC eighth edition-defined stage group, univariable Cox proportional hazard regression yielded HRs of 1. 25

DISCUSSION
HPV+ OPSCC has better clinical prognosis and treatment response than the alcohol-and tobacco-related HPV− OPSCC (23). Because of this, treatment deintensification to reduce therapy-related morbidity in "low-risk" HPV+ OPSCC patients is being considered. However, a recent phase II randomized controlled trial by Yom et al. reported that patients in whom therapy was deintensified based on traditional TNM staging information did not meet the goal of 2-year DFS >85% (7). This was likely due to the lack of biomarkers for identifying patients who would most benefit from therapy deintensification. Although the AJCC eighth edition staging system represents a  major improvement over the seventh edition, better and more reliable methods for pretreatment prognostication are needed for therapeutic decision-making. The radiomic biomarkers presented in this study aimed to identify those truly low-risk HPV+ patients within both the whole and AJCC eighth edition stage 1 population groups, for whom treatment deintensification should be considered.
In this work, we investigated the ability of both intratumoral and peritumoral radiomic biomarkers on CT scans to predict HPV status for a large cohort of 462 OPSCC patients. Additionally, we addressed the OPSCC prognosis prediction problem independently within the HPV+ and HPV− populations. The radiomic signature which was prognostic for the whole HPV+ population was also prognostic within the AJCC eighth edition stage I and stage II patients. Finally, the radiomic features were combined with pathologic staging factors to form a radiomic nomogram for individualized prognosis estimation for HPV+ OPSCC patients.
Currently, measuring p16 protein expression via immunohistochemistry is the recommended test for determining HPV status. However, distinct populations of patients exist in whom the tumors overexpress p16 but are in fact negative for HPV-DNA or mRNA expression and vice versa (24). Consequently, using p16 testing alone to determine HPV status results in some misclassified patients. Additional biomarkers are needed to complement the p16 testing. In the first experiment, we demonstrated that radiomic features from within the tumor and annular rings of 0-15 mm outside the tumor on CT imaging could reasonably predict HPV status of OPSCC, with an overall accuracy of 76%. A recent study by Leijenaar et al. designed a statistical framework for HPV status prediction, and they found that HPV+ tumors are more homogeneous in CT densities (25). This finding is in alignment with the result of this study. Specifically, HPV+ tumors possess a more homogeneous morphologic appearance in terms of CT texture patterns compared with HPV− tumors, which in turn is characterized by the Haralick correlation of information measured within the tumor (Figure 2A, second column). A higher value of the Haralick correlation indicates less pixel intensity disorders and decreased morphologic appearance heterogeneity for HPV+ tumors compared with HPV− tumors. Similar to the findings we report, Bagher-Ebadian et al. reported that HPV+ tumors have consistently lower energy components for seven frequency bands quantified by the DOST features (26). Although the discriminating textural features we identified are  different from previous studies, the interpretations of the features are similar. However, we also found higher gray-level intensity values of HPV+ tumors compared with HPV− tumors on CT, which appears to be at odds to the findings of the study by Leijenaar et al. (25). In addition, our study represents the first study for integrating CT radiomic both within and around the tumor for OPSCC HPV status prediction. We demonstrate the superior discriminability of peritumoral CT radiomic features, which appears to suggest discriminable differences of the microenvironment in the regions immediately outside the tumor. HPV− patients are best characterized by a combination of local intensity disorder and microscale heterogeneity in gradient orientation, particularly outside the tumor. Specifically, a higher peritumoral variation of the gradient orientation defined by the CoLlAGe sum of variance was observed in HPV− compared with HPV+ tumors. Reduced expression of Haralick information  capturing low level of correlation between adjacent pixels was also an important component of the peritumoral radiomic signature for HPV− patients. Furthermore, HPV− peritumoral regions were characterized by abrupt changes in edges, ripples, and intensity smoothness, as detected by the elevated expression of Laws features. These findings are consistent with previous studies, which showed that HPV+ tumors have less stroma overall, have smoother borders to the nests and leading edges, and are more homogeneously cellular usually without keratin production (27)(28)(29).
In the second experiment, we evaluated the prognostic ability of the HPV status-associated radiomics features found in experiment    (7). Their results showed that patients randomly assigned to the non-chemotherapy arm did not meet the goal of 2-year DFS >85%. The authors also reported that for these patients using a low-than-standard-dose radiotherapy resulted in a higher rate of locoregional failure. This indicates that risk stratification based on grouping of clinicopathologic factors alone is not robust and there is a clear unmet clinical need for developing more granular and more objective biomarkers to identify patients who could truly benefit from treatment deescalation. Clinically, stage I and stage II patients are the current target for treatment de-escalation. However, a subset of these patients still had poor survival outcome and would not benefit from treatment de-escalation (33). The CT radiomic risk score developed in this study represents a potential useful tool for guiding treatment intensities within the early-stage HPV+ OPSCC patients. By applying the threshold defined from the training set, the learned radiomic risk scores could further separate the stage I and stage II HPV+ patients into high-vs. low-risk groups based on disease-free survival. As such, the radiomic biomarker presented in this study could potentially help to distinguish patients within the current AJCC eighth edition definition of low risk as to which patients will benefit from treatment deintensification vs. those who will not. The findings of this study are consistent with a recently published study, where a histology-based imaging biomarker (MuNi) was found to be associated with survival for stage I and stage II patients (34).
Noticeably, the prognostic radiomic features identified in this study are mainly from peritumoral compartments. One possible explanation for this is that the peritumoral radiomic features are associated with tumor-infiltrating HPV-specific immune responses prior to treatment, which have been more commonly found floating around the tumor without actual penetration or action into the tumor core and are strongly associated with prognosis (35,36).
In the third experiment, we developed a radiomic nomogram M p+RRS for HPV+ DFS prediction. Currently, the conventional TNM and AJCC staging systems are routinely used for risk stratification and prognosis estimation. They reflect tumor size (T), lymph node status (N), and cancer metastasis (M). However, these staging factors could not capture the intratumor heterogeneity, which has been shown to be a significant prognostic factor. The radiomic nomogram combined the pathologic staging information with the radiomic features extracted from the entire tumor on CT scans, enabling for robust pretreatment survival estimation. Combining the RRS HPV+ from experiment 2 with the pathologic T stage resulted in a nomogram that leads to a more individualized prognosis prediction. With refinement and improvement, this type of radiomic approach might guide more tailored treatment for patients with better survival outcome. Compared with risk stratification using only the conventional staging factors, the radiomic nomogram had an improved DFS estimation. A previous study by Fakhry et al. showed that a nomogram integrating clinicopathologic factors (i.e., HPV status, T and N stages) could reliably predict progression-free survival (37), which is in alignment with our results. Based on our results, T3 stage is significantly associated with worse DFS for HPV+ patients, while N1 and N2 stages are significantly associated with worse DFS for HPV− patients in both univariable and multivariable analyses. With regard to the radiomic nomogram for head and neck cancer, Zhang et al. built multiparametric MRI-based radiomic nomograms for predicting nasopharyngeal carcinoma prognosis and obtained C-index of 0.776 for PFS prediction (38). Yuan et al. proved that a nomogram consisting of MRI radiomic signatures and TNM stage could better predict head and neck cancer prognosis with a C-index of 0.72 on the validation set (39). The prognostic performance difference between our new CT-based nomogram and the MRI-based nomogram may be on the higher image resolution offered by multiparametric MRI, although CT tends to be used more routinely compared with MRI for head and neck cancer. We also noted that both patient cohorts from Zhang et al. (38) and Yuan et al. (39) comprised a majority of advanced stage head and neck cancer patients (100% and 70.6%), while our training and validation cohorts consist of only 17% and 26.6% of stage III and no stage IV OPSCC patients. This may also influence the performance of the model since intratumor heterogeneity could be more easily captured within those aggressive tumors.
Our study does have several limitations. First, the prognostic biomarker validation on a single cohort was done in a retrospective manner. Second, we predicted neither benefit of existing treatments for the two populations nor treatment response within the individual AJCC eighth edition-defined stage groups. These aims will be part of our future study involving large multisite and multimodality evaluation of radiomic signatures in predicting treatment response for the two populations. Third, we acknowledge the limitation of our dataset from the TCIA in terms of HPV status based on p16 testing which may not accurately reflect the true transcriptionally active HPV status, at least for a small percentage of patients. Nonetheless, this study demonstrates that CT radiomic features could, in theory, complement the existing p16 testing method in distinguishing HPV status.
Despite the aforementioned limitations, this study is the first to show the role of combined intratumoral and peritumoral radiomic features in predicting HPV status of OPSCC patients. It is also the first study to incorporate both radiomic signatures and corresponding nomograms for prognosis prediction for HPV+ and HPV− patients. If confirmed in prospective clinical trials, this radiomic nomogram pipeline could enrich the existing AJCC eighth staging systems for risk-stratifying OPSCC patients. One can imagine a strategy where numerous sources of data go into predictive models for patient care. Especially attractive here is that all patients received pretreatment cross-sectional CT scans so the data are already garnered in digitized form, easily available for radiomics-based nomograms for prognosis prediction.

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 Cleveland Clinic. The patients/participants provided their written informed consent to participate in this study.

AUTHOR CONTRIBUTIONS
BS: conceptualization, data curation, formal analysis, methodology, investigation, writing-original draft, and writing -review and editing. KY: conceptualization, investigation, and writing-review and editing. JG: writing-review and editing. CL: writing-review and editing. LL: formal analysis and investigation. NB: writing-review and editing. CK: writingreview and editing. PT: resources and investigation. PF: conceptualization, investigation, and writing-review and editing. SK: resources, investigation, and writing-review and editing. JSL: resources, investigation, and writing-review and editing. AM: conceptualization, methodology, investigation, supervision, funding acquisition, project administration, and writing-review and editing. All authors contributed to the article and approved the submitted version.