Introduction

Bats are considered the natural reservoirs of a variety of zoonotic RNA viruses, such as Ebola viruses, Marburg viruses, and severe acute respiratory syndrome coronavirus [1,2,3]. Several DNA viruses, including adenoviruses, herpesviruses, and polyomaviruses (PyVs), have also been detected [4,5,6]. However, the pathogenic and zoonotic roles of these DNA viruses have been not clarified.

PyVs are small double-stranded DNA viruses with a circular genome of approximately 5 kbp. The viral genome consists of three regions: regulatory, early, and late regions. The regulatory region is responsible for transcription from both the early and late promoters and the initiation of viral DNA synthesis. The early region contains genes encoding the large T antigen (TAg) and small t antigen (tAg). The late region contains the structural proteins VP1, VP2, and VP3 [7]. Although PyV diversity in bat populations in North, Central, and South America, Africa, Indonesia, and New Zealand were investigated in previous studies [6, 8,9,10,11,12], the prevalence and genetic diversity of PyVs in Japanese bats remain unclear. The aims of this study were to (i) determine the presence of PyVs in Japanese bats, (ii) characterize the genomic structure of bat PyVs, and (iii) analyze the evolutionary relationships between the bat PyV detected in this study and other known bat PyVs.

Results and discussion

Eighteen bats (Miniopterus fuliginosus) were collected in Wakayama Prefecture, Japan. The pooled intestinal contents obtained from each bat (approximately 1 g/body) were prepared as a 10% suspension in sterilized phosphate-buffered saline (PBS). The suspension was centrifuged at 17,000×g for 30 min at 4 °C. The supernatant was filtered three times through Qualitative Filter Paper Grade I (GE Healthcare, Piscataway, NJ, USA) to remove large debris and then through a 0.45-μm bottle top filter (Thermo Scientific, Waltham, MA, USA) to remove bacterium-sized debris. The filtrate was centrifuged at 107,000×g for 3 h at 4 °C. After the supernatant was removed, the precipitate was suspended after overnight incubation with 500 μl of PBS at 4 °C. The enriched fraction containing the predicted viral capsid was obtained via cesium chloride (CsCl) density gradient (1.7, 1.5, 1.35, and 1.2 g/ml) centrifugation at 234,000×g for 18 h at 20 °C. Next, the obtained fraction (1.35–1.5 g/ml) was diluted tenfold with PBS and centrifuged at 107,000×g for 3 h at 4 °C to remove CsCl. After the supernatant was removed, the precipitate was suspended again in 100 μl of PBS and treated with Deoxyribonuclease (RT Grade) for Heat Stop (NIPPON GENE, Tokyo, Japan) and RNase ONE (Promega, Madison, WI, USA) to remove any extra-viral capsid nucleic acids. DNA extraction was performed using a QIAamp DNA Mini Kit (QIAGEN) according to the manufacturer’s instructions. The DNA was then randomly amplified using a Genome Plex Complete Whole Genome Amplification Kit (SIGMA) with random hexamer primers. Amplicons were purified through MobiSpin S-400 columns (MobiTec, Gottingen, Germany). Metagenomic analysis was performed using Illumina HiSeq (Roche, Basel, Switzerland). The genomic data were assembled using GS De Novo Assembler ver. 2.8 computer software (Roche).

A total of 10,136,210 reads and 5362 contigs were obtained in a pool of sample from 18 bats. To identify homologous sequences, the obtained genomic data were analyzed via a BLASTn search using the DNA Data Bank of Japan database in accordance with a previously reported method [13, 14]. Virus-related sequences were identified in 123 contigs. Of these, 14 contigs contained PyV-like sequences with high identities. Other contigs contained the sequences of Eel River basin pequenovirus, Montastraea cavernosa colony-associated virus, and Grapevine-associated totivirus-1. To determine the complete viral genome of these PyV-like sequences, PCR was performed using LA Taq DNA polymerase (Takara Bio, Otsu, Japan) in accordance with the manufacturer’s instructions. Specific PCR primers were designed on the basis of the sequences obtained from the contigs. The primer sequences were as follows: Pyv-F1 (sense, 5′-AAG TTT GCA GTA GTC TTT GAA GAT GTG AAG GGT C-3′), Pyv-R1 (antisense, 5′-CAC TCC TGG GCT TTC CTG CTC ATA TTT ATG CA-3′), Pyv-F2 (sense, 5′-CAT AAA CAG GGT CAA ACC AC-3′), and Pyv-R2 (antisense, 5′-AAG CAC TCC ACC AAA GGA AA-3′). DNA extracted from the pooled sample of bat intestinal contents was used as the template. The PCR products were visualized via electrophoresis on a 1% agarose gel stained with SYBR Safe (Life Technologies, Carlsbad, CA, USA). The full genomic DNA could be amplified by two independent PCR using the aforementioned primers. The amplified DNA was cloned by inserting the PCR product into the pCR2.1 TOPO vector (Life Technologies) in accordance with the manufacturer’s instructions. The obtained sequences were analyzed using the BigDye Terminator v3.1 Cycle Sequencing Kit (Life Technologies), and nucleotide sequences were assembled using ATGC computer software (Genetyx Corporation). A homology search was performed using NCBI BLAST.

The genome of Miniopterus fuliginosus polyomavirus (MfPyV) has a length of 4956 bp (accession number: LC529726). The genome organization includes an early region coding for tAg and TAg on one strand and a late region encoding the capsid proteins VP1, VP2, and VP3 on the opposite strand. A noncoding regulatory region (NCCR) was located between the start of the early region and that of the late region, in line with previous findings for bat PyVs (Fig. 1a and Supplementary Table 1) [9,10,11,12]. Interestingly, open-reading frames encoding VP2 and VP3 of MfPyV did not overlap with that of VP1. The stop codons of VP2 and VP3 are located at base positions 1184–1186, whereas the start codon of VP1 is located at base positions 1188–1190. A single nucleotide (guanine) at 1187 separates the VP2/3 and VP1 regions (Fig. 1a and Supplementary Fig. 1). Therefore, genomic composition of MfPyV is genetically different from those of typical PyVs in terms of non-overlapping VP regions. Figure 1b and c present the phylogenetic trees of VP1 and TAg of MfPyV and 28 other bat PyVs constructed using neighbor-joining analysis. Based on phylogenetic analyses of the VP1 and TAg amino acid sequences, both regions of MfPyV are closely related to those of other Miniopterus PyVs and group B bat PyVs (Fig. 1b and c). VP1 is a major PyV structural protein that is indispensable for entry of the virus into host cells [7]. MfPyV VP1 displayed less than 72% nucleotide sequence identity with other bat PyVs (Supplementary Table 2). TAg is a multifunctional protein that plays important roles in viral DNA replication and the regulation of viral and cellular gene expression [15,16,17]. The predicted MfPyV TAg exhibited low similarity (< 73%) with those of other PyVs (Supplementary Table 2). MfPyV TAg sequences contained features known to be conserved in TAgs of other bat PyVs, including the highly conserved DnaJ domain (HPDKGG), a retinoblastoma (Rb)-binding motif (LYCNE), and several functional motifs (Supplementary Fig. 2). According to a previous report, these elements work together to bind Rb and interrupt its interaction with the E2F transcription factor to promote viral replication and cell cycle progression [18]. TAg is generated via alternative splicing of the early mRNA transcript [11, 19]. In the early region of the MfPyV genome, conserved predicted splice donor sites are located at base positions 4729–4734 (CCTGAG/GTAAGG) and 4346–4351 (TTTCAG/GTCTTC) (Fig. 1a). In the deduced NCCR region of the MfPyV genome, several conserved elements were identified, including several copies of the consensus TAg binding site GAGGC and its reverse complement GCCTC (Supplementary Fig. 3). These elements are likely to comprise the core of the replication origin [20].

Fig. 1
figure 1

a Genome organization of Miniopterus fuliginosus polyomavirus (MfPyV). The entire dsDNA viral genome consisted of 4956 bp. The viral genome appeared to have a noncoding regulatory region, an upstream regulatory region, potential open-reading frames for the late proteins VP1, VP2, and VP3, and the early proteins small t antigen and large T antigen (TAg). b and c Phylogenetic relationships between the amino acid sequences of VP1 (b) and TAg (c) of bat polyomaviruses. MfPyV is indicated by bold text and representative clusters are shown as A-E in the phylogenetic trees. The sequences for the other polyomaviruses were obtained from GenBank. Phylogenetic analysis was performed using the neighbor-joining method with 1000 bootstrap replicates. Bars indicate 0.1 amino acid residue replacements

Full size image

Comparison of the full-length genome sequence of MfPyV with those of other bat PyVs revealed that MfPyV is most closely related to the KY156 strain with 70% nucleotide sequence identity (Supplementary Table 2). According to the Polyomaviridae Study Group of the International Committee on Taxonomy of Viruses, a novel PyV species should have < 81% nucleotide sequence identity to other known PyVs [7]. MfPyV exhibited less than 81% nucleotide sequence homology to the known reference PyVs including previously reported bat PyVs. In line with the nomenclature of the other bat PyVs, we propose the name MfPyV for the newly discovered virus.

For virus isolation, we attempted to propagate the MfPyV strain using the Tb1-Lu cell line derived from the lungs of the free-tailed bat Tadaria brasiliensis (ATCC #CCL-88). However, a cytopathic effect was not observed in the cells following serial passage of the cultures. Viral DNA replication was also not detected in the cells and supernatant collected at each passage. There is a need for additional research to identify efficient cell culture systems for bat PyVs to elucidate the viral infection/replication mechanisms and their pathogenicity.

In conclusion, we detected a novel PyV genome sequence in Japanese bats. Further epidemiological investigations are needed to determine the extent of PyV genetic variation in various bat species in Japan.