Apparent Diffusion Coefficient Histograms of Human Papillomavirus–Positive and Human Papillomavirus–Negative Head and Neck Squamous Cell Carcinoma: Assessment of Tumor Heterogeneity and Comparison with Histopathology

One hundred five consecutive patients with primary oropharyngeal and oral cavity head and neck squamous cell carcinoma underwent MR imaging with anatomic and DWI sequences. The collected tumor voxels from the contoured ROIs provided histograms from which position, dispersion, and form parameters were computed. Histogram data were correlated with histopathology, p16-immunohistochemistry, and polymerase chain reaction for human papillomavirus DNA. Diffusion phenotypes of human papillomavirus–positive and human papillomavirus–negative head and neck squamous cell carcinomas showed significant differences, which reflect their distinct degree of tumor heterogeneity. BACKGROUND AND PURPOSE: Head and neck squamous cell carcinoma associated with human papillomavirus infection represents a distinct tumor entity. We hypothesized that diffusion phenotypes based on the histogram analysis of ADC values reflect distinct degrees of tumor heterogeneity in human papillomavirus–positive and human papillomavirus–negative head and neck squamous cell carcinomas. MATERIALS AND METHODS: One hundred five consecutive patients (mean age, 64 years; range, 45–87 years) with primary oropharyngeal (n = 52) and oral cavity (n = 53) head and neck squamous cell carcinoma underwent MR imaging with anatomic and diffusion-weighted sequences (b = 0, b = 1000 s/mm2, monoexponential ADC calculation). The collected tumor voxels from the contoured ROIs provided histograms from which position, dispersion, and form parameters were computed. Histogram data were correlated with histopathology, p16-immunohistochemistry, and polymerase chain reaction for human papillomavirus DNA. RESULTS: There were 21 human papillomavirus–positive and 84 human papillomavirus–negative head and neck squamous cell carcinomas. At histopathology, human papillomavirus–positive cancers were more often nonkeratinizing (13/21, 62%) than human papillomavirus–negative cancers (19/84, 23%; P = .001), and their mitotic index was higher (71% versus 49%; P = .005). ROI-based mean and median ADCs were significantly lower in human papillomavirus–positive (1014 ± 178 × 10−6 mm2/s and 970 ± 187 × 10−6 mm2/s, respectively) than in human papillomavirus–negative tumors (1184 ± 168 × 10−6 mm2/s and 1161 ± 175 × 10−6 mm2/s, respectively; P < .001), whereas excess kurtosis and skewness were significantly higher in human papillomavirus–positive (1.934 ± 1.386 and 0.923 ± 0.510, respectively) than in human papillomavirus–negative tumors (0.643 ± 0.982 and 0.399 ± 0.516, respectively; P < .001). Human papillomavirus–negative head and neck squamous cell carcinoma had symmetric normally distributed ADC histograms, which corresponded histologically to heterogeneous tumors with variable cellularity, high stromal component, keratin pearls, and necrosis. Human papillomavirus–positive head and neck squamous cell carcinomas had leptokurtic skewed right histograms, which corresponded to homogeneous tumors with back-to-back densely packed cells, scant stromal component, and scattered comedonecrosis. CONCLUSIONS: Diffusion phenotypes of human papillomavirus–positive and human papillomavirus–negative head and neck squamous cell carcinomas show significant differences, which reflect their distinct degree of tumor heterogeneity.

I nfection with certain types of human papillomavirus (HPV), particularly HPV16 and HPV18, is a well-established cause for the development of head and neck squamous cell carcinoma (HNSCC) arising in the oropharynx and, to a lesser extent, also in other head and neck sites. During the infection process, the viral DNA is incorporated in the cell DNA and induces genetic modifications responsible for carcinogenesis. 1 The HPV oncoprotein E7 binds to and degrades the retinoblastoma protein, which leads to the accumulation and overexpression of p16 protein. 2 P16 immunohistochemistry appears to be the best stand-alone test for HPV testing of HNSCC because it is widely available, very sensitive (up to 100%), and specific for the presence of transcription-ally active HPV. Polymerase chain reaction is also very sensitive; however, its specificity is limited because the mere presence of HPV DNA does not convey whether the virus is transcriptionally active or simply a bystander. 2,3 In the oropharynx, HPV-and/or p16-positive (HPVϩ) HNSCC has a better prognosis than HPV-and/or p16-negative (HPVϪ) HNSCC, 4 whereas in the oral cavity and in other head and neck sites, reported results are contradictory. [4][5][6] Morphologic imaging techniques do not allow distinction between HPVϩ and HPVϪ cancers. Some authors found that HPVϩ HNSCCs tend to have cystic nodal metastases and extracapsular spread, whereas others reported that morphologic criteria (size, T stage, and cystic nodal disease) and sociodemographic factors were not reliable predictors of HPV status. 7,8 DWI is increasingly used for the pretherapeutic work-up and monitoring of treatment response in HNSCC. On DWI, HNSCC displays restricted diffusivity because of increased cellularity, and mean ADC values are in the range of 900 -1300 ϫ 10 Ϫ6 mm 2 /s. [9][10][11][12] These reported ADC values are based on a monoexponential fitting model with 2 b-values (b ϭ 0 and b ϭ 800 -1000 s/mm 2 ). [9][10][11]13 Mean ADC values are often calculated based on a single freehand ROI drawn over the largest axial tumor area, 14 or alternatively, they are based on elliptical ROIs drawn over the "most cellular" tumor regions, 11 the choice of these cellular regions being made visually. In contrast to simple ROI techniques, pixel-based calculations of ADC maps have the advantage of providing a detailed view of the entire area of interest in the form of parametric maps and histograms offering the possibility to assess tumor heterogeneity in a noninvasive fashion. 15 The assessment of tumor heterogeneity plays a major role in cancer research because it influences treatment outcome and has an impact on the development of targeted therapies. Some authors have reported the added value of ADC histograms as predictors of treatment response, whereas others have noted an association between low pretreatment mean ADC and favorable outcome after radiochemotherapy. 13,16,17 Because the HPV status was not included in the analysis, it is not known whether ADC characteristics associated with improved outcome could actually be attributed to the HPV status.
HPVϩ and HPVϪ HNSCCs are distinct tumors with different molecular pathogenesis. At histopathology, they diverge in terms of microstructural heterogeneity. 18, 19 We hypothesized that diffusion phenotypes based on a histogram analysis of ADC values reflect different degrees of tumor heterogeneity in HPVϩ and HPVϪ HNSCCs. The purpose of this study was to analyze the relationship between ADC histogram parameters and histopathologic features of HPVϩ and HPVϪ HNSCCs and to provide an understanding of the histologic underpinning of this quantitative tool.

Patient Selection and Inclusion Criteria
This retrospective study was approved by the Institutional Ethics and Research Committee of the University Hospitals Geneva and was performed in accordance with the guidelines of the Helsinki II declaration. Informed consent was waived. Our institution is a tertiary referral center for head and neck oncology and has a weekly multidisciplinary head and neck oncology tumor board.
The inclusion criteria for the current study were as follows: 1) adult patients seen at the multidisciplinary head and neck tumor board, 2) primary HNSCC of the oropharynx and oral cavity, 3) MR imaging examinations (standard of care for HNSCC in our institution) obtained with morphologic sequences and singleshot spin-echo echo-planar DWI with b ϭ 0 and b ϭ 1000 s/mm 2 , 4) unequivocal categorization of the HPV status, and 5) histopathology. We included both the oropharynx and the oral cavity in this series to avoid selection bias because though some tumors were classified as originating from the oral cavity, they showed involvement of both subsites and/or originated at the anatomic limit between the 2 subsites. Moreover, recent data suggest that ectopic tonsillar tissue may actually be the origin of HPVϩ HN-SCC in nonoropharyngeal sites (eg, the floor of the mouth). 18,19 In addition, because HPVϩ HNSCC tends to arise from the depths of tonsillar crypts and then spread submucosally, emerging at a distant site, the tumors can be misclassified as primary nonoropharyngeal HPVϩ HNSCC though they are oropharyngeal in origin. 18, 19 A computerized search of the medical records of the multidisciplinary head and neck tumor board over a period of 3 years yielded a total of 112 consecutive patients meeting the inclusion criteria. Three patients had to be excluded because of poor quality DWI caused by severe motion artifacts and metal implants. Four further patients were excluded because the lesions were not seen on diffusion-weighted MRI (2 in situ lesions, 2 T1 lesions according to the Union for International Cancer Control TNM classification). Therefore, a total of 105 consecutive patients formed the basis of the current study.

MR Imaging Protocol
The MR imaging examinations were performed on 1.5T (n ϭ 45) and 3T (n ϭ 60) MR imaging scanners provided by 2 different manufacturers. There was no statistically significant difference between the number of patients examined at 1.5T and at 3T (P ϭ .22). For both field strengths, our imaging protocols included axial T2, coronal STIR, axial T1 before and after intravenous injection of a gadolinium-based contrast agent, and fat-saturated contrast-enhanced T1 in the axial and coronal plane. A singleshot spin-echo echo-planar DWI sequence was performed in all patients before intravenous contrast material with 2 b-values (b ϭ 0 and b ϭ 1000 s/mm 2 ) and with similar parameters at 1.5T (TR, 7300 ms; TE, 75 ms) and 3T (TR, 7500 ms; TE, 75 ms), respectively. Fat suppression was used for all DWI sequences (spectral fat saturation or STIR with TI, 170 -180 ms at 1.5T and 230 ms at 3T). Depending on patient morphology, the FOV ranged from 220 -260 mm with corresponding acquisition matrices of 77-192 voxels in height and 112-208 voxels in width. The section thickness was 3-4 mm in 100 cases and 5 mm in 5 patients. All ADC maps were calculated by monoexponential fitting. [9][10][11][12][13]

Image Analysis, Data Processing, and Criteria for the Interpretation of Histograms
Two radiologists evaluated the MR imaging examinations independently. They were blinded to clinical data, histopathologic results, and HPV status. OsiriX MD version 3 (http:// www.osirixviewer.com) was used for image interpretation and measures. The readers were experienced radiologists with Ͼ8 years of expertise 2154 de Perrot Nov 2017 www.ajnr.org in oncologic head and neck MR imaging. Each tumor was identified by using morphologic MR imaging and DWI. The diagnosis of HNSCC was made according to established criteria (mass lesion with high signal on B1000 and low signal on a gray-scale ADC map). 10, 11 After scrolling through all the images, the 2 largest consecutive cross-sectional tumor areas were selected and contoured manually on axial B1000 and corresponding ADC maps by using freehand ROIs. Because DWI shows spatially nonlinear distortions, we did not contour tumor ROIs on morphologic MR imaging, but directly on DWI to avoid including nontumoral pixels in the analysis. Then, the freehand ROIs were saved individually. Although freehand ROIs of the entire tumor may be more representative of tumor heterogeneity, we deliberately contoured only the 2 largest axial tumor areas because we wanted to be certain about tumor margins. In our experience, the identification of tumor margins on the most cranial and most caudal sections may be fraught with errors (partial volume effect with nontumoral tissues and nonlinear image distortions), especially in smaller lesions. Because we were interested in tumor heterogeneity, we contoured the entire cross-sectional tumor area without excluding necrotic portions. The voxels within the contoured tumor areas were exported in an .xml file by using the OsiriX plugin "Export ROIs," and they were retrieved by the statistical software R (http://www.r-project.org/). 20 For each set of voxels belonging to a tumor, the software calculated a histogram. For each tumor, the following histogramderived ADC parameters were obtained: mean and median values (relevant for the histogram center); minimum; maximum; 25% quartile; 75% quartile; interquartile range (value spread around the center); skewness and kurtosis (histogram symmetry and peakedness). The interpretation of histogram parameters was done in accordance with the literature. 15 Kurtosis in a normal distribution is 3. Excess kurtosis (raw kurtosis minus 3) is 0 for a normal distribution. Positive numbers represent "narrow" leptokurtic histograms, and negative numbers correspond to "wide" platykurtic histograms. Skewness associated with a normal and symmetric distribution is 0. Positive skewness represents a histogram with a right longer tail (skewed right), whereas negative skewness represents a histogram with a left longer tail (skewed left). 15 If skewness is ϽϪ1 or Ͼϩ1, the distribution is highly skewed. If skewness is between Ϫ1 and Ϫ0.5 or between ϩ0.5 and ϩ1, the distribution is moderately skewed. If skewness is between Ϫ0.5 and ϩ0.5, the distribution is approximately symmetric.

Standard of Reference
For each patient, the HPV status was determined based on the tissue samples obtained by endoscopic biopsy (n ϭ 52) or surgical resection (n ϭ 53). A 2-step procedure consisting of p16 immunohistochemistry and polymerase chain reaction was used to differentiate between HPVϩ and HPVϪ HNSCC. 21 Tumors were considered p16ϩ in the presence of a strong nuclear and cytoplasmic staining of Ն70% of the neoplastic cells. 2 In the presence of p16 overexpression, polymerase chain reaction for the detection of HPV DNA was carried out. Cases were categorized as HPVϩ according to p16 immunohistochemistry and polymerase chain reaction results according to American Joint Committee on Cancer 2017 guidelines. 3 A degree of histologic differentiation was not assigned to HPVϩ HNSCCs based on American Joint Committee on Cancer 2017 recommendations. 3 A semiquantitative histologic score (On-line Table 1) was used to facilitate correlation between tumor microarchitecture and ADC histograms. This semiquantitative score took the following features into consideration: shape and size of tumor cells and nuclei, quantitative and qualitative stromal characteristics (ie, desmoplasia, lymphocytic infiltration), intratumoral keratinization, mitotic index estimated by MIB-1 immunohistochemistry, and the size and distribution of necrosis.

Statistical Analysis
All statistical analyses were performed by using the R statistical environment version 3.0.2. 20 For all tests, the significance level was fixed at .05. In the patient characteristics analysis, categoric variables were compared by using the 2 test with Yates continuity correction for data grouped in a 2 ϫ 2 contingency matrix. If the expected frequency in a cell was lower than 5, the Fisher exact test was applied. For continuous variables (eg, patient age), a nonparametric 2-sided unpaired Wilcoxon rank sum test with continuity correction was used. Multivariate logistic regression analysis was obtained to correct for possible confounding factors in the correlation between HPV status and ADC. For quantitative parameters derived from tumor histograms, a comparison between HPVϩ and HPVϪ HNSCCs was performed by using the Wilcoxon rank sum test mentioned above.

RESULTS
Patient and tumor characteristics are summarized in Table 1. The prevalence of HPVϩ HNSCCs in the oropharynx (18/ 52, 35%) was significantly higher (P Ͻ .001) than in the oral cavity (3/53, 6%). Based on the standard of reference, HPVϩ tumors were less often keratinized than HPVϪ tumors, and the tumor cells were more often basaloid (with limited cytoplasm) or undifferentiated. HPVϩ HNSCCs showed a higher mitotic index, a higher proportion of histologic score 4/5, and a higher proportion of nodal metastases than HPVϪ cancers. These differences were statistically significant (Table 1). Mean and median (95% CI) ROI size for ADC measurements were 7.7 cm 2 (range, 2.3-56 cm 2 ) and 5.72 cm 2 (range, 2.84 -10.58 cm 2 ), respectively. Comparison between HPVϩ and HPVϪ HNSCC at both 1.5T and at 3T ( Table 2) revealed significant differences between histogram parameters, in particular for mean and median ADC, skewness, and excess kurtosis (P Ͻ .05). Comparison of histogram parameters of HPVϪ HNSCCs at 1.5T versus 3T revealed no significant difference between the 2 field strengths (P Ͼ .05). Similarly, histogram parameters of HPVϩ HNSCCs showed no significant difference between 1.5T and 3T, respectively (P Ͼ .05). The comparison between HPVϩ and HPVϪ HNSCC for the pooled data irrespective of field strength and subsite is summarized in Table 2 and Fig 1. Comparison between HPVϩ HNSCC and HPVϪ HNSCC in the oropharynx revealed significant differences for all histogram parameters (P Ͻ .05), whereas in the oral cavity, only a trend was observed because of the small number of HPVϩ cases. Comparison between HPVϩHNSCC histograms in the oropharynx versus oral cavity revealed no significant differences (P Ͼ .05). Multivariate logistic regression analysis showed that the correlation between HPV status and mean ADC was independent of age, sex, tumor site, and T, N, and histologic parameters (On-line Table 2).      Figure). The observed leptokurtic and skewed right histogram in HPVϩ HNSCC corresponded histologically to homogeneous tumors with high cellularity, back-to-back homogeneous cell arrangement with little interstitial space, basaloid cells with scant cytoplasm, few to absent keratin pearls, and a higher mitotic index (Figs 3 and 4 and the On-line Figure). These histologic features of HPVϩ cancers generated a high number of pixels with low ADCs and, therefore, higher kurtosis. The significantly higher positive (right) skewness in HPVϩ cancers was caused by scattered comedonecrosis, which resulted in a small number of pixels with high ADCs (driving the mean upwards but not affecting the median of the histogram).

DISCUSSION
Although tobacco and alcohol consumption are the most common risk factors for the development of HNSCC, HPV has been recognized as the main etiologic factor in a subgroup of patients.
The incidence of HPVϩ HNSCC and the proportion of HPVϩ to HPVϪ HNSCCs depend on tumor localization, geographic and socioeconomic factors, smoking status, and sexual behavior. Most authors currently agree that HPVϩ HNSCC represents a distinct tumor entity with a more favorable prognosis than alcohol-or tobacco-related HNSCC. Recent publications have attempted to determine whether mean ADC values calculated by monoexponential fitting could be used to determine the HPV status in HNSCC. 12,14,22 In contrast to the intravoxel incoherent motion biexponential model, the monoexponential model does not differentiate between molecular movement caused by microperfusion and true diffusion caused by Brownian movement. 23 To obtain a robust fitting, biexponential models require the acquisition of 6 -12 b-values, leading to increased examination time and a higher percentage of motion artifacts. Biexponential models are, therefore, not widely used in clinical practice. Publications evaluating the association between mean ADC values with monoexponential fitting and the HPV status in HNSCC have yielded contradictory results; although some authors suggested that mean ADCs at 1.5T are lower in HPVϩ than in HPVϪ HNSCCs, others found no significant differences. 12,14,22 In these studies, the methodology of ROI contouring, the type of DWI sequence used, and the range and num-  22 we found that mean ADCs were significantly lower in HPVϩ than in HPVϪ cancers. However, our data are in contradiction to results reported by Schouten et al, 14 who did not find a statistically significant difference between ADCs of both tumor types. 14 A possible explanation is the fact that Schouten et al 14 used a PROPELLER sequence, whereas the other authors, including us, used EPI-based DWI.
Imaging-based first-order histogram analysis of ADC and of dynamic contrast-enhanced MR imaging parameters can be used to measure organ and tumor heterogeneity in a noninvasive fashion. 15 Histogram analysis has been applied for the discrimination between histologic tumor types in renal cell cancers, brain tumors and HNSCC. 15,16,24,25 Histogram analysis of pretreatment ADCs measured at b ϭ 2000 s/mm 2 was able to distinguish between well, moderately, and poorly differentiated HNSCCs, whereas analysis based on b ϭ 1000 s/mm 2 was not able to do so. 25 King et al 16 used ADC histogram analysis to evaluate the capability of DWI to predict treatment outcome after radiochemotherapy of advanced HNSCC, and Shukla-Dave et al 26 used histogram analysis of dynamic contrast-enhanced MRI parameters to predict outcome in HNSCC with nodal metastases.
As suggested by Just, kurtosis and skewness can be regarded as "quantitative surrogate markers of tumor heterogeneity." Kurtosis reflects distribution peakedness, a narrow distribution being linked to a low voxel value variation and, therefore, to a homogeneous tumor matrix. Skewness measures distribution asymmetry, the meaning of this metric being more difficult to interpret. 15 Data from the current study reveal that ADC histograms of HNSCC reflect the histologic microarchitecture unique to the HPV status. Histologically, in HPVϩ tumors, cells are organized in homogeneous clusters with little interstitial space. This finding explains the lower mean and median ADC in HPVϩ cancers because the reduced intercellular space restricts diffusion of water molecules. The homogeneous microstructure without major regional variability explains our observation that histograms of HPVϩ HNSCCs are more slender (leptokurtic) than histograms with a classic Gaussian distribution (as seen in HPVϪ cancers). In other terms, higher kurtosis reflects higher tumor matrix homo- geneity. The significantly higher positive (right) skewness in HPVϩ cancers is caused by a small amount of pixels with high ADCs and a large amount of pixels with low ADCs. This ADC distribution corresponds histologically to large tumor areas with densely packed cells and only minor amounts of scattered comedonecrosis. Kurtosis and skewness thus demonstrate that the ADC distribution in HPVϩ cancers tends to be non-Gaussian, whereas HPVϪ cancers tend to have Gaussian histograms. This difference can be attributed to the higher heterogeneity of HPVϪ cancers (more scores 1-3 in this series). Data of our study, therefore, support the hypothesis that the distinct histologic microarchitecture of HPVϩ and HPVϪ HNSCCs can be elegantly demonstrated in a noninvasive fashion by ADC histogram analysis. Furthermore, the observed differences between histograms in HPVϩ and HPVϪ cancers were independent of MR imaging field strength and manufacturer.
Commercially available software packages are already able to provide histogram analysis of ROIs contoured during routine image interpretation. ADC histograms could, therefore, be used as a noninvasive triage tool for HPV testing by identifying HPVϩ HNSCC located in nonoropharyngeal sites, where p16/polymerase chain reaction testing is not obtained routinely. High kurtosis, positive skewness, and low mean ADC would warrant p16/polymerase chain reaction testing. In addition, ADC histograms could also serve as surrogate indicators of HPVϩ HNSCC in cases with insufficient tissue remaining for p16/polymerase chain reaction testing and at the time of intraoperative frozen section evaluation.
HPVϩ HNSCC has been linked to a better outcome after radiochemotherapy. Some authors have shown that low mean ADC and/or histograms wherein Ͼ45% of voxels have ADCs Ͻ1.15 ϫ 10 Ϫ3 mm 2 /s indicate a favorable response to radiochemotherapy. 13,16,17 Although the HPV status was not analyzed in these studies, our data suggest that the reported predictive value of ADC histograms and lower mean ADC may in part be related to the HPV status. Therefore, histogram analysis used to monitor treatment response and target new therapeutic approaches should take into consideration the distinct diffusion phenotypes of HPVϩ and HPVϪ HNSCC as revealed by this study.
The current study has limitations inherent to its retrospective study design. Limitations include section thickness variability and the small number of HPVϩ HNSCC cases in the oral cavity. Further limitations concern freehand ROI contouring (as opposed to automatic segmentation) and ROI contouring of 2 tumor sections only (as opposed to contouring all tumor sections). Automatic tumor segmentation would in theory produce more robust lesion delineation and require less experienced readers. 27 Nevertheless, the method we used is similar to ROI-contouring methods published in other studies. 27

CONCLUSIONS
Quantitative ADC histogram analysis reveals significant differences between HPVϩ and HPVϪ HNSCCs irrespective of the MR imaging field strength used. The observed difference in ADC histograms reflects the distinct histopathologic features of these tumors, especially their different degree of microstructural heterogeneity. Because DWI is increasingly used to predict and monitor treatment outcome in HNSCC, diffusion phenotypes of HPVϩ and HPVϪ tumors as demonstrated in this study should be considered in future research protocols.