Broad-Spectrum Detection of HPV in Male Genital Samples Using Target-Enriched Whole-Genome Sequencing

Most human papillomavirus (HPV) surveillance studies target 30–50 of the more than 200 known types. We applied our recently described enriched whole-genome sequencing (eWGS) assay to demonstrate the impact of detecting all known and novel HPV types in male genital samples (n = 50). HPV was detected in nearly all (82%) samples, (mean number of types/samples 13.6; range 1–85), and nearly all HPV-positive samples included types in multiple genera (88%). A total of 560 HPV detections (237 unique HPV types: 46 alpha, 55 beta, 135 gamma, and 1 mu types) were made. The most frequently detected HPV types were alpha (HPV90, 43, and 74), beta (HPV115, 195, and 120), and gamma (HPV134, mSD2, and HPV50). High-risk alpha types (HPV16, 18, 31, 39, 52, and 58) were not common. A novel gamma type was identified (now officially HPV229) along with 90 unclassified types. This pilot study demonstrates the utility of the eWGS assay for broad-spectrum type detection and suggests a significantly higher type diversity in males compared to females that warrants further study.


Introduction
Human papillomaviruses (HPV) are DNA viruses in the family Papillomaviridae that includes 49 genera.More than 200 HPV types that infect cutaneous and mucosal surfaces have been found in five of these genera, with most in alpha, beta, and gamma genera.While most infections clear, persistent infection with some alpha types (HPV16, 18,31,33,35,39,45,51,52,56,58,59, 66, and 68) is a risk factor for anogenital and oropharyngeal cancers [1].Beta and gamma types have not been studied systematically.However, recent reports have detected alpha, beta, and gamma types in oral samples [2], cervical cancers [3], head and neck cancers [4], anal and genital areas of males [5], and juvenile-onset recurrent respiratory papillomatosis (JoRRP) [6].In addition, interactions between HPV types in different genera may have biologic and clinical significance, as reported for the differential interaction of HPV16 with beta and gamma HPVs in head and neck cancers [7], the association of alpha and beta HPV coinfection with the clinical severity of JoRRP [6], and the detection of a high diversity of alpha, beta, and gamma HPVs in penile samples from males infected with HIV [8].Mechanistic studies suggest that some beta 3 types, such as HPV49, 75, and 76, can immortalize primary keratinocytes such as those immortalized with high-risk alpha HPV16 [9].Beta types may act as co-carcinogens with UV radiation for nonmelanoma skin cancers (NMSC) [10].Individuals with epidermodysplasia verruciformis (EV), a rare autosomal recessive hereditary skin disorder, have a high susceptibility to beta HPV infection and NMSC [11].Gamma HPVs can induce cutaneous lesions, and gamma 11 and 12 species were recently shown to be associated with head and neck squamous cell carcinoma (HNSCC) [4].The widespread occurrence of alpha, beta, and gamma HPVs in various anatomical sites indicates the importance of studies to clarify their natural history and pathogenesis, and the potential impact of HPV vaccines on their prevalence.
Currently, most studies have used PCR assays with degenerate primers (FAP59/64, BGC-PCR, and CUT primers) to detect and type HPV within specific genera [12,13].Epidemiological studies in men frequently use self-collected external genital swabs that sample both mucosal and cutaneous epithelia, yet most HPV typing assays only target mucosal types that are mainly alpha types [14].Next-generation sequencing (NGS)-based assays are making advances with viral genome research, and a few NGS studies, such as FAP PCR followed by 454 sequencing or Illumina MiSeq, report detection of alpha, beta, and gamma types using amplicons of limited sizes [8,15,16].NGS following rolling-circle amplification (NGS-RCA) has also been reported with moderate success [17].There are no NGS-based reports of detecting multiple HPV types belonging to different genera from any anatomical site in a single reaction.We recently developed a universal HPV assay based on target-enriched whole-genome sequencing (eWGS) to detect and type all known and potentially novel HPVs in a single reaction [18].Anticipating that samples from keratinized (cutaneous) and non-keratinized (mucosal) surfaces of male genital tract would have high HPV diversity, we used eWGS to detect and identify HPV in a pilot study of 50 male genital samples.

Samples
We used residual anonymized DNA extracts (n = 50) from external genital swabs in specimen transport medium collected from US males between 2010 and 2011 (IRB determination-Exempt).The vaccination history was unknown, but collection was prior to the US recommendation of routine HPV vaccination for males.DNA extracts were prepared using a Chemagic gDNA kit with an automated Chemagic MSM1 extractor (Perkin Elmer, Waltham, MA, USA) and stored at −80 • C for 63-68 months prior to use in the eWGS assay.The DNA concentration of the thawed extracts was determined using a Qubit dsDNA HS assay (Life Technologies, Eugene, OR, USA).
Due to variation in the input DNA (10-100 ng/sample; mean ± SD = 27.3 ± 23.6), HPV-negative placental DNA was added to bring the amount of DNA/sample to 100 ng prior to shearing (Supplemental Table S1).Libraries from each of the 63 samples (13 controls and 50 extracts) were prepared and barcoded as previously described [18].The barcoded samples were pooled (15-16 per pool) at equal molar amounts (pre-capture pooling) and enriched through hybridization (65 • C for 24 h) with HPV bait.Each pool contained at least one HPV-negative and one -positive control.The pooled captured libraries were amplified using 14 cycles of PCR, followed by purification with Ampure beads (Beckman Coulter, Indianapolis, IN).Quantity and quality of the captured libraries were assessed by a Bioanalyzer 2100 (Agilent Technologies, Santa Clara, CA, USA) and qPCR using the Kapa library quantification kit (Roche, Indianapolis, IN, USA) on a Light Cycler 480 (Roche Diagnostics, Indianapolis, IN, USA).The four target enriched libraries were combined at equal molar concentrations to make 2 sequencing libraries that were paired end sequenced at 3.2 pM on a two-lane flow cell using a Illumina HiSeq Rapid SBS kit (200 cycles) on a Hiseq 2500 (Illumina, San Diego, CA, USA).

Bioinformatics
The sequence data were analyzed using CLC Genomics Workbench (CLC Bio, Waltham, MA, USA).Read de-multiplexing, quality control, alignment, and HPV type determination were performed as described previously [18], except that reads were mapped at 90% similarity to 445 reference genomes (including both classified and unclassified types) available at PAVE (https://pave.niaid.nih.gov/)(as of 24 February 2022).A sample was considered positive for an HPV type if there were at least 100 mapped reads that covered at least 20% of a classified or unclassified reference HPV genome.This cut off is more stringent than recommended (minimum 10 reads and 10% genome coverage) for determining positive signal in HPV whole-genome studies [19].

Results
Sequencing libraries 1 (samples 1-31) and 2 (samples 32-63) generated similar cluster densities of 962 K/mm 2 and 995 K/mm 2 , with a total of 30, 935, and 33,083 mega base (Mb) sequence data, respectively.All samples with DNA generated an average of 10.5 million reads/sample with a mean Q score of 36.2;93.4% of the bases had a Q score ≥ 30.Results from analysis of these data are presented below.The 13 control samples gave the expected results (Supplementary Table S2).
Among the 122 alpha HPV detections, 46 unique types were identified in 32 samples.The number of reads/alpha genome was 88,996 ± 29,619 and >80% of the genome was sequenced in over 75% (94/122) of these detections.Among the 179 beta HPV detections, 55 unique beta types were identified in 35 samples.The number of reads/beta genome was 17,020 ± 4769, and >80% of the genome was sequenced in over 57% (102/179) of the detections.Among the 257 gamma HPV detections, 135 unique types were identified in 36 samples.The number of reads/gamma genome was 36,585 ± 9738, and >80% of the genome was sequenced in nearly 42% (107/257) of the detections.The two mu papillomavirus signals were mapped from two samples (reads/genome = 529 ± 444, genome coverage = 0.3 ± 0.12) to the same unclassified novel type, HPV-md01c06.

Identification of Novel and Unclassified Types
Interestingly, of the 237 unique HPV types identified by eWGS, there were no RNA baits in the enrichment hybridization for ninety-seven types (fourteen beta, eighty-two

Identification of Novel and Unclassified Types
Interestingly, of the 237 unique HPV types identified by eWGS, there were no RNA baits in the enrichment hybridization for ninety-seven types (fourteen beta, eighty-two gamma, and one mu types) and eighty-three of these types were unclassified.The study also identified a novel HPV genome (GenBank Accession # MW535770), that has subsequently been officially classified by the International HPV Reference Center, Sweden, as gamma HPV type 229.

Discussion
This pilot study was designed to demonstrate the capacity of the eWGS to identify all known and unknown HPVs in samples from male external genital swabs that include both mucosal and cutaneous epithelial surfaces.The expected universal capacity of the eWGS assay to detect HPV was supported by the detection of two-hundred thirty-seven unique types from four of the five Papillomaviridae genera recognized to have HPV types (forty-six alpha, fifty-five beta, one-hundred thirty-five gamma, and one mu types) (Figure 6).For three of the genera, multiple species were detected (fourteen species of alpha, five species of beta and twenty-three species of gamma).The broad typing capacity of the assay resulted in a high number of multiple types detected per sample (mean 13.6, range 1-85).This extensive diversity and high prevalence of HPV in male genital samples suggests further studies are warranted to evaluate the potential for males to serve as a reservoir for a large number of HPVs.It is interesting to note that of the five most prevalent high-risk types in a meta-analysis of data from women with normal cervical cytology in the pre-vaccine era (HPV16, 18, 31, 52, and 58) [20], HPV 58 was not detected in the male samples in this study.Low-risk alpha types HPV90 (18%), 43 (14%), and 84 and 74 (each at 12%) were more prevalent than high-risk alpha types HPV16, 18, and 31 (each at 8%) in our study.Alpha 3 types (HPV62, 84, and 89) were the most prevalent low-risk types in male genital samples in the National Health and Nutrition Examination Survey of adult US men [21] and women [22].A higher prevalence of HPV84 than of HPV16 has been reported in prior studies of men from culturally and geographically diverse regions of the world [23].Also, in agreement with that report, we did not identify HPV 61, a common vaginal type that is also an alpha 3 species like HPV 84 [24].

Figure 2 .
Figure 2. Alpha HPV in Male Samples.Prevalence by type (A) and species (B).All alpha types detected were officially classified.

Figure 3 .
Figure 3. Beta HPV in Male Samples.Prevalence by type (A) and species (B).A total of 55 unique beta types were detected that included 11 unclassified types.

Figure 2 .
Figure 2. Alpha HPV in Male Samples.Prevalence by type (A) and species (B).All alpha types detected were officially classified.

Figure 2 .
Figure 2. Alpha HPV in Male Samples.Prevalence by type (A) and species (B).All alpha types detected were officially classified.

Figure 3 .
Figure 3. Beta HPV in Male Samples.Prevalence by type (A) and species (B).A total of 55 unique beta types were detected that included 11 unclassified types.

Figure 3 .
Figure 3. Beta HPV in Male Samples.Prevalence by type (A) and species (B).A total of 55 unique beta types were detected that included 11 unclassified types.

Viruses 2023 , 10 Figure 5 .
Figure 5. Overall HPV Prevalence in Male Samples.Alpha, beta, and gamma types are indicated by orange, blue, and purple bars, respectively.

Figure 5 .
Figure 5. Overall HPV Prevalence in Male Samples.Alpha, beta, and gamma types are indicated by orange, blue, and purple bars, respectively.

Figure 6 .
Figure 6.Phylogenetic Tree of HPVs Detected in Male Samples.The tree is based on the L1 sequence of HPVs detected in this study, constructed using the tools at the PAVE EPSTEIME database.HPVs are grouped within alpha, beta, gamma, and mu (in green are gamma types with official classification).Tree branches in gray under gamma and beta indicate unclassified types.No unclassified types under alpha.