The sampling procedures were approved by the Ethics Committee of the Hainan Medical University. A total of 108 fecal samples were collected from 108 different crocodiles (Crocodylus siamensis) which were fed with tilapia or chicken at three crocodile breeding sites located in the counties/cities of Wenchang and Lingshui of Hainan Province (Figure 1). Following the discovery of the tilapia parvovirus HMU-HKU, a total of 70 intestine samples were additionally collected from 70 different tilapia (Oreochromis niloticus) at four tilapia breeding sites located in four counties/cities (Ding’an, Dongfang, Haikou, and Qionghai) of Hainan Province, China (Supplementary Figure 1). The samples were immersed in maintenance medium in virus-sampling tubes (Yocon/Beijing, China, which is composed of Hank’s, gentamicin, fungal antibiotics, bovine serum albumin) and were transported to the laboratory within 24 h using cold chain transportation and stored at −80°C. The tilapia species was determined morphologically and confirmed by DNA sequence analysis of mitochondrial cytochrome b (Kiril’chik et al., 1995).

Eight fecal samples from crocodiles fed with tilapia in Lingshui County were pooled and processed using a virus-particle-protected, nucleic acid purification protocol and nucleic acid library preparation method, as described previously (Wu et al., 2018). Amplified viral nucleic acid libraries were analyzed on an Illumina HiSeq2500 sequencer, for a single read of 100 bp in length. The raw sequence reads were filtered using previously described criteria to obtain valid sequences (Yang et al., 2011).

Sequence-similarity-based taxonomic assignments were conducted as described previously (Yang et al., 2011). Briefly, each read was evaluated for viral origin by conducting alignments in the National Center for Biotechnology Information (NCBI) non-redundant nucleotide database (NT), and protein database (NR), using BLASTn and BLASTx (-e: Expected value < 10^–5, -F: Filter query sequence, default = T). Taxonomies of aligned reads with the best BLAST scores (E-value < 10^–5) from all lanes were parsed and exported using the MEGAN 6- MEtaGenome Analyzer (Tamura et al., 2013).

The molecular clues provided by metagenomic analyses were used to classify sequence reads into a virus family or genus using MEGAN 6. Viral nucleotide (including both viral genome DNAs and viral mRNAs) was isolated using a QIAmp MinElute Virus Spin Kit (Qiagen/United States) from the specimens and used as the template for genome sequencing. Gene-specific primers were designed based on representative identified reads of the novel parvovirus to cover the partial genomes by PCR amplification and Sanger sequencing. The remaining genomic sequences were determined by genome walking (Takara/Japan), 5’-race system version 2.0 combo (Invitrogen/United States) and 3’-Full RACE Core Set with PrimeScript RTase (Takara/Japan).

Using the complete genomic sequences of the viruses obtained by the ends amplification, we designed specific, semi-nest primers targeting the NS1 gene for PCR to screen for parvoviruses in the fecal samples from crocodiles and intestine samples of tilapia. PCRs were performed using 2 × Taq PCR Mastermix (Tiangen, China), and the first round PCR was primed with outer primers (TParv-F:5′-TAATAAAGCGGGCGTGGAACA-3′ and TParv-R:5′-CTATTACTCAATGCTGCTCGTTCA-3′) and the second round PCR was primed with inner primers (TParv-Fn:5′-ATTACAGAGAAGTGACCGTTTTAG-3′ and TParv-R). 2 μl of the first round PCR product was used as the template for the second round of PCR. The thermal cycling conditions for both PCRs were 94°C for 5 min, followed by 35 cycles of 94°C for 30 s, 57°C for 30 s, 72°C for 45 s, and a final elongation step at 72°C for 10 min. Finally, the PCR products were analyzed on 1.5% agarose gel electrophoresis ultraviolet imaging. Positive samples were determined with 872 bp amplified products. Approximately 30% of the PCR products were randomly selected for Sanger sequencing for confirmation.

The nucleotide sequences of the genomes and the amino acid sequences of the ORFs were deduced by comparing the sequences with those of other parvoviruses. Conserved protein families and domains were predicted using Pfam and InterProScan 5¹. Routine sequence alignments were performed using Clustal Omega, Needle² (Zhu et al., 2016).

MEGA6.0 was used to align nucleotide sequences and deduce amino acid sequences, using the MUSCLE package and default parameters. The phylogenetic tree was constructed with the maximum-likelihood method. The substitution models rtREV with Freqs(+F) model were carried out by model selection function of MEGA6.0 with 1000 bootstrap replicates. Amino acid pairwise alignment between Tilapia parvovirus HMU-HKU and other members of Parvoviridae were performed with NCBI Basic Local Alignment Search Tool. The p-distances were calculated using the “Compute Pariwise Distances” calculation in MEGA6 (Supplementary Table 1).

The near complete genome sequence of tilapia parvovirus HMU-HKU (HMU-HKU-1) and three complete ORFs of NS1 (HMU-HKU-2, HMU-HKU-3, and HMU-HKU-4) were deposited to GenBank under accession numbers MN162688 to MN162691. The Illumina HiSeq2500 sequence data of crocodilian feces samples were deposited into the NCBI sequence reads archive (SRA) under accession number PRJNA543153.

To identify novel viruses in crocodiles in Hainan, we first collected eight fecal samples from eight different crocodiles in a crocodile breeding site in the Lingshui County in the Southern part of Hainan (Supplementary Figure 1). The eight samples were combined into one pool for NGS. A total of 3.7 GB of nucleotide data (38,225,149 valid reads, 100 bp in length) were obtained. The reads that were classified as cellular organisms (bacteria, archaea, and eukaryotes) and reads with no significant similarities to any amino acid amino acid sequence in the NR database were removed. The remaining 197,800 reads were best-matched with viral proteins from the NR database (approximately 0.52% of the total sequence reads). There were 620 Parvovirinae-associated reads in the lane. They were assembled as one contig for 902 bp genome fragment.

Based on the NGS-generated fragment sequences, primers were designed to obtain the intervening portions of the genome. The terminal sequences were then acquired using a 5′- and 3′-RACE method. The nearly complete genome of the novel tilapia parvovirus HMU-HKU-1 was 3566 bp in length (Figure 2), with a G+C content of 43.89%. A putative ORF1 (nt 24–350) was found before the putative NS1. NS1 (nt 269–2125) was 1,857 bp in length and encoded a 618-amino acid non-structural protein. The NS1 in tilapia parvovirus HMU-HKU1 had less than 45.5% amino acid identity with the NS1’s of other known members of Parvoviridae. A putative nucleoprotein (NP) was identified to be overlapping with the NS1 at nt 1428–2039. This putative NP shared only 33.99% amino acid identity with that of simian parvo-like virus 3 (accession number: KT961660). The VP1 ORF was 1,260 bp in length and encoded a putative 419-amino acid capsid protein VP1. The capsid protein VP1 had less than 40.38% amino acid identity with that of other known members of Parvoviridae. Based on these findings, the predicted genome structure of the novel tilapia parvovirus HMU-HKU-1 is more similar to that of members of Chaphamaparvovirus than those of the other members of Parvoviridae.

At the nucleic acid level, the novel tilapia parvovirus HMU-HKU-1 genome sequence was barely comparable with other known parvoviruses in the NCBI nucleotide database. Conserved polyadenylation signals were absent, downstream of the putative NS1 and the VP protein-coding regions (Arthur et al., 2009). In the putative NS1 protein, the following motifs were present: a Walker loop motif, conserved domains for viral replication, a ³⁸⁰GASSSGKS³⁸⁷ [GxxxxGK(T/S)] domain (Walker et al., 1982), a catalytic domain for rolling circle replication, and a ⁴⁵HYHVLV⁵⁰ domain [HUHUUU, U being hydrophobic amino acid] (Ilyina and Koonin, 1992).

Phylogenetic analysis of the novel tilapia parvovirus HMU-HKU-1 was conducted with 62 reference genome sequences of viruses in the family of Parvoviridae available in GenBank (accessed on 1st August 2019). Evolutionary trees were constructed for the complete protein sequences of helicase domain (Figure 3). The topologies of the tree showed that the tilapia parvovirus HMU-HKU-1 formed a relatively independent branch in the newly proposed genus Chaphamaparvovirus in the subfamily Hamaparvovirinae. The Syngnathus scovelli ChVP, another fish parvovirus and the only member of the novel proposed genus Ichthamaparvovirus, clustered with Chaphamaparvovirus in the same branch, and represented a deeper rooted branch in the subfamily Hamaparvovirinae (Pénzes et al., 2019).

Pairwise comparisons were performed for the nucleotide and amino acid sequences of tilapia parvovirus HMU-HKU-1 with other parvoviruses. The results showed that the genome of tilapia parvovirus HMU-HKU-1 shared less than 8% similarity with other members of Parvoviridae at nucleotide level and was barely comparable with other known parvoviruses. At the amino acid level, tilapia parvovirus HMU-HKU-1 also exhibited low level of similarity with other members of Parvoviridae. The NS1 and VP1 of tilapia parvovirus HMU-HKU-1 exhibited ≤45.50 and ≤40.38% amino acid identity, respectively, with those of Chaphamaparvovirus. The NS1 of tilapia parvovirus HMU-HKU-1 exhibited 31.37% amino acid identity with that of the recently reported Syngnathus scovelli ChPV and exhibited the highest identity of 45.50% with that of porcine parvovirus 7 (accession number: KU563733). As the NS1 of tilapia parvovirus HMU-HKU-1 exhibited ≥35% amino acid identity with those of Chaphamaparvovirus, this novel virus likely represented a new species within the new genus Chaphamaparvovirus according to the ICTV’s most recent taxonomic criteria for Parvoviridae classification.

As the eight crocodiles we sampled in Lingshui County were exclusively fed with tilapia, we further investigated the differential prevalence rates of tilapia parvovirus HMU-HKU among crocodiles and tilapia in different regions of Hainan. To do this, we first designed a novel semi-nest PCR assay using specific primers targeting the viral NS1 gene. Next, we collected 50 fecal samples from 50 different crocodiles which were fed exclusively with tilapia, and another 50 fecal samples from 50 different crocodiles which were originally fed with tilapia until 2 weeks before sampling, after which they were exclusively fed with chickens. These two cohorts of crocodiles were bred in two different breeding sites in the Wenchang County at the northeastern part of Hainan. As shown in Supplementary Figure 1 and Table 1, the prevalence rate of tilapia parvovirus HMU-HKU in the crocodiles fed exclusively with tilapia was 86.0% (43/50), which was significantly higher (P < 0.001) than the one of the crocodiles fed with chicken (2/50, 4.0%).

To determine whether the novel tilapia parvovirus HMU-HKU was also found in the tilapia in other regions of Hainan, we then collected additional intestine samples from tilapia bred in other counties/cities in the province. As shown in Supplementary Figure 1, the virus is highly prevalent among tilapia bred in different locations, including Qionghai (18/20, 90.0%), Ding’an (16/20, 80.0%), Dongfang (8/10, 80%), and Haikou (8/20, 40.0%). The positive NS1 amplicon in feces samples from crocodile and tilapia were sequenced. On the basis of phylogenetic analysis and multiple sequence alignments of nucleotide sequences, there were four high identity parvovirus strains; HMU-HKU-1, HMU-HKU-2, HMU-HKU-3, and HMU-HKU-4, with 92.6–94.9% nucleotide identify with each other. The four strains formed an independent cluster with outgroup Chaphamaparvovirus (Figure 4). In multiple amino acid sequence alignments, in the Walker loop motif, [GxxxxGK(T/S)] domain, the first two amino acids are GA in strain 1 and 3, ST in strain 2 and DK in strain 4 (Supplementary Figure 2).