Molecular and Histological Identification of Bovine Papillomavirus 1, 2 and a Novel Genotype in Cutaneous Papillomas of Dairy Cattle in Taiwan

TL;DR

This study identifies BPV1, BPV2, and a new BPV type in skin tumors of dairy cattle in Taiwan, revealing broader tissue infection potential.

Contribution

Discovery of a novel BPV genotype with expanded tissue tropism in dairy cattle in Taiwan.

Findings

BPV1 and BPV2 were confirmed in cutaneous papillomas from dairy cattle in Taiwan.

A novel BPV genotype with 80.9% sequence identity to BPV38 was identified.

The novel BPV infected both epithelial and mesenchymal cells, indicating broader tissue tropism.

Abstract

Bovine papillomaviruses (BPVs) are host-specific and strongly epitheliotropic infectious agents that cause benign epithelial and mucosal proliferations, with potential for malignant transformation. However, BPV1, BPV2, and BPV5 are unique in their ability to infect both epithelial and connective tissues. While BPV infections had been documented globally, there was no disease information reported from Taiwan. To investigate whether BPVs are associated with the development of cutaneous papillomas in dairy cattle in Taiwan, in the present study, eight cutaneous papilloma samples from six dairy farms were collected and analyzed by using histopathology, immunohistochemical (IHC) staining, and molecular biology methods. BPV1 and BPV2 were identified, along with a novel BPV sharing 80.9% sequence identity with BPV38. This novel BPV, classified under Xipapillomavirus, was detected in both epithelial and mesenchymal cells through in situ hybridization (ISH), suggesting a broader tissue tropism than typical Xipapillomavirus infections. These findings provide new insights into BPV diversity and pathogenesis.

Figures

Full text

1. Introduction

Papillomavirus (PV) is an epitheliotrophic pathogen that causes hyperplastic and neoplastic lesions in skin and mucosa [1, 2]. PVs are widely distributed among vertebrates, infecting a broad spectrum of species in mammal, avian, reptile, and Osteichthyes. However, most PV types exhibit strict host-specificity, indicating that cross-species transmission is rare [3]. Although the PV-associated lesions are usually benign and may spontaneously regress along with the elevation of the host's immunity [4, 5], PVs can also lead to malignancies [6–9]. Cervical cancer [6], squamous cell carcinoma (SCC) [7], tricoblastoma [7], Bowenoid in situ carcinoma (BISC), and transitional cell carcinoma (TCC) [8, 9] are linked with specific types of PVs in human, canine, feline, and bovine. Each PV harbors a circular, double-stranded genomic DNA sized approximately 5.7–8.6 kilobase pair (kbp), which encodes seven to nine open reading frames (ORFs) and a noncoding long control region (LCR) [4]. These ORFs encode two capsid proteins (L1 and L2) and four to seven nonstructural proteins (E1, E2, E4, E5, E6, E7, and E8) [4].

According to the definitions provided by the International Committee on Taxonomy of Viruses (ICTV) of PVs, with the sequence identity of L1 ORF <45% could be considered as a different subfamily; <60% a different genus; <90% a different genotype; <98% a different subtype; 98%–99% a different variant [10, 11]. To date, there are 44 genotypes of BPV mainly classified into five genera: Delta-, Xi-, Epsilon-, Dyoxi-, and Dyokappa-PV [12]. While most PVs exhibit strict host specificity and limited tissue tropism [3], BPV1, 2, 5, and 13 have broader host ranges, which have been identified in bovine, horse, giraffe, and donkey [13, 14]. BPV is unique not only in the host range but also in tissue tropism. BPV1, 2, and 5 can infect and transform subepithelial mesenchymal cells, and the subsequent proliferation of both epithelial and subcutaneous connective tissues ultimately leads to the formation of a fibropapilloma [14]. Most BPVs cause benign hyperplastic lesions, however, BPV1 and BPV2 are linked to malignant urinary bladder cancers in cattle [8, 9] and invasive sarcoids in equine [13], while BPV2, 4, 13, and 44 have been associated with alimentary SCCs [15, 16].

It is evidenced that certain BPVs could present not only in the exfoliated keratins but also in body fluids [17]. BPV-induced lesions typically appear on the hairy skin, teats, and urogenital mucosa [14]. Lesions surround oral cavity may cause dysphagia and anorexia, while those on the teats and urogenital mucosa can reduce milk yield, increase mastitis risk, and impair reproduction [18]. Consequently, BPV infections may result in economic losses at the farm level. Moreover, although the PVs are generally host-restricted, specific BPV types with a broad host range could be an invisible threat to both farm animals and wildlife.

Given the importance of BPV in both domestic and wild animals, BPV infections had been reported across several Asian countries [19–21], but relevant reports from Taiwan are limited. In this study, skin lesions clinically diagnosed as papillomas in lactating cattle were collected for BPV detection. Histopathology and polymerase chain reaction (PCR) using two degenerate primer sets were performed, followed by full-genome sequencing and phylogenetic analysis. Immunohistochemical (IHC) and in situ hybridization (ISH) staining further confirmed the presence of BPV in tissue sections.

2. Material and Methods

2.1. Sample Collection

Eight samples (Case nos. 23-00 to 23-07) were collected from six dairy farms in Taiwan, using excisional (7/8) or scraping (1/8) biopsy. All samples were protruding growths collected from adult dairy cattle. All the masses collected were located on the hairy skin of the flank and were clinically diagnosed as papillomas by veterinarians. They varied in size (diameter of the masses ranged 1–5 cm), but shared a consistent papillomatous appearance. Each tissue sample was divided in half, one fixed in neutral formalin for the histopathological examination while the other was stored at −20°C for further DNA extraction and molecular analysis. However, the sample from Case no. 23-04, collected via scraping biopsy, was submitted only for molecular detection. All procedures involving sample collection from the cattle were performed in accordance with the guidelines of the Institutional Animal Care and Use Committee of National Chung Hsing University.

2.2. DNA Extraction and Sequencing Analysis

The nucleic acids were extracted from the fresh tissues by using DNeasy Blood and Tissue Kit (Qiagen, Hilden, Germany). To detect the existence of BPV, two consensus degenerative primer pairs for PVs were used: MY09/11 (forward: 5′-GCMCAGGGWCATAAYAATGG-3′ and reverse: 5′-CGTCCMARRGGAWACTGATC-3′) [22] and CP4/5 (forward: 5′-ATGGTACARTGGGCATWTGA-3′ and reverse: 5′-GAGGYTGCAACCAAAAMTGRCT-3′) [23]. For each PCR reaction, 2 µL of DNA template, 10 µL of amaR OnePCR Buffer (GeneDireX, Taiwan), 0.1–0.2 µM of each primer, and PCR-grade water were mixed to the final volume of 20 µL. The PCR condition was listed as follows: 94°C for 3 min, 35 cycles of 94°C for 30 s, 45°C (MY09/11) or 50°C (CP4/5) for 1 min, 72°C for 45 s, and a 5-min final extension. Electrophoresis was conducted with 1.5% (w/v) agarose gel and the gel was subsequently stained with SYBR Safe DNA Gel Stain (Invitrogen, Waltham, MA, USA). The expected size of the amplicon of each primer pair was approximately 450 base pair (bp). All the amplicons on expected size were sent to Tri-I Biotech, Inc (Taipei, Taiwan) for Sanger sequencing bidirectionally to confirm the existence of PVs. BigDye Terminator v3.1 Cycle Sequencing Kit (Applied Biosystems, USA) was used for Sanger sequencing and the signals were detected by Applied Biosystems 3730xl DNA Analyzer (Applied Biosystems).

2.3. PCR for Full Genome Assembly of BPV1 and BPV2

Six out of eight (6/8) tissue samples tested positive for BPV1 and one (1/8) for BPV2, using conventional PCR technique as described above. For full genome genetic analysis of identified BPV1 and BPV2 isolates, three primer pairs (Table A1) were, respectively, designed to cover the full genome of BPV1 and BPV2 with overlapping regions. The Platinum Green Hot Start PCR Master Mix (Invitrogen) was utilized for long sequence amplification. Each reaction contained 25 µL of Master Mix, 0.2 µM of each primer, 2 µL of DNA template, and PCR-grade water for a total volume of 50 µL. The PCR condition included 94°C for 2 min, 35 cycles of 94°C for 30 s, 50°C for 30 s, and 72°C for 3 min, followed by a final extension of 5 min. The complete genomic assemblies of each case were uploaded to the NCBI database.

2.4. De Novo Next-Generation Sequencing (NGS)

One sample was confirmed with PV infection by conventional PCR technique with two degenerative primer sets. However, the preliminary results showed the PV detected in this case had a low sequence identity to the known PVs. Therefore, de novo NGS was used to characterize the full genome of the novel BPV. The DNA extracted from fresh tissue was sent to Tri-I Biotech, Inc (Taipei, Taiwan) for de novo NGS. The DNA sample was fragmented using the enzymatic fragmentation DNA-Seq kit (Takara Bio, CA, USA) and sequenced with the Illumina MiSeq System (Illumina, CA, USA). The quality of raw reads was assessed using NanoPlot v1.28.1. Genome assembly was performed with SPAdes v3.14. The complete genome assembly of this novel BPV was also uploaded to the NCBI database.

2.5. IHC Staining

To assess the presence of viral antigen in the lesions, IHC staining was performed. Serial deparaffinization was performed by using the Non-Xylene solution (Muto Chemical, Japan), followed by a rehydration process with graded ethanol from 99% to 80%. Slides were boiled in Tris–EDTA buffer (10 mM Tris, 1 mM EDTA, 0.05% Tween 20, pH 9.0) for 20 min in a steamer to achieve complete antigen-retrieval. After three washes in TBST buffer, the slides were covered with Power Block solution (BioGenex, CA, USA) for 10 min to avoid nonspecific signals. Following three washes in TBST buffer, the mouse anti-HPV antibody (K1H8, ab245950; Abcam, Cambridge, UK) was 1:200 diluted in PBS and subsequently applied onto the slides for 1-h incubation under room temperature. The endogenous peroxidase activity was blocked by incubating the slides with UltraVision Hydrogen Peroxide Block (Thermo Scientific, MA, USA) for 10 min under room temperature. Following three washes in TBST buffer, the slides were covered with the goat anti-mouse IgG antibody conjugated with horseradish peroxidase (HRP; Dako, Agilent Technologies, CA, USA) for 1-h incubation under room temperature. The signals were colorized by using AEC + high sensitivity substrate chromogen (Dako) for 10 min under room temperature. The subtract was removed and the slides were washed in TBST buffer. The counterstaining was performed by immersing the slides into hematoxylin (Muto Chemical) for 1 min and then, washing the slides with running tap water. The signals were visible under microscopic examination.

2.6. ISH Staining

ISH staining was performed to better depict the novel BPV in the lesion. The ISH procedure was performed as described by Wang et al. [24], with some modifications. The amplicons produced from PCR with CP4/5 was used to generated the digoxigenin-labeled DNA probe. Subsequently, the embedded tissues were sectioned to obtain 4 µm slides. The slides were deparaffinated and rehydrated as described above in the IHC staining. Antigen retrieval was achieved by soaking the slides in 0.1 N HCl solution followed by enzymic digestion with proteinase K (Qiagen) at 37°C for 15 min. Digoxigenin-labeled probes was added to the slides, then, the slides were heated on the 95°C plate for 6 min, followed by a 16–18 h incubation at 42°C. The endogenous peroxidase activity in the tissue was halted by incubating the slides in 3% H_2_O_2_ for 6 min. The potential nonspecific bindings were blocked by 10% normal goat serum. Probes were detected using the mouse anti-digoxigenin antibody (catalog no.11333062910, Roche) and anti-mouse IgG HRP secondary antibody (Dako). Liquid DAB + substrate chromogen system (Dako) was used to visualize the signals. The tissues on the slides were counterstained with hematoxylin for a better depiction.

2.7. Phylogenetic Analysis

Following the guidelines proposed by the ICTV of PV, the L1 ORF, which encodes the major viral capsid protein, serves as the basis for classifying and genotyping PVs. The full-length sequences of the L1 ORF were input to perform the phylogenetic analysis. The genetic information of 44 different genotypes of BPV classified in five genera was included and aligned by MEGA11 software. The phylogenetic tree was generated accordingly using the maximum likelihood method based on the Tamura-Nei model. The sequence identity was calculated using DNAStar software.

3. Results

3.1. Histological Findings

Eight clinical samples diagnosed as skin papillomas by local veterinarians were collected from six dairy farms through excisional or scraping biopsy (Table 1). Histological examination revealed that five samples (Case nos. 23-00, -01, -03, -05, and -07) were fibropapillomas, characterized by notable expansion of the stratum spinosum, connective tissue proliferation, and rete ridge formation (Figure 1). The Case nos. 23-02 and -06 were histologically diagnosed as papilloma (Figure 1). Obvious hyperkeratosis was also noted in all the cases. Eosinophilic intranuclear inclusions with chromatin margination were observed only in Case no. 23-00.

3.2. Molecular Diagnosis of PVs by Degenerative Primers

To detect PVs in the sample, two sets of degenerative primers were utilized to amplify the viral genome by PCR technique. All the cases (8/8) exhibited positive bands for CP4/5 (Figure 2A), while seven cases (7/8), excluding Case no. 23-02, showed positive bands for MY9/11 (Figure 2B). Based on the preliminary BLAST analysis on the NCBI database, the PVs in Case nos. 23-00, -01, -03, -04, -05, and -07 were temporarily determined as BPV1, while Case no. 23-06 was classified as BPV2. Interestingly, Case no. 23-02 displayed significant divergence (>10%) from currently known PVs and was temporarily categorized as an “unclassified” genotype.

3.3. IHC Staining and ISH Staining

To assess the presence of PV antigens in lesions, IHC staining with anti-HPV antibody (K1H8) was conducted. As shown in Figure 3, the positive signals in Case no. 23-00 appeared in both keratinized epithelial cells and the underlying hyperplastic fibroblasts, with a higher quantity and intensity of signals compared to other cases. Case nos. 23-01, -03, -05, and -07 showed few positive signals in epithelial cells, but none in the connective tissue. Case nos. 23-02 and 23-06 displayed no detectable signals. Since the PV detected in Case no. 23-02 was distant from the known genotype, to further elucidate the causal relationship between the novel BPV and lesion in Case no. 23-02, ISH staining was performed. As shown in Figure 4, the presence of genomic material of BPV in the lesion was clearly depicted by ISH staining. The positive signals mainly appeared at the nucleus of epithelial cells, whereas a small proportion of mesenchymal cells also had intranuclear positive signals.

3.4. Phylogenetic Analysis of BPVs

The full genome assembly of BPV in each case was obtained by overlapping Sanger sequencing with genotype-specific primers (Table A1), while the full genome assembly of the “unclassified” BPV was obtained by performing de novo NGS. The genetic information was uploaded onto the GenBank and the accession numbers for each case were assigned (Table 1). The full-length L1 ORF was aligned with those of reference strains representing 44 distinct BPV genotypes for the phylogenetic analysis. As shown in Figures 5 and 6, Case nos. 23-00, -01, -03, -04, -05, and -07 displayed high sequence identity (>98%) with the BPV1 reference strain (LC510378.1), thus, being classified as the variants of BPV1. Nevertheless, Case no. 23-03 still displayed notable genetic variability (~1.6%) compared to the reference and local strains, suggesting the possibility of local evolution of BPV in herds. Additionally, a BPV2 variant with 99.9% sequence identity to the reference strain (KC878306.1) was also detected. Remarkably, the BPV from Case no. 23-02 showed low sequence identities (55.8–80.9%) with known genotypes and formed a distinct subbranch under the genus of Xipapillomavirus (Figure 5), thus, temporarily designated as an “unclassified” genotype. The most closely related strain was BPV38 (MW404259.1), with only 80.9% similarity in the L1 ORF. Following taxonomic rules [10], the unclassified BPV was proposed as a novel genotype (putative BPV genotype 45) within the genus of Xipapillomavirus.

3.5. Full Genome and Putative Protein Analysis of BPV1 and BPV2

The full genome of BPV1 and BPV2 from the clinical samples was sequenced and assembled. The putative amino acids of each ORF were analyzed and aligned with the reference strains. Although the L1 ORF sequence identities among variants exceeded 98%, several nonsynonymous mutations were still noticed (Table 2). Notably, the BPV1 detected in Case no. 23-03 had the highest variation compared to both the reference strain and other local variants. Case no. 23-07 had a 10-amino acid mutation at the N-terminus of L2. This might indicate that within the same geographic region, BPVs of the same genotypes might independently undergo adaptation and accumulate mutations within their respective cattle populations. Our findings also revealed that the E7 protein of BPV1 exhibited greater conservation, whereas the L2 protein displayed a higher degree of divergence. Notably, all six variants displayed three unique substitutions compared to the reference BPV1 strain: K493R in E1, Q22P in E4, and S452N in L2, indicating potential local evolution of BPV1 in Taiwan. For BPV2, only three amino acid substitutions were identified throughout all viral proteins, including two in E6 and one in L2.

3.6. Full Genome and Putative Protein Analysis of the Novel BPV

The full genome of the novel “unclassified” BPV identified in this study was obtained through de novo NGS. As listed in Table 3, the genome spanned 7298 bp and comprised seven putative ORFs: E8, E7, E1, E2, E4, L2, and L1, in referring to BPV38. Amino acid sequence comparisons revealed that the L1, L2, E1, E2, E7, and E8 ORFs of this novel BPV shared 87%, 80%, 82%, 69%, 70.7%, and 76% identity, respectively, with those of BPV38. Several conserved functional motifs of PVs present in BPV38 were also identified in the novel BPV. In the E7 ORF, the retinoblastoma tumor suppressor protein (pRb)-binding motif (L_x_C_x_E) was characterized as L_N_C_E_E, alongside a zinc-binding domain (C_xx_C_x29_C_xx_C) represented as C_YV_C_x29_C_AS_C. Both motifs were recognized regulators of the host's cell cycle. The E8 ORF of the novel BPV contained a hydrophobic transmembrane domain with conserved C-terminal residues (L_x_GWD), which had been linked to immune evasion. Additionally, a conserved ATP-binding site (G_x4_GKS) and an RXL cyclin-binding motif (KRRLL) were both identified in the E1 ORF. The L1 ORF featured two E1-binding sites (E1BS and AACAAT) and one modified E1BS (TAACAA). In the LCR, one modified E1BS (TAACAA) and a canonical E2-binding site (E2BS and ACCG_x4_CGGT) were also found. The ATP-binding site, E1BS, E2BS, and RXL motif were associated with viral genome replication. Notably, no leucine-zipper domain (L_x6_L_x6_L_x6_L_X6_L) was detected in the E1 and E2 ORFs.

4. Discussion

Bovine PV (BPV) has been strongly linked to benign and malignant proliferative cutaneous lesions in cattle [1, 2]. Unlike most PVs, which are strictly host-restricted and highly epitheliotrophic, BPVs are unique for interspecies transmission and broad tissue tropism [18, 25]. To investigate the BPV in lactating farms in Taiwan, eight cutaneous papilloma samples from six farms were collected and analyzed through histopathology, molecular diagnosis, and IHC staining as summarized in Table 4. The findings indicated that BPV1 was the most prevalent genotype causing bovine papillomas in Taiwan, although BPV2 was also detected. Notably, a novel BPV genotype (putative BPV45), sharing only 80.9% L1 identity with BPV38, was identified and classified within Xipapillomavirus based on full-genome sequencing and phylogenetic analysis. ISH further confirmed its involvement in the development of hyperplastic cutaneous lesions. These findings contribute valuable insights into BPV infections in Taiwan and expand the global understanding of BPVs.

According to published reports, infections with Delta-PV, particularly BPV1 and BPV2, often involved both epithelial cells and stromal cells, leading to fibropapilloma in cattle or sarcoid in horses [26, 27]. BPV5, an Epsilon-PV, could also cause mesenchymal transformation [28]. In contrast, Xipapillomavirus infections in bovine were typically restricted to epithelial tissues [27]. In this study, a novel BPV in the Xipapillomavirus was identified. Despite the lesion caused by this novel BPV was histologically diagnosed as papilloma, the genetic material of this novel BPV was detected in both proliferating epithelial cells and the subepithelial mesenchymal cells by ISH staining, suggesting an expanded tissue tropism for this novel BPV.

BPVs classified within the Xi- genus were known for lacking the E6 ORF [27]. Seven putative ORFs and a noncoding LCR were annotated for the novel BPV. The presence of the pRb-binding motif on the E7 suggested its potential to transform epithelial cells [27, 29]. Interestingly, according to the published reports, the absence of the pRb-binding motif in BPV1, 2, and 5 resulted in the lesion manifesting as a fibropapilloma rather than a papilloma [30]. Consistent with this, all BPV1-positive cases in this study were diagnosed as fibropapillomas, while the BPV2-positive case was diagnosed as a papilloma due to minimal stromal proliferation. Although the novel BPV was detected in mesenchymal cells by ISH staining, no significant stromal proliferation was observed, possibly reflecting the early disease stage and supporting previous findings linking the functional pRb-binding motif with reduced fibropapilloma incidence [30]. Additionally, a transmembrane domain and the conserved hydrophilic tail (L_x_GWD) were identified in E8, potentially contributing to immune evasion by modulating MHC expression [31, 32].

The prevalent BPV genotypes vary across countries and even among cattle breeds. Large-scale surveys had shown BPV1 or BPV2 to be predominant in the dairy sector in Egypt, Italy, Turkey, Iraq, Brazil, China, and Japan [19, 21, 33–37], while BPV8 and BPV10 were most common in Germany [38] and BPV7 and BPV8 dominated in Poland [39]. These regional differences in genotype distribution demonstrated the need for local epidemiological investigations. Our study contributed to filling this gap by providing the first molecular and pathological evidence of BPVs in Taiwan. Additionally, many recent reports have highlighted coinfections with different BPV genotypes in lesions [12, 18]. In the present study, BPV1 and BPV2 were detected in different cases by conventional PCR technique with two consensus degenerative primer sets for PVs. By performing conventional PCR with these primers, the possibility of coinfections cannot be completely excluded. Nevertheless, previous research has shown that BPV1 and BPV2 can independently cause hyperplastic or neoplastic lesions in cattle [25, 40, 41], which was consistent with our findings. In Case no. 23-02, the novel BPV was initially detected using degenerative primers, subsequently, de novo NGS was employed to assemble the complete genome and also exclude the possibility of coinfection, further confirming that this novel BPV alone could induce epithelial proliferation leading to papilloma.

To further confirm the presence of BPV in the lesions, IHC staining was performed with the anti-HPV monoclonal antibody (K1H8). Although this antibody exhibited some cross-reactivity [42], positive signals were detected only in the BPV1-positive cases, with no signals observed in BPV2 and the novel BPV cases. Among the IHC-positive tissues, only Case no. 23-00 presented strong positive signals in both epithelial cells and mesenchymal cells, while the other BPV1-positive cases exhibited only a limited number of positive signals. The target of this anti-HPV antibody is the viral capsid protein encoded by the late transcription gene [42]. In early infection stages, viral capsid proteins may not be produced or accumulated in significant quantities, making the detection challenging. Given the lack of a specific antibody, ISH staining with genotype-specific probes was a useful technique to confirm the presence of the BPV in the lesions, highlighting the complementary role of ISH in diagnosing BPV infections.

In summary, a novel BPV genotype, which was suggested to be BPV genotype 45, was identified and characterized from a papilloma in a lactating Holstein. While BPV1 was also detected from hyperplastic lesions from different individuals in the same farm, coinfection of BPVs in the same lesion with the novel BPV was excluded, confirming the independent transforming potential of this novel BPV. Furthermore, IHC and ISH staining revealed that both BPV1 and the novel BPV exhibited tropism for epithelial and connective tissues, expanding current knowledge on BPV tissue specificity. This study also highlighted the circulation of BPV1, BPV2, and the novel BPV in Taiwan, with BPV1 being the most prevalent in skin papillomas.