Nasal virome of dogs with respiratory infection signs include novel taupapillomaviruses.

Using viral metagenomics, we characterized the mammalian virome of nasal swabs from 57 dogs with unexplained signs of respiratory infection showing mostly negative results using the IDEXX Canine Respiratory Disease RealPCR™ Panel. We identified canine parainfluenza virus 5, canine respiratory coronavirus, carnivore bocaparvovirus 3, canine circovirus and canine papillomavirus 9. Novel canine taupapillomaviruses (CPV21-23) were also identified in 3 dogs and their complete genome sequenced showing L1 nucleotide identity ranging from 68.4 to 70.3% to their closest taupapillomavirus relative. Taupapillomavirus were the only mammalian viral nucleic acids detected in two affected dogs, while a third dog was coinfected with low levels of canine parainfluenza 5. A role for these taupapillomavirues in canine respiratory disease remains to be determined.

Introduction

Canine infectious respiratory disease (CIRD) is a major cause of morbidity and an important welfare issue for kennelled dog populations worldwide. The canine respiratory disease is generally characterized by coughing, snorting, wheezing, sneezing, nasal and eye discharge, fever, and dry and itchy nose. In some cases, the disease can be lethal if bronchopneumonia develops [1]. Viral causes for respiratory disease in dogs include canine adenovirus [2], canine respiratory coronavirus, canine herpes virus, distemper virus [3], influenza virus [4], canine parainfluenza virus [1] and a recently described canine protoparvovirus [5].

The Papillomaviridae family currently consists of at least 53 genera formed by 133 papillomavirus (PV) species that have been identified in humans, non-human mammals, birds, reptiles and fish [6]. Transmission of warts between dogs was first shown in 1898 [7] and cell-free extracts were shown to transmit canine papillomavirus (CPV) in 1959 [8]. Canine oral papillomas are typically benign and rapidly regress rarely causing respiratory problems [8] except in immuno-suppressed dogs which can develop severe oral papillomatosis [9-11].

The first CPV genome was sequenced in 1986 [12]. To date a total of 20 CPV have been described that can be classified into 3 different genera [13-16]. Most recently canine papillomavirus type 18 genome was characterized in pigmented plaques [17], canine papillomavirus 19 was sequenced from an oral papilloma coinfected with CPV1 and CPV2 [16], and the genome of canine papillomavirus 20 was released in 2016 (GenBank accession: KT901797). Of these 20 canine papillomavirus genomes, five are classified in the Taupapillomavirus genus (CPV 2,7,13,17,19), 2 in the lambda genera (CPV 1,6) and 13 in the Chi genera (CPV3,4,5,8,9,10,11,12,14,15,16,18,20). Lambdapapillomavirus CPV1 has been shown to cause oral papillomatosis [18]. CPV2, 7, 13, 17 and 19 within the taupapillomavirus genus were reported in endophytic and exophytic papillomas, squamous cell carcinoma, oral papillomatosis and oral squamous cell [16, 17, 19, 20].

Here we characterized the eukaryotic virome of nasal swabs from dogs with unexplained respiratory infection signs and characterize three novel taupapillomavirus genomes.

Materials and methods

Canine samples

Nasal swab samples were collected from dogs with unexplained respiratory disease. With a few exceptions discussed below, these samples had pre-tested negative for IDEXX Canine Respiratory Disease RealPCR™ Panel. A total of 57 nasal swab samples from different dogs were first analyzed in 11 pools of 1–7 samples using a viral metagenomic approach.

Viral metagenomics

Nasal swab samples were diluted with 750 µl Dulbecco’s phosphate-buffered saline, pooled, clarified by 15,000 g centrifugation for 5 min, and supernatants filtered using a 0.45-µm filter (Millipore). Free nucleic acids in the filtrates (not protected in viral capsids) were digested using DNAse and RNAse enzymes to enrich for viral nucleic acids. Nucleic acids were then extracted (MagMAX Viral RNA Isolation Kit, Ambion, Inc, Austin, Tx, USA) [21] and amplified by random RT-PCR followed by use of the Nextera™ XT Sample Preparation Kit (Illumina) to generate a library for Illumina MiSeq (2 × 250 bases) with dual barcoding as previously described [22].

Following de novo assembly using the Ensemble program [23], both contigs and singlets viral sequences were then analyzed using translated protein sequence similarity search (BLASTx v.2.2.7) to all annotated viral proteins available in GenBank. Geneious R10 was then used to align reads and contigs to reference viral genomes from GenBank the program.

Generation of full genomes of taupapillomaviruses

Back to back (inverse) PCR primers were designed based on metagenomics-derived sequences matching papillomaviruses (Supplementary material 1). PCR was carried out in a total volume of 25 µl. PCR mixes contained 3 µl of extracted DNA, 1.25 units Takara La Taq DNA polymerase (Takara), 10 pmol of each primer sets dNTP mixture (Takara) and the supplied 10X La PCR buffer II (Takara). First round PCR products were used as template in a 25 µl nested PCR. After an initial denaturation at 94 °C for 1 min, amplification was performed for 40 cycles consisting of 20 s at 94 °C, 30 s at 55 °C and 2 to 7 min (depending on PCR product size) at 72 °C, followed by a final extension for 10 min at 72 °C. Long PCR amplicons were sequenced using the Nextera library making kit and the MiSeq Illumina sequencer. PCR primers used to identify individual taupapilllomavirus positive samples in the pool samples are described (Supplementary material 1).

Phylogenetic analysis

The E1-E2-L1-L2 protein sequences were aligned using MAFFT [24, 25] in Geneious v10.1.3., and phylogenetic tree was derived using the maximum likelihood method based on the general time reversible model [26]. The percentage of trees in which the associated taxa clustered together is shown next to the branch points. Initial tree(s) for the heuristic search were obtained by applying the neighbor-joining method to a matrix of pairwise distances estimated using the maximum composite likelihood (MCL) approach. The tree is drawn to scale, with branch lengths measured in the number of substitutions per site. The analysis involved 27 nucleotide sequences of all currently known canine papillomaviruses. Positions containing gaps and missing data were eliminated. There were a total of 5051 positions in the final dataset. Evolutionary analyses were conducted in MEGA6 [27].

The L1 gene nucleotide sequences were aligned using MAFFT [24, 25] in Geneious v10.1.3., and phylogenetic tree was derived using the neighbor-joining method [28]. The tree is drawn to scale, with branch lengths in the same units as those of the evolutionary distances used to infer the phylogenetic tree. The evolutionary distances were computed using the maximum composite likelihood method [29] and are in the units of the number of base substitutions per site. The analysis involved 27 nucleotide sequences. There were a total of 1423 positions in the final dataset. Evolutionary analyses were conducted in MEGA6 [27].

Results

Viral metagenomic

Mammalian viruses detected with BLASTx-translated protein matches (E score < 10−30) were: rubulavirus 5 (also known as canine parainfluenza virus 5 or PIV5) in 3 different pools, carnivore bocaparvovirus 3 in two pools, canine betacoronavirus in one pool, canine circovirus in one pool, and canine papillomavirus 9 (Chi genus) in one pool (Table 1). Three pools also contained reads showing E score < 10−30 to taupapillomaviruses (Table 1). PCR based on these taupapillomavirus sequences were first used to identify the individual positive samples within each of the three pools (Supplementary material 1). These three samples were then deep sequenced individually (rather than in pools) to confirm the presence of the taupapillomavirus, acquire more of their genome sequences and test for coinfecting viral genomes.

Table 1

Distribution of sequence reads to different viruses in 11 nasal swab pools and 3 individual taupapillomavirus positive swab samples

	No. of sample in pool	Total reads	Circoviridae	Paramyxoviridae	Papillomaviridae		Parvoviridae	Coronaviridae
			Canine circovirus	Parainfluenza virus 5	Taupapillomavirus	Chi papillomavirus	Carnivore bocaparvovirus 3	Betacoronavirus 1
Pool1	5	497,208		20	29,029
Pool2	5	1,516,776				131	2
Pool3	5	523,486					2
Pool4	5	2,493,284
Pool5	5	667,534
Pool6	6	606,050		8
Pool7	7	1,125,636
Pool8	7	673,820
Pool9	7	674,772			287
Pool10	1	1,107,958			665
Pool11	4	930,598	4278	31				11
Case 1 (individual sample in Pool1)		202,354		8	3084
Case 2 (individual sample in Pool9)		353,980			2424
Case 3 (individual sample in Pool10)		399,948			1812

Taupapillomavirus positive animal

Clinical case 1 (in pool 1) shedding taupapillomavirus CPV21 was a 12-month-old intact male Labrador retriever with respiratory infection signs. The nasal swab of case 1 was tested by the IDEXX Canine Respiratory Disease RealPCR™ Panel consisting of real-time PCRs for Bordetella bronchiseptica, canine adenovirus type 2, canine distemper virus, canine herpesvirus type 1, canine parainfluenza virus, canine pneumovirus, canine respiratory coronavirus, H3N2 canine influenza virus, other influenza A virus (H1N1, H3N8, H7N2), Mycoplasma cynos, and Streptococcus equi subsp. zooepidemicus. Case 1 was positive only for Mycoplasma cynos. Two other dogs from the same household also exhibited similar respiratory infection signs. One of these dogs was also tested with the IDEXX Canine Respiratory Disease RealPCR™ Panel with negative results for all agents. Swabs from both of these related dogs tested negative for taupapillomavirus CPV21 by both metagenomic sequencing and PCR.

Clinical case 2 (in pool 9) infected with CPV22 was an eleven months old, intact female Miniature Australian Shepherd with sneezing and coughing since boarding at a clinic. All tests on the IDEXX Canine Respiratory Disease RealPCR™ Panel were negative.

Clinical case 3 (in pool 10) infected with CPV23 was a 1-year 2-month-old male Pembroke Welsh Corgi. This animal was part of an outbreak of tracheobronchitis in up to 20 dogs in a canine daycare center with upper respiratory infection signs of coughing with occasional clear phlegm/saliva. No ocular or nasal discharges were reported. All tests on the IDEXX Canine Respiratory Disease RealPCR™ Panel on case 3 were negative. The same respiratory panel was also used on two other patients in same outbreak which were positive for Mycoplasma cynos only. Neither of these two nor other animals from that outbreak were available for metagenomics or taupapillomavirus PCR.

Results of full genomes generation

Two to four hundred thousand sequence reads were generated for each of the three individual respiratory samples (Table 1). The raw sequence data for each of these 3 nasal swab samples are available at NCBI’s Short Reads Archive under GenBank accession number SRP133969.

The most commonly detected eukaryotic viral reads in these three nasal swabs belonged to the Taupapillomavirus genus in the Papillomaviridae family. The frequency of reads belonging to the assembled genomes ranged from 0.45 to 1.5% of the total read numbers. One swab from case 1 also contained nucleic acids from another canine virus, canine parainfluenza virus 5 (PIV5) at the low read frequency of 0.0032% (8 out of 202,354 reads).

Taupapillomavirus genomes

Gaps in the assembled genome sequences were filled by long-range PCR using the primers described (Supplementary material 1). The full genome of the canine taupapillomaviruses was 8140 (CPV23), 8300 (CPV22) and 8225 bp (CPV21) in length, with 46.8 to 48.40% GC content coding for five early proteins (E1, E2, E4, E6 and E7) and two late proteins (L1 and L2) (GenBank accession numbers; MH285952 for CPV21, MH285953 for CPV22, MH285954 for CPV23, Supplementary material 2). The predicted ORFs and characteristics of their putative protein products are shown in supplementary material 3. As customary the ‘‘A’’ in the first start codon of the ORF E6 was assigned position 1. The ORF genetic organization of CPV21-23 was identical to that of other taupapillomaviruses.

All E6 and E7 proteins showed conserved zinc-binding domains (two domain CXXC-X29-CXXC for E6 and single domain CXXC-X30-CXXC for E7) but lacked for PDZ-binding motif (ETQL) and retinoblastoma (pRb) protein-binding site that is present in high-risk oncogenic PVs [30] although it has been shown for canine papillomavirus 2 that a distinct region is responsible for pRb binding and degradation [31]. The long control region [32] between E2 and L2 was 648, 714 and 842 bases long. Conserved helicase domain, the conserved ATP-binding site of the ATP-dependent helicase (GXXXXGK[T/S]), cyclin A interaction motifs and the leucine–zipper domain were identified within the hypothetical proteins (Supplementary material 4).

Maximum likelihood (ML) tree based on the concatenated alignment of E1, E2, L1 and L2 from new taupapillomavirus tree with bootstrap support for the position of the branches is shown in Fig. 1a. The novel taupapillomaviruses clustered with the Taupapillomavirus genus with high bootstrap support.

Fig. 1

Phylogenetic analysis of all currently sequenced canine papillomaviruses. Genomes from this study are marked with solid triangles. a Maximum likelihood (ML) tree based on the concatenated alignment of E1, E2, L1 and L2 from new taupapilloma virus and 24 other PV types of different species and genera. b The L1 gene nucleotide sequences were aligned using MAFFT in Geneious v10.1.3., and phylogenetic tree was derived using the neighbor-joining method. Abbreviations used include canine papillomavirus, CPV; Felis catus (cat) papillomavirus, FcaPV; Mustela putorius papillomavirus (ferret), MpPV; Vulpes vulpes papillomavirus (fox), VvPV

The three L1 ORF showed nucleotide identity of 67.3 to 69.6% to each other and from 68.4 to 70.3% with their closest relative, canine papillomavirus 13 isolate Zurich/2011 (CPV13, GenBank accession number JX141478) [20] and 52 to 63.5% similarity with other canine papillomavirus L1 (Table 2).

Table 2

Percent identity of L1 nucleotide sequences between the taupapilomaviruses described here and all previously reported canine papillomaviruses

Accession number_name_genus	CPV21	CPV22	CPV23
D55633_CPV1_Lambda	54.516	55.838	54.664
AY722648_CPV2_Tau	60.947	59.698	59.711
DQ295066_CPV3_Chi	53.746	54.022	53.399
EF584537_CPV4_Chi	52.047	52.838	52.577
FJ492743_CPV5_Chi	52.789	54.622	53.906
FJ492744_CPV6_Lambda	55.858	54.386	54.646
FJ492742_CPV7_Tau	63.57	62.187	62.788
HQ262536_CPV8_Chi	54.902	54.134	54.659
JF800656_CPV9_Chi	54.142	54.879	54.42
JF800657_CPV10_Chi	55.708	54.944	54.486
JF800658_CPV11_Chi	54.142	54.42	53.373
JQ754321_CPV12_Chi	53.49	54.551	54.158
JX141478_CPV13_Tau	68.445	69.048	70.363
JQ701802_CPV14_Chi	54.729	54.028	53.635
JX899359_CPV15_Chi	55.033	54.724	56.168
KP099966_CPV16_Chi	54.929	54.581	52.381
KT272399_CPV17_Tau	62.467	61.286	61.166
KT326919_CPV18_Chi	52.964	53.172	52.745
KX599536_CPV19_Tau	62.846	62.582	61.012
KT901797_CPV20_Chi	53.947	53.569	53.111
MH285952_CPV21_Tau	100	69.625	69.226
MH285953_CPV22_Tau	69.625	100	67.324
MH285954_CPV23_Tau	69.226	67.324	100

According to the phylogenetic tree constructed with the L1 nucleotide region, the papillomavirus genomes in nasal swab clustered with the taupapillomaviruses, most closely with CPV13 (described in 2011 from a dog with papillomatosis) [16, 17, 19, 20] (Fig. 1b). The suggested ICTV criteria for naming a new genera require an L1 nucleotide identity of approximately < 60% [33-35]. A distribution of pair-wise distances shows that L1 sequences classified with the same species generally showed identity of 70-100% [33, 34]. Given these criteria, the latest ICTV release for the Papillomaviridae family (June 2018) list four species in the Taupapillomavirus genus. Taupapillomaviruses 1 and 2 consist of only canine viruses. Taupapillomavirus 3 is represented by feline viruses (FcaPV3, FcaPV4) [36, 37] and Taupapillomavirus 4 by a ferret papillomavirus (MpPV1) [38, 39] (Fig. 1). The L1 nucleotide identity between the new canine taupapillomaviruses described here and their closest relative, CPV13 in taupapillomavirus species 2, ranged from 68.4 to 70.3% (Table 2). The new taupapillomaviruses described here (CPV21-23) can therefore most conservatively be considered new types within Taupapillomavirus 2 species.

Discussion

We analyzed the virome of nasal swabs from fifty-seven dogs with unexplained respiratory infection signs. Reads from three novel taupapillomaviruses were identified and their complete genomes sequenced (CPV21-23). Genomes of canine parainfluenza 5, carnivore bocaparvovirus 3, canine betacoronavirus, canine circovirus, and canine papillomavirus 9 were also detected. All of these viruses except canine circovirus and canine papillomavirus 9 have already been associated with canine respiratory disease although with different level of pathogenicity [40].

It is not known whether papillomatosis was present in the 3 taupapillomavirus containing clinical cases described here. Papillomavirus-induced signs such as epithelial hyperplasia may not be readily visible unless worsened to the level of papillomatosis due to immune deficiency [41, 42]. For clinical case 1, the taupapillomavirus could only be detected by PCR in one of three dogs from the same household showing respiratory infection signs indicating that the taupapillomavirus may be unrelated to this small respiratory disease outbreak or that the samples from the other two dogs were collected after CPV21 was cleared. Clinical case 1 also yielded canine parainfluenza 5 at low level (based on low % of reads and lack of detection using the commercial IDEXX PCR assay). This paramyxovirus is a well-documented canine respiratory pathogen [43].

Case 2 was a sporadic case of respiratory infection signs where the only mammalian virus detected by metagenomics was taupapillomavirus. All the commercial IDEXX PCR assays were negative.

Clinical case 3 was part of a larger respiratory outbreak of 20 dogs, three of which were tested on the IDEXX respiratory panel and were negative except for one case of Mycoplasma cynos, an infection also commonly seen in asymptomatic dogs [43]. Taupapillomavirus CPV23-infected clinical case 3 was the only dog from that outbreak available for this study.

Evidence for a pathogenic role of taupapillomavirus infection in these 3 dogs with respiratory infection signs therefore remains limited to their detection in nasal swabs in the absence of known respiratory pathogens in 2 of these 3 dogs. Detection of taupapillomavirus DNA may also reflect frequent asymptomatic or chronic infections unrelated to these animals’ respiratory disease. Non-viral pathogens such as untested-for bacteria or fungi may have also been involved in these dog’s respiratory signs. Coinfections with parainfluenza 5 and/or Mycoplasma cynos may have also worsened the outcome of taupapillomavirus infection as could have immunodeficiency as shown by induction of papillomatosis by other taupapilomaviruses [41, 42].

Further studies will be required to determine the role of taupapillomaviruses in canine respiratory diseases. Such studies may include localizing the precise sites of viral replication in the canine respiratory track using in situ RNA hybridization and of their association with papillomatosis or other respiratory disease signs using matched case-control studies.

Below is the link to the electronic supplementary material.

Supplementary material 1. Primers were used in this study. PCR primers used to identify taupapilllomavirus positive samples in the pool samples and complete the genomes (XLSX 12 KB)

Supplementary material 2. ORF maps of CPV21, 22, 23 and CPV13. GenBank accession numbers; MH285952, MH285953, MH285954 (PPTX 424 KB)

Supplementary material 3. Predicted open reading frames (ORFs) in the genome of CPV21, 22, 23. The predicted ORFs and characteristics of their putative protein products (XLSX 10 KB)

Supplementary material 4. Genome motifs of CPV21, 22, 23. Table include ; Genome length, GC content, the long control region, zinc binding domains, ATP binding site, leucine-zipper domain, Cycling A interaction motif, Polyadenylation signal, E1 C-terminal, E1 N-terminal (DOCX 15 KB)