Cutaneous Human Papillomaviruses and the Risk of Keratinocyte Carcinomas

Keratinocyte carcinomas (KC), composed primarily of cutaneous squamous cell carcinoma (cuSCC) and basal cell carcinoma (BCC), are the most common skin malignancies in the United States, with approximately 5.4 million KC cases diagnosed annually. Although KC-associated mortality is low, patients are often diagnosed with multiple KCs, resulting in significant morbidity and treatment costs of approximately $8 billion annually (1–3). Established risk factors for cuSCC and/or BCC include age, male sex, light-colored skin, hair, and eyes, immunosuppression, and past ultraviolet radiation (UVR) exposure (4), many of which have been incorporated into KC risk prediction models (5–7). Additional biomarkers of recent, potentially modifiable exposures are needed to inform novel primary prevention strategies and personalized skin cancer screening and treatment regimens.

Infections with cutaneous human papillomavirus (cuHPV) have been associated with KC in several epidemiologic studies, particularly types in genus β with cuSCC (8–19). β-HPV detection could be a biomarker of future KC risk with clinical utility for personalizing KC prevention. However, data needed to inform such risk models are limited, as most previous epidemiologic studies of β-HPV and KC have been retrospective or cross-sectional in design, many relying on markers of past infection (serology) (9–13, 15–17, 19), with fewer incorporating measures of recent infection (viral DNA-based markers; refs. 8, 9, 15, 16, 18). Furthermore, no previous epidemiologic studies have examined quantitative measures of recent UVR exposure as predictors of KC risk in the context of cuHPV, and limited studies have investigated KC risk associated with the growing number of cuHPV types identified in genus γ (11, 13, 16, 19, 20).

The Viruses in Skin Cancer (VIRUSCAN) Study was specifically designed to investigate KC risk associated with cuHPV in an immunocompetent population, using a cohort design and multiple biomarkers of past and recent infection with β- and γ-HPV, as well as quantitative measures of recent UVR exposure (21). Importantly, the VIRUSCAN Study was uniquely poised to directly assess the development of virus-positive KC following the detection of cuHPV at baseline, analogous to the seminal studies of high-risk HPV and cervical cancer (22, 23).

Methods are summarized here with details provided in Supplementary Materials and Methods.

The VIRUSCAN Study design and cohort participants have been described previously (21, 24, 25). Briefly, individuals undergoing routine skin cancer screening exams were recruited from the University of South Florida Dermatology Clinic (USFDC) in Tampa, Florida, in July 2014 to August 2017. Eligible patients were age ≥60 who had not had a history of both cuSCC and BCC at the time of study enrollment based on medical record review. A total of 1,179 patients consented to participate and were confirmed not to have a prevalent KC at the time of enrollment (21). Patients underwent total body skin examinations (TBSE) every 6 to 12 months as part of their routine care at USFDC, and VIRUSCAN Study visits were coupled with these routine TBSE through September 2018 (21). Six patients (0.5%) withdrew, and 165 (14.0%) did not return to the clinic during the follow-up period, resulting in 1,008 participants contributing 2,157 person-years of follow-up time. The mean follow-up time was 788 days, with a median of 792 days (interquartile range = 508–1,081 days). All patients provided written informed consent in accordance with the Belmont Report and study methods were approved by the USF Institutional Review Board.

At baseline, a blood sample, eyebrow hair (EBH), and forearm skin swab (SSW) were obtained from participants (21, 25). Pathology reports corresponding to biopsies conducted during follow-up were routinely reviewed for KC diagnoses. Recuts of formalin-fixed, paraffin-embedded KC tumors were obtained from the USF Pathology Department. The study dermatopathologist confirmed all KC diagnoses, and in the adjacent normal tissue, graded the degree of solar elastosis, a measure of cumulative UVR exposure. Solar elastosis was graded on a 4-point scale ranging from no evidence of UV-radiation induced damage (0) to the highest degree of solar elastosis (3): 0: no presence of elastotic fibers or very occasional elastotic fibers only discernible with 20× power; 1: scattered elastotic fibers lying as individual units between collagen bundles; 2: densely scattered elastotic fibers distributed predominantly as bushels rather than individual units; 3: dense aggregates of elastotic fibers forming amorphous deposits of blue-gray material with lost fiber texture, as described previously (21, 26). Incident cuSCC were graded for the histologic features of papillomatosis, hypergranulosis, and crateriform architecture (present vs. absent).

Baseline serum samples were analyzed at the German Cancer Research Center (21, 27) using multiplex detection methods for assessment of antibodies to the major protein L1 corresponding to 17 β and 7 γ-HPV types. The cuHPV types selected for serological testing were those previously shown to be associated with KC in the literature and/or the most prevalent types detected in KC tumors in our previous case-control study (13, 18). A Luminex analyzer (Luminex Corp.) was used to identify the internal bead color and quantify the fluorescence expressed as median fluorescence intensity (MFI). At the International Agency for Cancer Research, EBH, SSW, and tumor recuts underwent DNA extraction and cuHPV detection (21, 25, 27) for 46 β and 52 γ types using multiplex polymerase chain reaction/Luminex assays (28, 29).

Demographic and skin cancer risk factors were assessed through self-administered questionnaires at baseline (21). Quantitative estimates of recent UVR exposure were obtained by measuring the difference in skin pigmentation of the inner arm (one's natural skin tone) and top of the forearm (sun-exposed skin) using the Konica Minolta CM-600 spectrophotometer (21, 24, 27). This “UVR indicator” has been shown previously to correlate with sun-exposure (21).

A power analysis was conducted during the design phase of the VIRUSCAN Study, which was based on published incidence rates of KC (30) and cuHPV prevalence, as described previously (21). Analyses were conducted separately for cuSCC and BCC. Baseline demographic characteristics including age, sex, race/ethnicity, and other skin cancer risk factors were described for cohort participants who did or did not develop KC. Log-rank P values were calculated to compare time to KC incidence between categories of each skin cancer risk factor. Potential confounders of the associations between cuHPV and KC were identified through backward stepwise elimination.

cuHPV prevalence was calculated for each baseline biomarker including serology, EBH, and SSW, across two viral DNA-based biomarkers (EBH and SSW), and across all three biomarkers. The number of cuHPV types detected in each viral biomarker was calculated. Cumulative KC incidences were plotted and stratified by any β- and γ-HPV as measured by each biomarker, with log-rank P values calculated to compare the incidence between virus-positive and virus-negative groups. Two sets of Cox proportional hazards models were constructed to estimate the risk of KC associated with baseline cuHPV prevalence as follows, separately for cuSCC and BCC: (i) minimum model, adjusting for age and sex; (ii) full model, adjusting for all potential confounders previously identified. UVR exposure and spectrophotometer-based natural skin tone were assessed as potential effect modifiers of the association between cuHPV and KC development through: (i) stratification of the minimum model by higher versus lower UVR indicator, and darker versus lighter natural skin tone; (ii) inclusion of interaction terms in the minimum model. All models were examined for violation of the proportional hazards assumption. Model fitness was compared using analysis of deviance.

KC risks associated with detection of cuHPV at baseline with specific β-HPV types in specific phylogenetic branches were also assessed. To avoid multiple comparisons, a parallel gate-keeping strategy was designed a priori to analyze the association between cuSCC and β-HPV at the branch- and type-level.

Pie charts were used to visualize and compare the β-HPV type distributions detected in each DNA marker. Prevalence of cuHPV in incident KC tumors was calculated and reported. Detection of cuHPV at baseline was considered concordant with a subsequent virus-positive tumor if the tumor was positive for the same cuHPV type. MFI values were plotted for those that did or did not develop a concordant tumor and compared across groups using the Wilcoxon rank-sum test. All tests were two-sided using R Statistical package version 3.5.0 from November 2018 to December 2020.

Among the 1,008 participants followed, a total of 149 cuSCC and 132 BCC diagnoses were observed during follow-up, corresponding to incidence rates of 7.52 and 6.58 per 100 person-years for cuSCC and BCC, respectively. Participants who were male, and those self-reporting a history of occupational sun exposure, greater sun sensitivity, lighter phenotypic characteristics such as lighter skin color and hair color, and a history of KC at study enrollment were significantly more likely to develop cuSCC, whereas participants who were older, male, had a greater number of moles on the body, and had a self-reported history of KC at study enrollment were more likely to develop BCC (Table 1).

CuSCC incidence was higher among participants who were β-HPV DNA-positive in EBH or SSW at baseline compared with DNA-negative participants, whereas no clear separation of curves was observed by cuHPV baseline serostatus (Fig. 1A). BCC incidence patterns were similar (Fig. 1B). No differences in cuSCC or BCC incidence were observed by γ-HPV serostatus or presence/absence of γ-HPV DNA in EBH (Supplementary Figs. S1A and S1B); BCC incidence was slightly elevated among those with versus without γ-HPV DNA in SSW (Supplementary Fig. S1B).

Based on corresponding hazard ratios (HR), no statistically significant risks of KC were observed in association with baseline cuHPV seropositivity or viral DNA in EBH (Table 2). Individuals with any β-HPV DNA in their baseline SSW were four times more likely to develop cuSCC as β-HPV DNA SSW-negative individuals, after adjustment for age and sex [HR = 4.32; 95% confidence interval (CI), 1.00–18.66; P = 0.05]. Similarly, presence of any β-HPV in species 1 was associated with significant risk of cuSCC (HR = 1.92; 95% CI, 1.15–3.21; P = 0.01). HRs were slightly attenuated after additional adjustment for self-reported skin color, hair color, and history of KC. No associations were observed between β-HPV species 2 or γ-HPV in SSW and risk of cuSCC. Associations between cuHPV in SSW and cuSCC did not significantly differ between individuals with lighter versus darker natural skin tone (Supplementary Table S1).

A positive trend was observed between the number of cuHPV types detected in baseline SSW and risk of cuSCC, specifically for β-HPV species 1 (Table 3), with those testing positive for 6 to 11 types having twice the risk of subsequent cuSCC compared with those with no β-HPV species 1 in SSW (HR = 2.35; 95% CI, 1.22–4.53; P = 0.01). A significant trend was also observed for γ-HPV and cuSCC, although no clear pattern of dose-response was observed. No significant trends were observed between numbers of cuHPV types in SSW and BCC (Table 3).

KC risks associated with presence of baseline cuHPV across multiple biomarkers are presented in Supplementary Table S2. Individuals with at least one β-HPV species 1 type detected in both their EBH and SSW were more than twice as likely to develop cuSCC compared with those who had no β-HPV species 1 types detected in either EBH or SSW (HR = 2.13; 95% CI, 1.19–3.84; P = 0.01; Supplementary Table S2). No associations were observed between presence of viral DNA in both EBH and SSW and risk of BCC. Finally, no additional information was gained by adding serology into the prediction models (Supplementary Table S3). When comparing the goodness of fit across models, the models that included β-HPV species 1 detected in EBH&SSW, with (P = 0.77) or without seropositivity (P = 0.17), did not exhibit significantly superior fitness to the model that included β-HPV species 1 detected in SSW alone.

cuHPV species 1 branch- and type-specific cuSCC risks were calculated using a gate-keeping strategy (31) to minimize false discovery rates (Supplementary Figs. S2A and S2B). Branch 1A and 1D were the two most prevalent branches in β-HPV species 1, with only Branch A being significantly associated with cuSCC (HR = 1.74; 95% CI, 1.23–2.46; P_adjusted = 0.003). Within Branch 1A, HPV 24 and HPV 98 were the two most prevalent types, with only HPV 24 being significantly associated with cuSCC (HR = 1.63; 95% CI, 1.14–2.35; P_adjusted = 0.04). HPV types 5 and 8 were most prevalent in Branch 1D, and HPV 5 was associated with cuSCC at a significance level of 0.058 prior to multiple comparisons adjustment. Neither of the two most common species 2 branches (2D and 2E) were significantly associated with cuSCC and, therefore, cuSCC associations with specific types within those branches were not assessed. The branch-level association observed for cuHPV in SSW was not observed in EBH (HR_{branch 1A} = 1.20; 95% CI, 0.79–1.84; P_adjusted = 0.39).

The β-HPV type distribution was similar across SSW and EBH, with the most common individual types detected comprising 1% to 7% of all β-HPV detections, and all other types combined comprising 29% to 31% of all β-HPV detections. In contrast, specific β-HPV types comprised higher proportions of β-HPV detections among cuSCC (HPV 111, 145, 174) and BCC (HPV 110, 111, 174), with only 13% representing all other types combined in cuSCC (Fig. 2). The overall prevalence of any β-HPV DNA in incident KC tumors was 22.0% in cuSCC (Fig. 3A) and 8.6% in BCC tumors (Supplementary Fig. S3A). The percentage of baseline cuHPV types detected in EBH that were present in subsequently developed cuSCC was 1.67% (Fig. 3B), whereas only 0.62% of baseline β-HPV types detected in SSW were also present in subsequently developed cuSCC (Fig. 3B). Even lower proportions were observed for BCC (Supplementary Fig. S3B). Of note, baseline β-HPV types detected in EBH and SSW with higher MFI levels were more likely to be present in subsequent cuSCC as compared with those with lower MFI levels (Supplementary Fig. S4A–S4C); results were similar for EBH and BCC (Supplementary Fig. S5A–S5C). Nine participants developed more than one cuHPV-positive tumor, of whom eight had at least two tumors positive for the same cuHPV type and concordant with the cuHPV type detected in their baseline EBH or SSW (Supplementary Table S4).

No significant interactions between recent UVR exposure and β-HPV DNA in SSW at baseline were observed in relation to subsequent cuSCC. For example, the association between β-HPV species 1 in SSW and cuSCC was similar among those with lower UVR exposure (HR = 2.08; 95% CI, 0.71–6.09) versus higher UVR exposure (HR = 1.86; 95% CI, 1.04–3.33; P_interaction = 0.95). However, the highest degree of solar elastosis, a marker of cumulative sun damage in skin tissue, was significantly associated with the presence of β-HPV species 1 DNA in adjacent cuSCC tumor tissues (Table 4). Furthermore, β-HPV DNA was more likely to be present in cuSCC tumors exhibiting papillomatosis, a feature associated with viral warts, and crateriform architecture, a feature associated with keratoacanthoma, a variant of cuSCC in which HPV has been previously identified (32–34).

Here we report findings of a prospective study, the VIRUSCAN Study, that evaluates the incidence of cuSCC and BCC using different biomarkers for a broad spectrum of cuHPV types. Baseline β-HPV detection in SSW was a significant predictor of cuSCC risk. Overall, no clear patterns were observed between baseline β-HPV and subsequent development of BCC, or between baseline γ-HPV and either cuSCC or BCC.

The VIRUSCAN Study included viral DNA testing for 52 β-HPV and 46 γ-HPV types, which, to our knowledge, is more than any other study has included to date. Leveraging a novel gate-keeping strategy enabled us to examine cuSCC risk associated with cuHPV at the type- and branch-levels, while minimizing false discovery rates. Our finding of increased cuSCC risk associated with β-HPV types in branch A in SSW, and specifically with HPV24, is consistent with a meta-analysis that reported significant positive associations with type 24, as well as with types 5, 8, 15, 17, 20, 24, 36, and 38 (35). The lack of association between β-HPV positive EBH and cuSCC in this study is also consistent with some (36, 37), but not all (8, 9, 18), previous studies. Of note, no previous studies of immunocompetent individuals have examined KC risk associated with cuHPV in SSW, thus our observation that models incorporating SSW as a single biomarker of cuHPV infection were the most informative and parsimonious for prediction of subsequent cuSCC is novel.

We observed that past β-HPV infection measured by serology was not associated with subsequent cuSCC. These findings differ from positive associations reported by several case–control studies (9–11, 13, 15, 17) and one cohort study (12), but are consistent with a second cohort study (38). Serological assays included a smaller number of cuHPV types than the DNA-based assays in this study, although seropositivity for “any β-HPV types” can reflect exposure to the cuHPV types directly measured, as well as potentially cross-reactive types.

Several lines of evidence suggest UVR exposure and cuHPV act synergistically in the development of cuSCC, including experimental murine studies (39, 40) and observational human studies relying on recall of past UVR exposures (9, 17, 41). The VIRUSCAN Study leveraged a spectrophotometer to assess skin pigmentation levels as an indicator of recent UVR exposure, and this spectrophotometer-based UVR indicator was positively associated with all three markers of baseline β-HPV among VIRUSCAN participants (27). Here, we demonstrated that β-HPV was more likely to be detected in cuSCC tumors arising in the setting of cumulative UVR exposure, as indicated by higher levels of solar elastosis in skin tissue adjacent to β-HPV-positive tumors. The co-presence of solar elastosis and β-HPV DNA in cuSCC may reflect their cooperation in the carcinogenic process, as β-HPV E6 and E7 oncoproteins have been shown to inhibit cellular responses to UV radiation, including cell-cycle arrest and DNA repair, as well as apoptosis (42). Alternatively, chronic sun exposure leading to cuSCC could simply provide an optimal microenvironment for viral infection, with the virus playing no role in carcinogenesis. Assessment of β-HPV oncoprotein expression in the cuSCC tissues was not possible in the VIRUSCAN Study, due to the extremely small size of the FFPE screen-detected lesions ascertained.

We are the first to report on the concordance between cuHPV DNA detected at baseline and in subsequent skin tumors. Less than 5% of baseline β-HPV types detected in SSW were present in subsequent cuSCC tumors, with higher mean MFI values in SSW being associated with subsequent β-HPV-positive tumors, and β-HPV-positive tumors more often exhibiting papillomatosis (43) and crateriform architecture (44) two characteristics of HPV-associated skin lesions. Therefore, although a small percentage of cuHPV may be etiologically linked to cuSCC, our results do not support a direct role for β-HPV in the etiology and/or in the maintenance of most cuSCC. Previous studies have shown that β-HPV infections may be involved at an early stage of skin carcinogenesis and are dispensable for the maintenance of the cancer phenotype (45), with the diversity of β-HPV types detected gradually decreasing from healthy skin to premalignant lesions, actinic keratosis, and cuSCC (46, 47). Other studies reported the lack of viral RNA in cuSCC, despite the detection of β-HPV DNA, indicating that these viruses are not active in cuSCC (48). Our study was not designed to evaluate the possibility of a more indirect role for β-HPV in cuSCC etiology.

Study limitations include follow-up time, statistical power to detect cuHPV type-specific associations, information on actinic keratosis (AK) and generalizability of the cohort. Although this study did not examine mucosal HPV types in genus alpha, which have been previously detected in cutaneous skin (49), the measurement of a broad spectrum of cuHPV types in multiple genera (β- and γ-HPV) provided a global view of cuHPV in relation to both cuSCC and BCC. We were unable to examine cuHPV detection in AK, a precursor to cuSCC, because the standard of care at the USFDC is to treat AK lesions with cryotherapy, leaving no tissue for cuHPV analysis. However, no associations were observed between self-reported history of AK at baseline and β-HPV seropositivity (Chi-square P = 0.46) or β-HPV DNA detection in SSW (P = 1.00) or EBH (P = 0.98). VIRUSCAN Study participants are at higher risk of skin cancer than the general Southern U.S. population based on comparisons of skin cancer risk factor profiles (21) and KC incidence rates (50, 51), therefore comprising an ideal population in which to study novel skin cancer screening regimens/cuSCC prevention strategies incorporating cuHPV.

In conclusion, we demonstrated a positive association between presence of β-HPV and subsequent incidence of cuSCC. Importantly, the use of cuHPV detection in SSW as a biomarker of subsequent KC was novel, leading to the finding that cuHPV detected in SSW alone predicted subsequent cuSCC just as well as the combination of SSW and EBH detection, or SSW and EBH detection with serology. Although the positive predictive value of β-HPV in SSW was low in the study population overall (15%), it doubled among VIRUSCAN participants with a previous history of cuSCC, suggesting future studies should examine SSW β-HPV positivity in combination with other clinically relevant risk factors. Taken together, cuHPV detection could inform individual risk stratification for use in KC prevention and screening efforts, thus minimizing the incidence, morbidity, and/or healthcare costs associated with this most common cancer in the United States.

D.E. Rollison reports grants from NCI during the conduct of the study and other support from NanoString Technologies, Inc. outside the submitted work. R.P. Amorrortu reports grants from NCI during the conduct of the study. Y. Zhao reports grants from NCI at the NIH during the conduct of the study. M.J. Schell reports grants from NCI during the conduct of the study. N.A. Fenske reports grants from NIH during the conduct of the study. A.R. Guiliano reports grants and personal fees from Merck & CO, Inc., and personal fees from Immunomic Therapeutics outside the submitted work. V.K. Sondak reports personal fees from Merck, BMS, Eisai, Novartis, Regeneron, Replimune, and Alkermes, and other support from Neogene Therapeutics outside the submitted work. No disclosures were reported by the other authors.

D.E. Rollison: Conceptualization, resources, supervision, funding acquisition, methodology, writing–original draft, project administration. R.P. Amorrortu: Validation, visualization, writing–original draft, project administration. Y. Zhao: Data curation, formal analysis, validation, visualization, methodology, writing–original draft. J.L. Messina: Conceptualization, formal analysis, investigation, methodology, writing–review and editing. M.J. Schell: Conceptualization, methodology, writing–review and editing. N.A. Fenske: Resources, writing–review and editing. B.S. Cherpelis: Resources, writing–review and editing. A.R. Giuliano: Conceptualization, methodology, writing–review and editing. V.K. Sondak: Conceptualization, writing–review and editing. M. Pawlita: Conceptualization, resources, methodology, writing–review and editing. S. McKay-Chopin: Investigation, writing–review and editing. T. Gheit: Conceptualization, resources, investigation, methodology, writing–review and editing. T. Waterboer: Conceptualization, resources, investigation, methodology, writing–review and editing. M. Tommasino: Conceptualization, resources, investigation, methodology, writing–review and editing.

This work was supported by the NCI at the NIH [R01-CA177586] awarded to D.E. Rollison. This work was also supported in part by the Tissue Core and the Participant Research, Interventions, and Measurement (PRISM) Core at the H. Lee Moffitt Cancer Center and Research Institute, a comprehensive cancer center designated by the NCI and funded in part by Moffitt's Cancer Center Support Grant [P30-CA076292]. The funders were not involved in the study design, data collection/analysis, or manuscript preparation. The authors thank past VIRUSCAN Study team members including Juliana Balliu, Shalaka Hampras, Yessica C. Martinez, Rhianna Reed, and Laxmi Vijayan for their contribution to participant recruitment/retention, data collection, and study management. They would also like to thank Mary-Katherine Haver for her assistance with the literature search. Finally, the authors sincerely thank the USF Dermatology Clinic staff and patients who participated in the VIRUSCAN Study. Where authors are identified as personnel of the International Agency for Research on Cancer/WHO, the authors alone are responsible for the views expressed in this article and they do not necessarily represent the decisions, policy, or views of the International Agency for Research on Cancer/World Health Organization.

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.