Identification of two novel papillomaviruses in belugas

Introduction: Papillomaviruses (PVs) can cause hyperplasia in the skin and mucous membranes of humans, mammals, and non-mammalian animals, and are a significant risk factor for cervical and genital cancers.

Methods: Using next-generation sequencing (NGS), we identified two novel strains of papillomavirus, PV-HMU-1 and PV-HMU-2, in swabs taken from belugas (Delphinapterus leucas) at Polar Ocean Parks in Qingdao and Dalian.

Results: We amplified the complete genomes of both strains and screened ten belugas and one false killer whale (Pseudorca crassidens) for the late gene (L1) to determine the infection rate. In Qingdao, 50% of the two sampled belugas were infected with PV-HMU-1, while the false killer whale was negative. In Dalian, 71% of the eight sampled belugas were infected with PV-HMU-2. In their L1 genes, PV-HMU-1 and PV-HMU-2 showed 64.99 and 68.12% amino acid identity, respectively, with other members of Papillomaviridae. Phylogenetic analysis of combinatorial amino acid sequences revealed that PV-HMU-1 and PV-HMU-2 clustered with other known dolphin PVs but formed distinct branches. PVs carried by belugas were proposed as novel species under Firstpapillomavirinae.

Conclusion: The discovery of these two novel PVs enhances our understanding of the genetic diversity of papillomaviruses and their impact on the beluga population.

1. Introduction

Papillomavirus (PV) can cause persistent infections in various parts of the skin and mucous membranes in humans and mammals, resulting in epidermal proliferative lesions (Cheng et al., 2023; Yang et al., 2023). Human genital PV (HPV) types are primarily responsible for cervical, vulvar, penile, and other anogenital cancers (Pino et al., 2023; Seyoum et al., 2023). Although papillomaviruses have been detected in all vertebrate populations except amphibians, most of them have been identified in humans and other mammals (Van Doorslaer et al., 2017). However, novel papillomaviruses have been recently detected in the cartilage of Fulmarus glacialis (Gaynor et al., 2015) and the liver of Hemidactylus frenatus (Agius et al., 2019) using next-generation sequencing (NGS). The increase in the number of viruses is associated with viruses that are typically found in terrestrial mammals, being found in marine mammals as well. For instance, Mäkelä et al. (2015) found respiratory diseases in marine mammals caused by influenza A and B virus infections (Mäkelä et al., 2015). Fouchier et al. (2004) found that H7N7 infected seals, which then infected human beings causing conjunctivitis (Fouchier et al., 2004). In addition, the US Centers for Disease Control and Prevention revealed that from 2007 to 2009, various belugas in North America were infected with the seal pox virus, which caused various diseases (Roess et al., 2011) and the individuals who were in close contact with belugas also tested positive for the virus. A study at the University of Florida found that the San Miguel sea lion virus (SMSV) was indistinguishable from the vesicular exanthema of swine virus (VESV), which caused vesicular inflammation in California pigs (Waltzek et al., 2012). These studies suggest the possibility of zoonotic viruses being carried by marine mammals and the risk of cross-species transmission. Cui et al. (2019) and Zhou et al. (2020) found that severe acute respiratory syndrome coronavirus (SARS-CoV) and Middle East respiratory syndrome coronavirus (MERS-CoV), two highly infectious and pathogenic viruses, might have originated from bats. The genetic recombination ability of coronaviruses makes them more adaptable to their hosts and increases their survival (Su et al., 2016). Therefore, the growing global virus diversity and potential spillovers require attention (Cui et al., 2019).

Papillomaviruses (PVs) are a diverse group of small, non-enveloped viruses belonging to the family Papillomaviridae, with double-stranded DNA genomes ranging in length from 5 to 8 kb (Zhang et al., 2022). The structurally conserved circular genome of PVs is organized into several partially overlapping open reading frames (ORFs) (Arman and Munger, 2022). These include the major genes present in all PVs: the early genes (E1 and E2), which are responsible for replication (Kraberger et al., 2022; Prabhakar et al., 2023), regulation of the life cycle, viral gene expression. The late genes (L1 and L2) are controlling the formation of viral capsids. Many PVs also carry other early genes, viral oncogenes E6 and E7 that are responsible for their ability to interfere with tumor suppressor genes such as p53 and the retinoblastoma protein (pRB) family (Nadile et al., 2022; Sofiani et al., 2023; Zhi et al., 2023). The upstream regulatory region or long control region (URR or LCR) located between the L1 and E6 ORFs is responsible for the origin of virus replication, as well as virus and cell transcription factor binding sites (TFBS) (Warburton et al., 2018; Coursey et al., 2021).

While the understanding of non-human PV types has increased in recent years, knowledge about the diversity of animal PVs remains limited (Bernard et al., 2010). The discovery of bovine and human papillomavirus sequences in cats suggests the possibility of cross-species transmission (Egberink et al., 2013). Recently, research on invertebrates and vertebrates in marine ecosystems has focused on the origin, evolution, and emergence of terrestrial viruses. Papillomaviruses were first reported in porpoises (Phocoena phocoena) (Gottschling et al., 2011) and subsequently found in many other small Odontoceti, including common dolphins (Delphinus delphis) (Robles-Sikisaka et al., 2012), common bottlenose dolphins (Tursiops truncatus) (Bossart et al., 2019), and burmester’s porpoises (Phocoena spinipinnis) (Van Bressem et al., 2007; Brimer et al., 2017). In order to expand our understanding of the diversity of marine animal PVs, we conducted a marine mammal virus survey program and discovered two novel papillomaviruses (PV-HMU-1 and PV-HMU-2) from nose, throat, and anal swab samples of belugas (Delphinapterus leucas). We obtained and characterized the complete genome sequences of PV-HMU-1 and PV-HMU-2. The alignment results of the novel PVs revealed low amino acid identities with other known papillomaviruses and clustered with dolphin papillomaviruses in the evolutionary tree. The findings of this study deepen our understanding of the host and genetic diversity of papillomaviruses and provide essential data for the possible zoonotic transmission caused by animal-origin PVs.

2. Materials and methods

2.1. Swab sample collection

We collected a total of 21 nose, throat, and anal swabs from ten belugas and one false killer whale (Pseudorca crassidens), housed in aquaria at Qingdao Polar Haichang Ocean Park and Dalian Sun Asia Tourism Holding Co., China in April 2018 (Supplementary Table 1). The nine throat, nose, and anal swabs from mammals in Qingdao constituted pool 49 (Table 1), while the 12 throat and anal swabs from Dalian constituted pool 41 (Table 2). The sampling procedures were approved by the Ethics Committee of the Hainan Medical University (Approval number: HMUEC20180059). The collected samples were quickly immersed in the maintenance medium in virus-sampling tubes (Yocon Biology, Beijing, China), to ensure sample quality, and transported to the laboratory within 24 h using an ice-cold dry chain (Wu et al., 2018a). The samples were then stored at −80°C (Wu et al., 2018a).

TABLE 1

Table 1. Detection rate of PV-HMU-1 among belugas from Qingdao Ocean Park.

TABLE 2

Table 2. Detection rate of PV-HMU-2 among belugas from Dalian Ocean Park.

2.2. Viral nucleic acid library construction

Samples from each species were combined by adding 1 mL from each sample into one fresh sample tube. Based on the swab sample locations, the 21 samples were divided into two pools and passed through 0.45 μm filters (Millipore Sigma, Burlington, MA, USA) to remove eukaryotic cells and large bacteria. The filtrate was then ultracentrifuged at 100,000 × g at 4°C for 3 h. The precipitates collected from the two pool samples were resuspended in 100 μL Hank’s balanced salt solution and digested using DNase, which consisted of 14 U of TURBO DNase (Ambion, USA), 20 U of benzonase (Novagen, Germany), and 20 U of RNase One (Promega, USA) at 37°C for 2 h in 1 × DNase buffer to decompose and remove unprotected nucleic acids. Viral DNA and RNA were extracted using QIAamp Viral RNA Mini Kit (Qiagen, Hilden, Germany), according to the manufacturer’s instructions. First-strand cDNA was generated using Superscript III Reverse Transcriptase (Invitrogen, Thermo Fisher Scientific), as previously described (Wu et al., 2018b). Double-stranded DNA was synthesized using Klenow fragment enzyme combined with K-8N (GACCATCTAGCGACCTCCACNNNNNNNN). To eliminate excess enzymes and PCR product primers, a PCR purification kit was used, followed by the use of magnetic beads to adsorb 300–2,000 bp fragments and dissolve them in nuclease free water. The nucleic acid library was then constructed using the DNA Library Prep Kit (Invitrogen Collibri, USA).

2.3. Next-generation sequencing

The amplified viral nucleic acid libraries were sequenced with an Illumina HiSeq 2500 sequencer using the 150 bp paired-end method (Illumina Inc., San Diego, CA, USA). To maximize the available length and total output of the raw data, each pool sample library was sequenced with separate channels and indexing. The raw sequence reads were filtered, and valid sequences were obtained using previously described criteria. The initial image analysis and base calling were performed using the GAPipeline program with default parameters. Sequences that include the primer K and adaptor are extracted using the Illumina filter (Wu et al., 2012). The sequence data were deposited into the National Center for Biotechnology Information sequence read archive under the accession number PRJNA650224.

2.4. Taxonomic assignment

Raw sequence reads were filtered to obtain valid sequences using previously described criteria (Yang et al., 2011). Assembly of quality-controlled reads of each sample into contigs was performed with Trinity (Version 2.5.1). Diamond (Version 0.9.14.115) software was used to translate contigs into amino acid sequences. Sequence similarity-based taxonomic assignments were conducted as described previously (Yang et al., 2011). Briefly, each sequence was aligned with the NCBI non-redundant protein databases (NR) using the basic local alignment search tool (Blast)¹ Blast X (E-value < 10^–5, −F: filter query sequence, default = T) to determine its viral origin. All Blast X results of contigs were annotated with taxonomy IDs MegaHit and other information obtained from the NCBI Entrez server (Tamura et al., 2013).

2.5. Genome sequencing of papillomavirus

Molecular clues obtained from metagenomic analyses were utilized to classify sequence reads into viral families or genera using MEGAN.² Representative reads of novel papillomaviruses were selected for PCR amplification and sequencing verification, and the complete novel PVs sequences were amplified. To amplify and sequence parts of the genome, specific nested PCR primers were designed. cDNA was generated using Superscript III Reverse Transcriptase (Invitrogen/United States) and random primers. All primer sequences were based on newly obtained reads and amplified sequences, and the primer sequences used are shown in Supplementary Table 2.

2.6. Genome annotation

The viral nucleotide sequences of the genomes and the amino acid sequences of the ORFs were deduced by comparing them against those of other viral families. Conserved protein families and domains were predicted using Pfam,³ Blast p (see text footnote 1), and InterProScan 5 (see text footnote 3). The genome structure diagram was constructed using IBS 1.0.3 software, and routine sequence alignment was performed using Clustal Omega.⁴ TFBS was predicted with “BDGP”⁵ and “Softberry.”⁶

2.7. Calculation of viral prevalence

Specific primers targeting the non-structural gene were designed based on the partial viral genomic sequences obtained from NGS for PV screening in each individual sample (Supplementary Table 2). PCR was performed using GoTaq Colorless Master Mix (Promega, Madison, WI, USA). The first-round PCR product (2 μl) as the template for the second round of PCR was used to achieve high specificity and sensitivity. The thermal cycling conditions for PCRs were as follows: 94°C for 5 min, followed by 35 cycles of 94°C for 30 s, 57°C for 35 s, 72°C for 30 s, and a final elongation step at 72°C for 10 min. The PCR products were analyzed using 1.5% agarose gel electrophoresis and ultraviolet imaging (Bio-Rad, USA).

2.8. Phylogenetic analysis

An amino acid alignment was performed between the novel beluga PV and other reference sequences using the Blast tool on the NCBI to identify closely related sequences. The phylogenetic tree included both mammalian and non-mammalian papillomaviruses, and representative species from each genus were selected as reference sequences. According to a previous study, each PV reference sequence is concatenated independently (Agius et al., 2019). A total of 61 PVs reference sequences were combined with the novel PVs sequences in this study to construct a phylogenetic analysis tree. MEGAX was used to align the amino acid sequences using the MUSCLE package with default parameters. The phylogenetic tree was constructed using the maximum likelihood method. The best substitution models were the Poisson model with uniform rates among sites, and this model was replicated with 1,000 bootstrap replicates to construct the phylogenetic tree.

3. Results

3.1. Discovery of two novel PVs in belugas using NGS

In this study, we obtained 2.74 GB (214,466 valid reads, 150 bp in length) and 1.39 GB (102,834 valid reads, 150 bp in length) of nucleotide data from sample pools 41 and 49, respectively. Sequences of less than 50 bp and duplicate sequences were eliminated. The archaea, bacteria, microbial eukaryotes, and those that had no significant similarity with any amino acid sequences in the viral NR protein database were eliminated. In pool 49, 14,055 virus-associated reads, accounting for 6.55% of the total reads, were identified, while in pool sample 41, 2,404 virus-associated reads, accounting for 2.33% of the total reads were identified. The viromes of both pools consisted of double-stranded (ds) DNA, single-stranded (ss) DNA, ssRNA, and retro-transcribing viruses. The virus-associated reads of two pool samples were assigned to families such as Picornaviridae, Herpesviridae, Myoviridae, Podoviridae, Siphoviridae, Microviridae, Retroviridae, and unclassified viruses (Figure 1). Moreover, swab sample pools from Qingdao and Dalian contained 418 and 54 papillomavirus reads, respectively, with PVs accounting for 2.97 and 2.24% of the total viruses in pools 49 and 41, respectively. After PCR verification, two novel PV sequences were obtained and compared with the original sequencing PV reads, which showed that the PV reads in the metagenomic sequencing results belonged to the two novel PVs. In addition, molecular clues related to mammalian Picornaviridae and Herpesviridae were found in the whale swab samples from Dalian, with 491 reads related to herpesviruses showing less than 50.68% aa identity with other known herpesviruses and being related to Beluga alpha herpes virus 1. Furthermore, the sequencing results showed 56 reads related to picornavirus, with less than 40.24% aa identity with other known picornaviruses, and verification of this virus was completed by our laboratory (Wang et al., 2021).

FIGURE 1

Figure 1. Viral reads from two sample pools using next-generation sequencing (NGS). (A) Identification of viral groups in swab samples from Qingdao belugas and false killer whale. (B) Viral content in swab samples from Dalian belugas.

3.2. Genomic characterization of novel PVs

Specific primers were designed based on the NGS-generated fragment sequences to obtain the complete genome sequences of two novel PVs from samples obtained from belugas in Qingdao and Dalian Ocean Parks. The novel PVs were named PV-HMU-1 (GenBank accession: OP748391) and PV-HMU-2 (GenBank accession: OP748392), respectively. The complete genome length of PV-HMU-1 and PV-HMU-2 were 7,624 and 7,179 bp, respectively, with C + G contents of 58.77 and 54.37% (Figure 2 and Supplementary Table 1). The Pfam software predicted that the observed mammalian PV genomic structures included a core set of early (E6, E1, and E2) and late (L2 and L1) open-reading frames (ORFs) (Figure 2). ORF prediction software revealed that PV-HMU-1 E6 (184 aa) was 552 nt (1–552), E1 (625 aa) was 1,875 nt (687–2,561), early E2 (393 aa) was 1,179 nt (2,506–3,684), late L2 (535 aa) was 1,605 nt (3,974–5,578), and late L1 (582 aa) was 1,746 nt (5,337–7,082) in length. For PV-HMU-2, E6 (233 aa) was 699 nt (1–699), E1 (622 aa) was 1,866 nt (755–2,620), early E2 (372 aa) was 1,116 nt (2,556–3,671), late L2 (450 aa) was 1,350 nt (3,787–5,136), and late L1 (496 aa) was 1,488 nt (5,120–6,607) in length.

FIGURE 2

Figure 2. Sequence and motif annotation of PV-HMU-1 and PV-HMU-2. (A) Putative and confirmed early and late ORFs identified in PV-HMU-1 and PV-HMU-2 genomes using ORF predictor. We used IBS 1.0.3 software to construct the genomic structure diagram of PV-HMU-1, PV-HMU-2, Phocoena phocoena papillomavirus 4 (PphPV4) and Lagenorhynchus acutus papillomavirus (GU117624). (B) The E6 conserved protein domains of PV-HMU-1 and PV-HMU2 (light blue shade) was aligned with other papillomaviruses. The C-terminal zinc-binding domains for PV-HMU-1, PV-HMU-2, HPV16, BPV1, TtPV2, and PphPV4 are indicated.

A study by Schulz et al. (2012) used TFBS prediction to identify potential cis-regulatory elements in the upstream regulatory region (URR) of rodent PVs. The URR sequences of PV-HMU-1 contained important TFBS, such as activator protein 1 (AP1), which were located directly upstream of the TFIID BS (TATA box binding protein). However, the TFBS of PV-HMU-2 was not clearly predicted. Both PV-HMU-1 and PV-HMU-2 included a zinc-finger motif (C-X2-C-X29-C-X2-C) in E6, while E7 was absent (Figure 2). An adenosine triphosphate/guanosine triphosphate (ATP/GTP) binding site (G-X4-GKS) was present within the helicase domain of E1 in both sequences. The N-terminal, DNA-binding domain, and helicase were located at positions 689–1,058 bp| 755–1,124 bp, 1,235–1,652 bp| 1,273–1,700 bp, and 1,707–2,504 bp| 1,712–2,589 bp of the E1 region for PV-HMU-1 and PV-HMU-2, respectively. Comparing the L2 and L1 amino acid sequences of PV-HMU-1 and PV-HMU-2 with those of other papillomaviruses, PV-HMU-1 showed 41.97 and 64.99% identity to Phocoena papillomavirus 4 (NC 018076), respectively, while PV-HMU-2 showed 49.19 and 68.12% identity to Sus scrofa papillomavirus 1 (MT078972), respectively.

3.3. Phylogenetic analyses

Phylogenetic analyses were conducted using the complete E1-E2-L2-L1 amino acid sequences of PV-HMU-1 and PV-HMU-2, as well as 61 representative references related to Papillomaviridae in GenBank (Figure 3). Based on the results of Blast and preliminary evolutionary analysis, the novel PVs were attributed to the mammalian PV clade, so we mainly selected mammalian PVs as references, with a few Sauropsid and Fish PV clades selected as out-groups. More than 50 genera of PVs existed, and we selected the representative reference of each genus and the virus strains with the highest identities relative to the novel PVs to construct the evolutionary tree. The phylogenetic analysis showed that the amino acid combination of PV-HMU-1 clustered with Tursiops truncatus papillomavirus 2 (TtPV2) (NC 008184) and formed a relatively independent clade within the genus Upsilonpapillomavirus. PV-HMU-2 clustered with Phocoena phocoena papillomavirus 4 (PphPV4) (GU117623), the only species of Dyopipapillomavirus. Sus scrofa papillomavirus 1 is a member of the genus Dyoomikronpapillomavirus and was clustered with PV-HMU-2. PV-HMU-1 and PV-HMU-2 did not cluster on the same branch on the evolutionary tree and were far from human-associated papillomaviruses. The genome sequences of PV-HMU-1 and PV-HMU-2 carried by belugas were difficult to compare with those of other known papillomaviruses at the nucleic acid level in the NCBI nucleotide database. Genome identity and phylogenetic analyses suggested that these beluga PVs represented novel species.

FIGURE 3

Figure 3. Phylogenetic tree of concatenated early (E1, E2) and late (L1) papillomavirus proteins. The evolutionary history was inferred from 61 amino acid sequences, including the novel PV-HMU-1 and PV-HMU-2, using a maximum likelihood approach. A 1,000 bootstrap replication were used to construct the tree. PV-HMU-1 and PV-HMU-2 are marked with red circles under mammalian-related PVs. The virus strains close to the novel identified PVs are marked with pink background, Sauropsid PV, and Fish PV clades are marked with green and blue backgrounds, respectively.

3.4. Prevalence of PV infection in belugas

Gene-specific nested primers targeting the L1 gene region for PCR were designed based on the genomic sequences of novel PVs to screen for their presence in the throat, nose, and anal swabs. The throat and nose swabs of belugas in both locations showed a high infection rate of PVs, while the infection rate in anal swabs was low. Thirty-three percent (1/3) of throat swabs and nose swabs, but no anal swabs were positive for PV-HMU-1 (Table 1). Seventy-one percent (5/7) of throat swabs but no anal swabs were positive for PV-HMU-2 (Table 2). All samples collected in Qingdao and Dalian were collected three times in April, September, and December 2018 to observe the virus-carrying rate in different quarters. The throat swabs of three animals in Dalian consistently tested positive for PV-HMU-2 at all three sampling times, indicating persistent viral presence. To confirm the presence of the same virus infection in different months, a set of L primers was used to screen the samples. After amplification, fragments of appropriate size were obtained and sent for sequencing. The sequencing results of two positive samples collected at different times showed an aa identity of > 99%. The positive samples were identified by amplification products of 540 bp for PV-HMU-1 and 522 bp for PV-HMU-2.

4. Discussion

The increasing contact between humans and marine mammals due to the exploration of the oceans and interest in these animals has led to a higher risk of zoonotic diseases that can endanger human health. While efforts have been made to characterize novel PV types in different animal hosts, knowledge about PV diversity in marine mammals remains limited. In recent years, many viruses, including picornaviruses, coronaviruses, morbilliviruses, and poxviruses, have been detected in dolphins (Robles-Sikisaka et al., 2012), seals (Fouchier et al., 2004), and sea lions (Waltzek et al., 2012). Due to the continuous expansion of the host range of viruses in marine mammals, it is suggested that they may serve as important hosts for many viruses, and there may still be many unknown viruses in these hosts. Reports of marine mammals acting as reservoirs for viruses are increasing. In Groth et al. (2014) obtained a complete sequence of influenza A virus (H13N2) from a sick pilot whale, while Sacristán et al. (2019) and Vargas-Castro et al. (2021) found that dolphins in the South Atlantic and whales in Spain were highly infected with herpesvirus (Sacristán et al., 2019; Vargas-Castro et al., 2021). Beluga coronavirus was detected in a study at the Washington University School of Medicine (Mihindukulasuriya et al., 2008), and a novel bottlenose dolphin coronavirus was detected at the University of Hong Kong (Woo et al., 2014). Changes in living habits, social behaviors, and migration patterns are significant factors in the transmission and spread of viruses among cetaceans. Monitoring the viruses carried by marine mammals is crucial for preventing potential infectious diseases of cetacean origin. In this study, we detected PV-HMU-1 and PV-HMU-2 in belugas from Qingdao and the Dalian Polar Ocean Parks and obtained their complete genomes. The genomes of PV-HMU-1 and PV-HMU-2 contain five ORFs (E6, E1, E2, L2, and L1) of mammal PVs but do not have the E7 ORF typical of Alphapapillomavirus (Figure 2). The E6 protein of the novel virus contains a zinc-finger motif that is crucial for the virus to invade the cell environment. AP1 is present at the PV-HMU-1 TFBS site and is used to regulate α-HPV expression. The two novel PVs contain an ATP/GTP binding site, an N-terminal, DNA binding region, and a helix in the E1 region. Papillomaviruses are generally transmitted through sex, close contact, and mother-to-child transmission. The 50% infection rate of PV-HMU-2 infection in individual samples reflected the possibility of a high incidence of PV infection in belugas. Belugas that tested positive for PV-HMU-1 and PV-HMU-2 were positive three times within a year, indicating that the infection was persistent. It is important to note that some of the negative results may be due to inadequate sampling.

In 2019, the ICTV report on the classification of genera highlighted that, based on multiple sequences or comparative pairs of L1 genes, most PV genus types have less than 60% sequence identity with other types, but this is not the only criterion (de Villiers et al., 2004). The L gene aa identity between PV-HMU-1 and PV-HMU-2 was 65.21%, although they were found in the same beluga host. The two novel PVs showed less than 68.12% aa identity with the L gene region of known papillomaviruses. The new viruses mainly cluster with viruses whose hosts are small toothed whales in the phylogenetic tree, which are closely related to Upsilonpapillomavirus and Dyopipapillomavirus. The clustering of PV-HMU-2 and pig PV suggests that the evolution of PV in belugas and terrestrial animals is related.

Several phylogenetic studies have shown that a complex evolutionary scenario drives PV diversity, such as recombination, co-divergence with the host, interspecies transmission, and adaptation to new environments (Geoghegan et al., 2017). Although the belugas in this study did not exhibit obvious tumors or other clinical manifestations, the high infection rate suggests a possible association between papillomavirus infection and belugas. While this study amplified the complete genome sequences of two PVs and investigated their classification, it is merely a small step toward understanding the complexities of papillomaviruses in marine mammal populations. More work is necessary to examine the evolution, transmission, and biological role of papillomaviruses in marine mammal populations, isolate and record them, and provide fundamental data for future research on the evolution, classification, and traceability of mammalian viruses.

In this study, we reported the discovery of two novel papillomaviruses in the throat and nasal swabs collected from belugas in Qingdao and Dalian. To our knowledge, belugas have never been reported as hosts for papillomaviruses before. These findings increase our understanding of the mammal-borne virus community in Qingdao and Dalian and provide preliminary data for future research on the evolution of PVs. However, further investigation is necessary to determine the pathogenicity of these two novel viruses and their potential impact on both humans and animals.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in this article/Supplementary material.

Author contributions

FY, JD, and PZ designed the study. PZ, SB, XS, and XL collected the specimens. YL, ZL, MX, YZ, GW, RP, HS, and XH performed the experiments. CT and GL analyzed the data. YL wrote the draft manuscript. JD and PZ edited the manuscript. All authors read and approved the final manuscript.

Funding

This work was supported by the National Natural Science Foundation of China (No. 32060015), the National Key R&D Program of China (grant No. 2022YFE0210300), the Key Research and Development Plan of Hainan Province (Nos. ZDYF2020150 and ZDYF2022SHFZ085), the CAMS Innovation Fund for Medical Sciences (CIFMS) (grant No. 2021-I2M-1-039), the Major Science and Technology Program of Hainan Province (No. ZDKJ202003), the Hainan Provincial Natural Science Foundation of China (Nos. 820RC628 and 823RC598), the Research Project of the Hainan Academician Innovation Platform (No. YSPTZX202004), the Hainan Talent Development Project (No. SRC200003), and the Graduate student innovation grant of Hainan Medical University (Nos. HYYS2020-33, Hys2020-360, Hys2020-376, and Qhys2021-356).

Conflict of interest

XS was employed by Dalian Sun Asia Tourism Holding Co., Ltd.

The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

The reviewer SS declared a past co-authorship with SB and PZ to the handling editor MH.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2023.1165839/full#supplementary-material

Footnotes

1. ^ https://blast.ncbi.nlm.nih.gov
2. ^ https://megasoftware.net/show_eua
3. ^ http://www.ebi.ac.uk/services/proteins
4. ^ http://www.ebi.ac.uk/tools/
5. ^ http://www.fruitfly.org/seq_tools/promoter.html
6. ^ http://www.softberry.com/

References

Agius, J. E., Phalen, D. N., Rose, K., and Eden, J. S. (2019). New insights into Sauropsid Papillomaviridae evolution and epizootiology: Discovery of two novel papillomaviruses in native and invasive Island geckos. Virus Evol. 5:2. doi: 10.1093/ve/vez051

PubMed Abstract | CrossRef Full Text | Google Scholar

Arman, W., and Munger, K. (2022). Mechanistic contributions of lncRNAs to cellular signaling pathways crucial to the lifecycle of human papillomaviruses. Viruses 14:2439. doi: 10.3390/v14112439

PubMed Abstract | CrossRef Full Text | Google Scholar

Bernard, H. U., Burk, R. D., Chen, Z., van Doorslaer, K., zur Hausen, H., and de Villiers, E. M. (2010). Classification of papillomaviruses (PVs) based on 189 PV types and proposal of taxonomic amendments. Virology 401, 70–79. doi: 10.1016/j.virol.2010.02.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Bossart, G. D., Romano, T. A., Peden-Adams, M. M., Schaefer, A. M., Rice, C. D., Fair, P. A., et al. (2019). Comparative innate and adaptive immune responses in Atlantic bottlenose dolphins (Tursiops truncatus) with viral, bacterial, and fungal infections. Front. Immunol. 10:1125. doi: 10.3389/fimmu.2019.01125

PubMed Abstract | CrossRef Full Text | Google Scholar

Brimer, N., Drews, C. M., and Vande Pol, S. B. (2017). Association of papillomavirus E6 proteins with either MAML1 or E6AP clusters E6 proteins by structure, function, and evolutionary relatedness. PLoS Pathog. 13:e1006781. doi: 10.1371/journal.ppat.1006781

PubMed Abstract | CrossRef Full Text | Google Scholar

Cheng, D., Qiu, K., Rao, Y., Mao, M., Li, L., Wang, Y., et al. (2023). Proliferative exhausted CD8(+) T cells exacerbate long-lasting anti-tumor effects in human papillomavirus-positive head and neck squamous cell carcinoma. Elife 12:e82705. doi: 10.7554/eLife.82705

PubMed Abstract | CrossRef Full Text | Google Scholar

Coursey, T. L., Van Doorslaer, K., and McBride, A. A. (2021). Regulation of Human Papillomavirus 18 genome replication, establishment, and persistence by sequences in the viral upstream regulatory region. J. Virol. 95:e0068621. doi: 10.1128/jvi.00686-21

PubMed Abstract | CrossRef Full Text | Google Scholar

Cui, J., Li, F., and Shi, Z. L. (2019). Origin and evolution of pathogenic coronaviruses. Nat. Rev. Microbiol. 17, 181–192. doi: 10.1038/s41579-018-0118-9

PubMed Abstract | CrossRef Full Text | Google Scholar

de Villiers, E. M., Fauquet, C., Broker, T. R., Bernard, H. U., and zur Hausen, H. (2004). Classification of papillomaviruses. Virology 324, 17–27. doi: 10.1016/j.virol.2004.03.033

PubMed Abstract | CrossRef Full Text | Google Scholar

Egberink, H., Thiry, E., Möstl, K., Addie, D., Belák, S., Boucraut-Baralon, C., et al. (2013). Feline viral papillomatosis: ABCD guidelines on prevention and management. J. Feline Med. Surg. 15, 560–562. doi: 10.1177/1098612x13489213

PubMed Abstract | CrossRef Full Text | Google Scholar

Fouchier, R. A., Schneeberger, P. M., Rozendaal, F. W., Broekman, J. M., Kemink, S. A., Munster, V., et al. (2004). Avian influenza A virus (H7N7) associated with human conjunctivitis and a fatal case of acute respiratory distress syndrome. Proc. Natl. Acad. Sci. U.S.A. 101, 1356–1361. doi: 10.1073/pnas.0308352100

PubMed Abstract | CrossRef Full Text | Google Scholar

Gaynor, A. M., Fish, S., Duerr, R. S., Cruz, F. N. Jr., and Pesavento, P. A. (2015). Identification of a novel papillomavirus in a Northern Fulmar (Fulmarus glacialis) with viral production in cartilage. Vet. Pathol. 52, 553–561. doi: 10.1177/0300985814542812

PubMed Abstract | CrossRef Full Text | Google Scholar

Geoghegan, J. L., Duchêne, S., and Holmes, E. C. (2017). Comparative analysis estimates the relative frequencies of co-divergence and cross-species transmission within viral families. PLoS Pathog. 13:e1006215. doi: 10.1371/journal.ppat.1006215

PubMed Abstract | CrossRef Full Text | Google Scholar

Gottschling, M., Bravo, I. G., Schulz, E., Bracho, M. A., Deaville, R., Jepson, P. D., et al. (2011). Modular organizations of novel cetacean papillomaviruses. Mol. Phylogenet. Evol. 59, 34–42. doi: 10.1016/j.ympev.2010.12.013

PubMed Abstract | CrossRef Full Text | Google Scholar

Groth, M., Lange, J., Kanrai, P., Pleschka, S., Scholtissek, C., Krumbholz, A., et al. (2014). The genome of an influenza virus from a pilot whale: Relation to influenza viruses of gulls and marine mammals. Infect. Genet. Evol. 24, 183–186. doi: 10.1016/j.meegid.2014.03.026

PubMed Abstract | CrossRef Full Text | Google Scholar

Kraberger, S., Austin, C., Farkas, K., Desvignes, T., Postlethwait, J. H., Fontenele, R. S., et al. (2022). Discovery of novel fish papillomaviruses: From the Antarctic to the commercial fish market. Virology 565, 65–72. doi: 10.1016/j.virol.2021.10.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Mäkelä, S. M., Österlund, P., Westenius, V., Latvala, S., Diamond, M. S., and Gale, M. Jr., et al. (2015). RIG-I signaling is essential for influenza B virus-induced rapid interferon gene expression. J. Virol. 89, 12014–12025. doi: 10.1128/jvi.01576-15

PubMed Abstract | CrossRef Full Text | Google Scholar

Mihindukulasuriya, K. A., Wu, G., St Leger, J., Nordhausen, R. W., and Wang, D. (2008). Identification of a novel coronavirus from a beluga whale by using a panviral microarray. J. Virol. 82, 5084–5088. doi: 10.1128/jvi.02722-07

PubMed Abstract | CrossRef Full Text | Google Scholar

Nadile, M., Retsidou, M. I., Gioti, K., Beloukas, A., and Tsiani, E. (2022). Resveratrol against cervical cancer: Evidence from in vitro and in vivo studies. Nutrients 14:5273. doi: 10.3390/nu14245273

PubMed Abstract | CrossRef Full Text | Google Scholar

Pino, M. D., Matas, I., Carrillo, P., Martí, C., Glickman, A., Carreras-Dieguez, N., et al. (2023). Natural history of anal HPV infection in women treated for cervical intraepithelial neoplasia. Cancers 15:1147. doi: 10.3390/cancers15041147

PubMed Abstract | CrossRef Full Text | Google Scholar

Prabhakar, A. T., James, C. D., Fontan, C. T., Otoa, R., Wang, X., Bristol, M. L., et al. (2023). Human papillomavirus 16 E2 interaction with TopBP1 is required for E2 and viral genome stability during the viral life cycle. J. Virol. 97:e0006323. doi: 10.1128/jvi.00063-23

PubMed Abstract | CrossRef Full Text | Google Scholar

Robles-Sikisaka, R., Rivera, R., Nollens, H. H., St Leger, J., Durden, W. N., Stolen, M., et al. (2012). Evidence of recombination and positive selection in cetacean papillomaviruses. Virology 427, 189–197. doi: 10.1016/j.virol.2012.01.039

PubMed Abstract | CrossRef Full Text | Google Scholar

Roess, A. A., Levine, R. S., Barth, L., Monroe, B. P., Carroll, D. S., Damon, I. K., et al. (2011). Sealpox virus in marine mammal rehabilitation facilities, North America, 2007-2009. Emerg. Infect. Dis. 17, 2203–2208. doi: 10.3201/eid1712.101945

PubMed Abstract | CrossRef Full Text | Google Scholar

Sacristán, C., Esperón, F., Ewbank, A. C., Díaz-Delgado, J., Ferreira-Machado, E., Costa-Silva, S., et al. (2019). Novel herpesviruses in riverine and marine cetaceans from South America. Acta Trop. 190, 220–227. doi: 10.1016/j.actatropica.2018.11.021

PubMed Abstract | CrossRef Full Text | Google Scholar

Schulz, E., Gottschling, M., Ulrich, R. G., Richter, D., Stockfleth, E., and Nindl, I. (2012). Isolation of three novel rat and mouse papillomaviruses and their genomic characterization. PLoS One 7:e47164. doi: 10.1371/journal.pone.0047164

PubMed Abstract | CrossRef Full Text | Google Scholar

Seyoum, A., Seyoum, B., Gure, T., Alemu, A., Belachew, A., Abeje, D., et al. (2023). Genotype heterogeneity of high-risk human papillomavirus infection in Ethiopia. Front. Microbiol. 14:1116685. doi: 10.3389/fmicb.2023.1116685

PubMed Abstract | CrossRef Full Text | Google Scholar

Sofiani, V. H., Veisi, P., Rukerd, M. R. Z., Ghazi, R., and Nakhaie, M. (2023). The complexity of human papilloma virus in cancers: A narrative review. Infect. Agent Cancer 18:13. doi: 10.1186/s13027-023-00488-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Su, S., Wong, G., Shi, W., Liu, J., Lai, A. C. K., Zhou, J., et al. (2016). Epidemiology, genetic recombination, and pathogenesis of coronaviruses. Trends Microbiol. 24, 490–502. doi: 10.1016/j.tim.2016.03.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Tamura, K., Stecher, G., Peterson, D., Filipski, A., and Kumar, S. (2013). MEGA6: Molecular evolutionary genetics analysis version 6.0. Mol. Biol. Evol. 30, 2725–2729. doi: 10.1093/molbev/mst197

PubMed Abstract | CrossRef Full Text | Google Scholar

Van Bressem, M. F., Cassonnet, P., Rector, A., Desaintes, C., Van Waerebeek, K., Alfaro-Shigueto, J., et al. (2007). Genital warts in Burmeister’s porpoises: Characterization of Phocoena spinipinnis papillomavirus type 1 (PsPV-1) and evidence for a second, distantly related PsPV. J. Gen. Virol. 88, 1928–1933. doi: 10.1099/vir.0.82694-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Van Doorslaer, K., Ruoppolo, V., Schmidt, A., Lescroël, A., Jongsomjit, D., Elrod, M., et al. (2017). Unique genome organization of non-mammalian papillomaviruses provides insights into the evolution of viral early proteins. Virus Evol. 3:vex027. doi: 10.1093/ve/vex027

PubMed Abstract | CrossRef Full Text | Google Scholar

Vargas-Castro, I., Melero, M., Crespo-Picazo, J. L., Jiménez, M., Sierra, E., Rubio-Guerri, C., et al. (2021). Systematic determination of herpesvirus in free-ranging cetaceans stranded in the Western Mediterranean: Tissue tropism and associated lesions. Viruses 13:2180. doi: 10.3390/v13112180

PubMed Abstract | CrossRef Full Text | Google Scholar

Waltzek, T. B., Cortés-Hinojosa, G., Wellehan, J. F. Jr., and Gray, G. C. (2012). Marine mammal zoonoses: A review of disease manifestations. Zoonoses Public Health 59, 521–535. doi: 10.1111/j.1863-2378.2012.01492.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, G., Huang, Y., Zhang, W., Peng, R., Luo, J., Liu, S., et al. (2021). Identification and genome analysis of a novel picornavirus from captive belugas (Delphinapterus leucas) in China. Sci. Rep. 11:21018. doi: 10.1038/s41598-021-00605-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Warburton, A., Redmond, C. J., Dooley, K. E., Fu, H., Gillison, M. L., Akagi, K., et al. (2018). HPV integration hijacks and multimerizes a cellular enhancer to generate a viral-cellular super-enhancer that drives high viral oncogene expression. PLoS Genet. 14:e1007179. doi: 10.1371/journal.pgen.1007179

PubMed Abstract | CrossRef Full Text | Google Scholar

Woo, P. C., Lau, S. K., Lam, C. S., Tsang, A. K., Hui, S. W., Fan, R. Y., et al. (2014). Discovery of a novel bottlenose dolphin coronavirus reveals a distinct species of marine mammal coronavirus in Gammacoronavirus. J. Virol. 88, 1318–1331.

Google Scholar

Wu, Z., Du, J., Lu, L., Yang, L., Dong, J., Sun, L., et al. (2018a). Detection of hantaviruses and arenaviruzses in three-toed jerboas from the inner mongolia autonomous region, China. Emerg. Microbes Infect. 7:35. doi: 10.1038/s41426-018-0036-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, Z., Liu, B., Du, J., Zhang, J., Lu, L., Zhu, G., et al. (2018b). Discovery of diverse rodent and bat pestiviruses with distinct genomic and phylogenetic characteristics in several Chinese provinces. Front. Microbiol. 9:2562. doi: 10.3389/fmicb.2018.02562

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, Z., Ren, X., Yang, L., Hu, Y., Yang, J., He, G., et al. (2012). Virome analysis for identification of novel mammalian viruses in bat species from Chinese provinces. J. Virol. 86, 10999–11012. doi: 10.1128/jvi.01394-12

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, J., Yang, F., Ren, L., Xiong, Z., Wu, Z., Dong, J., et al. (2011). Unbiased parallel detection of viral pathogens in clinical samples by use of a metagenomic approach. J. Clin. Microbiol. 49, 3463–3469. doi: 10.1128/jcm.00273-11

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, X., Li, Y., Tang, Y., Li, Z., Wang, S., Luo, X., et al. (2023). Cervical HPV infection in Guangzhou, China: An epidemiological study of 198,111 women from 2015 to 2021. Emerg. Microbes Infect. 12:e2176009. doi: 10.1080/22221751.2023.2176009

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, M., Kiyono, T., Aoki, K., Goshima, N., Kobayashi, S., Hiranuma, K., et al. (2022). Development of an in vitro carcinogenesis model of human papillomavirus-induced cervical adenocarcinoma. Cancer Sci. 113, 904–915. doi: 10.1111/cas.15246

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhi, W., Wei, Y., Lazare, C., Meng, Y., Wu, P., Gao, P., et al. (2023). HPV-CCDC106 integration promotes cervical cancer progression by facilitating the high expression of CCDC106 after HPV E6 splicing. J. Med. Virol. 95:e28009. doi: 10.1002/jmv.28009

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, P., Yang, X. L., Wang, X. G., Hu, B., Zhang, L., Zhang, W., et al. (2020). A pneumonia outbreak associated with a new coronavirus of probable bat origin. Nature 579, 270–273. doi: 10.1038/s41586-020-2012-7

PubMed Abstract | CrossRef Full Text | Google Scholar