Molecular detection and genomic characteristics of bovine kobuvirus from dairy calves in China.

In this study, 96 diarrheic and 77 non-diarrheic fecal samples from dairy calves were collected from 14 dairy farms in 4 provinces to investigate the molecular prevalence and genomic characteristics of Bovine Kobuvirus (BKoV) in China. The results showed that the BKoV positive rate for the diarrheic feces (35.42%) was significantly higher than that for the non-diarrheic feces (11.69%, p < 0.001). Interestingly, three potential novel VP1 lineages were identified from 15 complete VP1 sequences, and a unique triple nucleotide insertion which can result in an aa insertion, was first observed in the 11/12 VP0 fragments with 660 bp long in this study, compared with known BKoV VP0 sequences. Moreover, the first Chinese BKoV genome was successfully obtained from a diarrheic fecal sample, named CHZ/CHINA. The open reading frame (ORF) of the genome from strain CHZ/China shares 87.4%-88.3% nucleotide (nt) and 93.7%-96.4% amino acid (aa) identity, compared with the three known genomes of BKoV. Interestingly, phylogenetic tree based on aa sequences of these genomes showed that CHZ/CHINA was clustered into an independent branch, suggesting the strain may represent a novel BKoV strain. The findings contribute to better understanding the molecular characteristics and evolution of BKoV.

Bovine kobuvirus (BKoV) is a member of Aichivirus B, and another member in Aichivirus B is the sheep kobuvirus (Khamrin et al., 2014; Reuter et al., 2010). Since BKoV was first identified in Japan in 2003 (Yamashita et al., 2003), this virus has been detected in bovine with and without diarrhea symptomatology in 9 countries (Liu et al., 2013). It had been suggested that BKoV may be associated with diarrhea in calves (Candido et al., 2017), but the pathogenicity of BKoV still needs to be determined. Recently, BKoV was also detected in the spinal fluid from the brain of an 11 day-old calf where the animal had a history of diarrhea and neurological disease (Moreira et al., 2017), which indicates that this virus can cause systemic infections.

The BKoV genome is approximately 8.2–8.4 kb long and has the typical kobuvirus genome organization comprising a leader (L) protein, followed by structural (capsid proteins VP0, VP3, and VP1) and nonstructural (2A–2C and 3A–3D) proteins. In BKoV, VP1 is the most variable immune determinant protein (Yamashita et al., 2003), making it appropriate for genetic typing (Oh et al., 2006; Pham et al., 2008; Shi et al., 2013). The function of VP0 and VP3 of BKoV protein remains unclear, but VP0 (residues 50–63) in Aichi virus (AIV) may be involved in cellular receptor recognition (Zhu et al., 2016). The VP3 viral protein in porcine kobuvirus (PKV) may play a immune evasion role via the IFN signaling pathway (Peng et al., 2017). The aim of this study was to further investigate the molecular prevalence and genomic characteristics of BKoV in China.

From January 2018 to April 2018, a total of 96 diarrheic and 77 non-diarrheic were collected from 14 dairy farms across four the Chinese provinces of Liaoning (three farms), Henan (three farms), Shandong (three farms) and Shanxi (five farms). The ages of the tested calves ranged from 2 days old to 4 months old. All samples were shipped on ice and stored at −80 °C in sterile 50-ml centrifuge tubes. The fecal samples were fully resuspended in phosphate-buffered saline (1:5) and centrifuged at 10,000 ×g for 10 min, followed by filtration through a 0.45-μm filter. Viral RNA was extracted from 300 μl of each fecal suspension using RNAios Plus (TaKaRa Bio Inc., Japan) according to the manufacturer's instructions. The cDNA was synthesized using the PrimeScript™ RT Reagent kit according to the manufacturer's instructions (TaKaRa Bio Inc.) and then stored at −20 °C. BKoV was detected by an RT-PCR assay targeting a 631 bp fragment of the 3D gene according to previous report (Jeoung et al., 2011). To screen for the presence of co-infections with bovine rotavirus (BRV), bovine coronavirus (BCoV), and bovine viral diarrhea virus (BVDV), all the BKoV-positive diarrheic samples were subjected to specific RT-PCR assays for these viruses (Guo et al., 2019; Zheng et al., 2014).

The complete VP1 (801 nt) sequences were amplified from the BKoV-positive samples according to previous report (Liu et al., 2013). A pair of primers was designed based on known BKoV VP0 sequences, located at positions 1012–1671 in CHZ/China genome sequence. Moreover, 13 pairs of primers (Table S1) were used to amplify the genome sequence of the BKoV CHZ/China strain. All PCR products were purified using the Omega Gel kit (Omega), cloned into the pMD19-T simple vector (TaKaRa Bio Inc.), and then sequenced (Sangon Biotech) in both directions. Sequences were assembled using SeqMan software (version 7.0; DNASTAR Inc., WI, USA). Single open reading frames (ORFs) were identified using the online ORF finder at https://www.ncbi.nlm.nih.gov/orffinder/. Nt and deduced aa sequence homologies were determined using the MegAlign program in DNASTAR 7.0 software (DNASTAR Inc.). Phylogenetic trees were constructed using the maximum likelihood method with the Jukes-Cantor model, 1000 bootstrap replicates and default parameters in MEGA 7. The BKoV detection rates for the diarrhea and healthy samples were calculated and compared using the Epi Info statistical program (version 7.0), and values of p < 0.01 were regarded as statistically significant.

Among the 96 diarrheic samples, 34 were detected as BKoV positive (35.42%), and 11 out of 77 non-diarrheic samples were detected as BKoV positive (11.69%), as shown in Table 1 . The BKoV detection rate for the diarrheic samples was significantly higher than that for the non-diarrheic samples (p < 0.001), suggesting that the virus may be associated with diarrhea in calves. However, of the 34 BKoV positive diarrheic samples, 28 were confirmed to be co-infection positive for other viruses, as shown in Table S2. Therefore, further investigations are needed to better understand the role of BKoV infection in cattle with diarrhea.It may be that calf challenge experiments will be the only way to answer the question as to whether BKoV causes diarrhea.

Table 1

Sample collection and BKoV detection information.

Region	The number of diarrhea samples and test results		The number of non-diarrhea samples and test results		Number of dairy farms
	Number of samples	Positive rate	Number of samples	Positive rate
Liaoning	18	83.33%(15/18)	0	0	3
Henan	27	44.44% (12/27)	25	32.00% (8/25)	3
Shandong	11	63.64% (7/11)	12	8.33% (1/12)	3
Shanxi	40	0	40	0	5
Total	96	35.42% (34/96)	77	11.69% (9/77)	14

15 complete VP1 sequences were successfully obtained for fecal samples (GenBank accession numbers MK080238–MK080252), which share 74.2%–100.0% nt sequence identity (80.1%–100.0% aa identity) with each other, and share 72.9%–94.9% nt sequence identity (80.5%–97.8% aa identity) with VP1 sequences in the GenBank database. A phylogenetic tree analysis based on all 103 available complete sequences for BKoV VP1 and the 15 VP1 sequences from the present study shows that all the VP1 sequences fall into 8 distinct branches (Fig. 1 ). 7/15 VP1 sequences from our study were clustered into the known lineage 2. Interestingly, the remaining 8/15 VP1 sequences from our study were clustered into 3 independent branches distinct from the 5 known lineages, suggesting that these branches may represent 3 novel lineages. The potential lineages 6, 7, 8 share 73.5%–76.7%, 74.0%–78.8% and 74.0%–85.9% nt sequence identity with all complete BKoV VP1, including the VP1 sequences from the present study, respectively. Furthermore, potential lineage 6 contains only 1 strain of B5/HN/China strains; potential lineage 7 contains 2 strains of 9/LN/China and 14/LN/China, which share 99.8% nt sequence identity with each other; potential lineage 8 contains 5 strains of 3/LN/China, 7/LN/China, A2/SD/China, A3/SD/China and B21/HN/China strains, which share 90.9%–100% nt sequence identity with each other. Further analysis found that the strain in potential lineage 6 has 5 out of 267 unique aa mutations (P35S, T40S, R98K, T114A, and Y254L), 2 strains in potential lineage 7 share 2 out of 267 unique aa mutations (R43Q and I107V), and 5 strains potential lineage 8 share 1 out of 267 unique aa mutations (E18D), compared to all known complete VP1 sequences of BKoV. A previous study reported that a strain with a nt identity <85.1% should be considered a novel lineage (Liu et al., 2013), but the value of 85.1% may be biased, due to the limitation of sequence numbers of BKoV VP1 at that time. Currently, the precise biological function of BKoV VP1 remains unclear. However, information from AIV shows that VP1 is the structural protein to the enteric receptor recognition and may be involved in viral pathogenesis (Zhu et al., 2016). As VP1 is the most variable immune determinant protein in kobuvirus, VP1 in BKoV, AIV and PKV had been divided into different lineages (Liu et al., 2013; Oh et al., 2006; Pham et al., 2008; Shi et al., 2013). Further investigations are needed to better understanding the antigenicity of different VP1 lineages.

Fig. 1

Maximum likelihood phylogenetic tree based on complete VP1 nt sequences. Black circles denote isolates from the present study and hollow circles denote isolates from a previous Chinese study (Liu et al., 2013). Bootstrap values based on 1000 replicates are shown on the nodes.

12 VP0 fragments were successfully obtained from fecal samples (GenBank accession numbers MK080253–MK080264). The 12 sequences were found to share 78.5%–100% nt and 84.1%–100% aa sequence identity with each other, and 80.2%–84.7% nt and 85.0%–91.8% aa sequence identity with the only 3 BKoV VP0 sequences in the GenBank database. A phylogenetic tree based on the three known VP0 aa sequences and the 12 sequences from this study revealed that 11 of the BKoV strains from the present study clustered on an independent branch. The remaining 1 strain from this study and 3 known strains clustered on an independent branch (Fig. 2 ). Further analysis found that the 11/12 sequences in an independent branch shared a unique aa insertion in two forms (S63 or N63) in the VP0 region, showing that the VP0 type was the high frequency in Chinese BKoV. Interestingly, a unique aa insertion within the VP0 region also found in AIV strains (Pham et al., 2008). Hence, whether the variation in VP0 is a potential evolutionary characteristic in kobuvirus needs further to be investigated.

Fig. 2

Maximum likelihood phylogenetic tree based on 220 aa sequence alignment of the VP0 protein from BKoV strains. Black circles denote isolates from the present study. Bootstrap values based on 1000 replicates are shown on the nodes.

Until now, there are 3 complete BKoV genomes (U-1, SC1, EGY-1) in the GenBank database, contributing to understanding the genetic characteristics of BKoV. In this study, we added to a nearly complete BKoV strain (GenBank accession No. MK080265) genome of 7907 nt in length which contains the 7395 bp complete ORF, which is the first BKoV genome from China. Compared with 3 known BKoV, with the exception that the CHZ/CHINA strain's VP0 sequence is 3 nt longer, the lengths of the other CHZ/CHINA genes are identical to those of the other three genomes (Table S3) and shares 87.4%–88.3% nt and 93.7%–96.4% aa identity. Further phylogenetic analysis based on genomic sequences revealed that CHZ/CHINA clusters on an independent branch, with the three VP0, VP3, VP1 protein aa sequences generating the same result (Fig. 3 ), showing that CHZ/CHINA displays a larger genetic distance from the other three genomes and indicating that CHZ/CHINA may represent a novel BKoV strain. Moreover, the most significant difference between CHZ/CHINA and other BKoV strains is the VP0 protein. And VP0 from the CHZ/CHINA strain contains 19 unique aa mutations, and a unique triple nt insertion which can result in an aa insertion. The function of BKoV VP0 remains unclear. But, both VP0 and VP1 in AIV may be involved in cellular receptor recognition (Zhu et al., 2016) and viral pathogenesis (Adzhubei et al., 2013). Thus, it is worth studying the functional effects of this unique VP0 aa mutation in BKoV strains.

Fig. 3

Maximum likelihood phylogenetic tree based on sequence alignments for the ORF, VP0, VP3, and VP1 protein encoding domains from BKoV strains. Black circles denote isolates from this study. Bootstrap values based on 1000 replicates are shown on the nodes.

In conclusion, the results of this study showed that three potential novel VP1 lineages in BKoV were identified and a unique BKoV VP0 sequence type was found in diarrheic feces. The first nearly complete BKoV genome was obtained and phylogenetic analysis shows that this strain may represent a novel BKoV strain. These data contribute to further understanding of the molecular characteristics and genetic evolution of BKoV.

Depositories

More information about sequences is in the GenBank database:MK080201–MK080265.

Ethical statement

This study did not involve animal experiments besides the fecal sampling of diarrhea calves that visited farm for clinical treatment.

Declaration of Competing Interest

None.