key: cord-0882279-03hrlogv
authors: Ito, Mika; Tsuchiaka, Shinobu; Naoi, Yuki; Otomaru, Konosuke; Sato, Mitsuo; Masuda, Tsuneyuki; Haga, Kei; Oka, Tomoichiro; Yamasato, Hiroshi; Omatsu, Tsutomu; Sugimura, Satoshi; Aoki, Hiroshi; Furuya, Tetsuya; Katayama, Yukie; Oba, Mami; Shirai, Junsuke; Katayama, Kazuhiko; Mizutani, Tetsuya; Nagai, Makoto
title: Whole genome analysis of Japanese bovine toroviruses reveals natural recombination between porcine and bovine toroviruses
date: 2016-03-31
journal: Infection, Genetics and Evolution
DOI: 10.1016/j.meegid.2015.12.013
sha: 45cea3e36f83adaa347b9b323647d0bf0913ad75
doc_id: 882279
cord_uid: 03hrlogv

Abstract Bovine toroviruses (BToVs), belong to the subfamily Toroviridae within the family Coronaviridae, and are pathogens, causing enteric disease in cattle. In Japan, BToVs are distributed throughout the country and cause gastrointestinal infection of calves and cows. In the present study, complete genome sequences of two Japanese BToVs and partial genome sequences of two Japanese BToVs and one porcine torovirus (PToV) from distant regions in Japan were determined and genetic analyses were performed. Pairwise nucleotide comparison and phylogenetic analyses revealed that Japanese BToVs shared high identity with each other and showed high similarities with BToV Breda1 strain in S, M, and HE coding regions. Japanese BToVs showed high similarities with porcine toroviruses in ORF1a, ORF1b, and N coding regions and the 5′ and 3′ untranslated regions, suggestive of a natural recombination event. Recombination analyses mapped the putative recombinant breakpoints to the 3′ ends of the ORF1b and HE regions. These findings suggest that the interspecies recombinant nature of Japanese BToVs resulted in a closer relationship between BToV Breda1 and PToVs.

Toroviruses (ToVs), belonging to the subfamily Toroviridae, are members of the family Coronaviridae, order Nidovirales (Cavanagh and Horzinek, 1993) . ToVs have a single-stranded positive-sense RNA genome of nearly 28 kb, containing two large open reading frames (ORFs) encoding nonstructural proteins, ORF1a and ORF1ab, and four structural proteins, spike glycoprotein (S), membrane glycoprotein (M), hemagglutinin esterase (HE), and nucleocapsid phosphoprotein (N) (Draker et al., 2006; Sun et al., 2014) .

ToVs are identified in various animals and humans and are thought to cause diarrheic and respiratory diseases (Duckmanton et al., 1997; Ito et al., 2009; Kroneman et al., 1998; Uziel et al., 1999; Vanopdenbosch et al., 1991; Woode et al., 1982) . Berne virus (EToV Berne), the first isolated ToV, was isolated from a horse with diarrhea in 1972 in Switzerland (Weiss et al., 1983) . In 1979, a bovine ToV (BToV), named Breda virus, was detected in calves with diarrhea in the United States (US). BToV Breda causes enteric disease in gnotobiotic reared calves (Woode et al., 1982) . BToVs are found throughout the world including North America (Hoet et al., 2002 (Hoet et al., , 2003 , Central America (Pérez et al., 1998) , Europe (Haschek et al., 2006; Koopmans et al., 1991; Matiz et al., 2002) , Asia (Aita et al., 2012; Ito et al., 2010; Park et al., 2008) and South Africa (Penrith and Gerdes, 1992) . BToVs have been detected in 2.9-36.4% of fecal samples obtained from cattle with diarrhea (Duckmanton et al., 1998; Hoet et al., 2003; Ito et al., 2007; Kirisawa et al., 2007; Nogueira et al., 2013; Park et al., 2008) . Porcine ToV (PToV) is also prevalent in piglets worldwide; however their pathogenicity in swine remains unclear (Anbalagan et al., 2014; Shin et al., 2010; Sun et al., 2014) .

So far, only three whole genome sequences of ToV, namely, BToV Breda1 (AY427798) (Draker et al., 2006) and PToV NPL/2014 (KM403390), which was identified in the US in 2014 (Anbalagan et al., 2014) and PToV SH1 (JQ860350), which was identified in China in 2010 (Sun et al., 2014) , are available on the database. To gain more information about genetic diversity, relationship, and evolution of ToVs, we performed whole genome analysis of Japanese BToVs and PToV. Our data showed natural interspecies recombination events of Japanese BToVs, which originated from genetic recombination of BToV Breda1 and PToV strains.

Four BToVs and one PToV were studied. BToV Ishikawa/2010 (BToV Ishi) was isolated using HRT-18-Aich cells (Aita et al., 2012; Kuwabara et al., 2007 ) from a fecal sample of a cow with diarrhea in Ishikawa Prefecture in 2010 (Ito et al., 2012) . BToV Kagoshima/2014 (BToV Kago), BToV Tochigi/2013 (BToV Tochi), and BToV Tokyo/2014 (BToV Tokyo) were detected in the course of metagenomics of fecal samples obtained from 18-, 12-, and 16-day-old calves with diarrhea in Kagoshima, Tochigi, and Tokyo Prefecture in 2014 , respectively. PToV Tottori/2015 was identified in the course of metagenomics in fecal samples collected from a healthy two-month-old pig in 2015 in Tottori Prefecture.

Since viruses could not be isolated from samples using HRT-18-Aich cells by repeating passage thrice, fecal suspensions (20% v/v in sterile phosphate-buffered saline) of BToV Kago, BToV Tochi, BToV Tokyo, and PToV Tottori were used for RNA extraction. Viral RNA was extracted from 0.25 mL supernatant of BToV Ishi culture (10 5.3 TCID 50 /mL) or 0.25 mL fecal suspensions by using TRIzol® LS Reagent (Life Technologies, Carlsbad, CA, USA), followed by treatment of the RNA with DNase I (TaKaRa Bio Inc., Shiga, Japan). cDNA library was constructed using NEBNext® Ultra RNA Library Prep Kit for Illumina version 2.0 (New England Biolabs, Ipswich, MA, USA), as described previously (Nagai et al., 2015) . The libraries obtained were loaded onto a MiSeq cartridge (MiSeq Reagent Kit V2 (300 cycles); Illumina, San Diego, CA, USA) and sequenced using a MiSeq bench-top sequencer (Illumina) with 151 paired-end reads. To obtain complete genome sequence of BToV Ishi and BToV Kago, rapid amplification of cDNA end method (RACE) (5′-Full RACE Core Set and 3′-Full RACE Core Set; TaKaRa Bio, Otsu, Japan) was employed.

The sequence data were collected using the MiSeq Reporter program (Illumina) to generate reads in FASTQ format. Collected reads were trimmed and assembled into contigs by de novo assembly using CLC Genomics Workbench 6.5.1 (CLC bio, Aarhus, Denmark). The sequences were aligned using ClustalW in MEGA5.22 (Tamura et al., 2011) . Pairwise sequence identity calculations were performed using CLC Genomics Workbench 6.5.1 (CLC bio). The phylogenetic tree was constructed by the maximum likelihood method statistically supported by bootstrapping with 1000 replicates by using MEGA5.22. Recombination analysis was performed using SimPlot software v. 3.5.1 (Lole et al., 1999) and the Recombination Detection Program (RDP) v. 4.58 (Martin and Rybicki, 2000; Martin et al., 2005) .

We performed deep sequencing using the Illumina MiSeq sequencing system. The total ToV sequence read counts (percentage of ToV sequence reads: ToV sequence reads/total reads) of BToV Ishi, BToV Kago, BToV Tochi, BToV Tokyo, and PToV Tottori samples were 198,526 (13.2%), 31,685 (2.4%), 492 (0.04%), 1527 (0.4%) and 1805 (0.05%), respectively. An approximately 28-kb contig was obtained from BToV Ishi and BToV Kago samples with sufficient average sequence read depth of 1007 (maximum depth: 2384) and 162 (maximum depth: 320), respectively; however, large contigs were not obtained from BToV Tochi, BToV Tokyo, and PToV Tottori samples. Since nearly complete sequences of contigs were obtained, complete genome length of BToV Ishi and BToV Kago were determined using the 5′ and 3′ RACE method. The complete genome length of BToV Ishi and BToV Kago excluding poly(A) were 28,309 nucleotides (nt) and 28,301 nt, respectively, which were similar to the length of published sequences of ToVs; BToV Breda1 (28,475 nt), PToV NPL/2014 (28,305 nt), and PToV SH1 (28,301 nt). The nucleotide sequence identities among complete genome sequences of BToV Ishi and BToV Breda1, PToV NPL/2014, PToV SH1, and BToV Kago were 82.3%, 86.3%, 85.3%, and 97.6%, respectively. The sequences of BToV Ishi and BToV Kago were deposited in the DNA Data Bank of Japan, DDBJ/EMBL/GenBank database under the accession numbers LC088094 and LC088095, respectively.

Pairwise alignment for comparing nucleotide sequences of BToV Ishi and BToV Kago to other BToV and PToVs was performed using the whole genomic region of the 5′ untranslated region (UTR), ORF1a, ORF1b, S, M, HE, N, and the 3′UTR (Table 1) . BToV Ishi shared high sequence identities (96.8-99.4%) with BToV Kago for all genomic regions analyzed. Japanese BToVs showed high identities with PToVs in 5′UTR, ORF1a, ORF1b, N and 3′UTR (86.9-95.1%) and showed low identities with PToVs in S, M, and HE (70.2-80.3%). On the other hand, Japanese BToVs showed high identities with BToV Breda1 in S, M, and HE (87.9-95.7%) and showed low identities with BToV Breda1 in 5′UTR, ORF1a, ORF1b, N and 3′UTR (68.9-83.0%). These results suggest the occurrence of recombination events between BToV Breda1 and PToVs.

To investigate the recombination events, the complete genomes of BToV Ishi, BToV Kago, PToV NPL/2014, PToV SH1, and BToV Breda1 were aligned using ClustalW program in MEGA5.22 and standard similarity plot analysis was performed using SimPlot software v. 3.5.1 with BToV Ishi (Fig. 1) and BToV Kagoshima/2014 sequences as separate queries. Both SimPlot graph indicated that the sequences of Japanese BToV had high similarity with those of PToVs except in the 3′ end of ORF1b, S, M, and most of the HE coding regions, which had high similarity with that of BToV Breda1. The bootscanning analysis using RDP v. 4.58 was performed to identify the presumed recombinant breakpoints. Schema of ToV genome structure is shown in Fig. 2. A. Beginning and end breakpoint positions were mapped near the 3′ of ORF1b and HE regions, respectively (Fig. 2B) . At the starting putative recombination breakpoint, 25-nt sequences of Japanese BToVs corresponding to nt 20,138-20,162 of BToV Ishi showed low identities (72-76%) with that of BToV Breda1, whereas 25 nt sequences of Japanese BToVs corresponding to nt 20,161-20,185 showed high identities (100%) with that of BToV Breda1. At the end putative recombination breakpoint, 25 nt sequences of Japanese BToV corresponding to nt 27,404-27,428 of BToV Ishi showed low identities (84-88%) with those of PToVs; however, the 25-nt Japanese BToV sequences corresponding to nt 27,418 to 27,442 showed high identities (96-100%) with those of PToVs (Fig. 2C) .

For further confirmation of the recombination event, phylogenetic trees of nine genomic regions were constructed using nucleotide sequences of BToVs and PToV determined in this study together with ToV and PToV genome sequences available in GenBank. Although the genome sequences of BToV Tochi, BToV Tokyo, and PToV Tottori had gaps owing to insufficient sequence read counts, sequence regions with no gap and corresponding to BToV Ishi were used: 5′UTR In addition, Japanese BToVs in this study together with other BToV strains from Japan and other countries selected from GenBank, clustered with BToV Breda1 but not with PToVs in the S, M, and HE phylogenetic trees (Fig. 3D-G) . Ito et al. classified Japanese BToVs using 5′ portion of S coding region into three genotypes, Clusters 1-3 (Ito et al., 2007 (Ito et al., , 2009 (Ito et al., , and 2010 . In the phylogenetic tree, the 5′ portion of S coding region of BToVs in this study showed that BToV Ishi and BToV Tochi clustered with Cluster-1 and Cluster-2 strains, respectively, while BToV Kago and BToV Tokyo clustered with Cluster-3 strain, though enough bootstrap support could not be obtained. (Fig. 3D) . BToV Tochi separately clustered with BToV Ishi, BToV Kago, and BToV Tokyo in the S (21,693-25,548), M, and HE phylogenetic trees (Fig. 3E-G) .

In the present study, we performed whole genome sequencing of Japanese BToVs and PToV by deep sequencing using supernatants of virus cultures and fecal suspensions. We could obtain a sufficient number of sequence reads and contigs of nearly 28 kb from 0.25 ml of 10 5.3 TCID 50 /mL BToV Ishi cell culture supernatant. Only one sample among three 0.25-mL fecal suspension samples gave a sufficient number of sequence reads and large contigs, though the number of sequence reads from this sample was lower than that from the viral culture supernatant. Unfortunately, we could not isolate BToV from fecal sample of BToV Kago by using HRT-18-Aich cells, which are useful for BToV isolation (Aita et al., 2012; Ito et al., 2012; Kuwabara et al., 2007) , possibly because the samples could not be preserved at −20°C at the veterinary clinic. Although numerous sequence reads originating from bacterial species and hosts might interfere and reduce the number of sequence reads from fecal suspensions, the whole genome sequence of BToV could be obtained from this sample by using our deepsequencing strategy.

The results of sequence comparison through pairwise alignment, similarity plot, bootscanning, and phylogenetic tree analyses of Japanese BToVs with BToVs and PToVs sequences from database were indicative of interspecies recombination between BToV Breda1 and PToV at the 3′ end of ORF1b and the end of HE coding region. Smits et al. provided evidence for recombination of ToV sequences between BToV Breda1 and PToV in the HE coding region (Smits et al., 2003) . Italian and Dutch BToVs identified before the year 2000 were reported to possess a single recombination breakpoint located approximately 1100 nucleotides from the ATG initiation codon of the HE coding region. The putative recombinant breakpoint in the HE coding region of Japanese BToVs was quite identical to that of European BToVs. In the phylogenetic tree analyses, Japanese BToVs were found to be closely related to European BToVs in S, M, HE, and N coding regions. These findings suggest that Japanese BToVs and European BToVs have a common ancestor. However, owing to the lack of sequence information of 5′UTR, ORF1a and ORF 1b of European BToVs, we could not further compare Japanese BToVs with European BToVs. The HE coding region of ToV shares sequence homology to that of coronavirus and influenza C virus, however it is suggested that HE region of these viruses were introduced through independent heterologous recombination events (Smits et al., 2005) . Homologous recombination of HE coding region was also reported in PToVs (Cong et al., 2013; Smits et al., 2003) , suggesting that the HE coding region may be a hot spot for recombination events.

To our present knowledge, this is the first report to present the identification of a natural recombination event between BToV Breda1 and PToV at ORF1b. Approximately 76% of the genomes of BToV Ishi and BToV Kago are similar to that of PToV. Natural recombination events are common in coronaviruses, resulting in evolution and emergence of new variants that are pathogenic to poultry, cats, pigs, and humans (Graham and Baric, 2010; Herrewegh et al., 1998; Hewson et al., 2014; Lau et al., 2015; Song et al., 2013; Terada et al., 2014; Tian et al., 2014; Wang et al., 2015; T. Zhang et al., 2015; Y. Zhang et al., 2015; Zhao et al., 2013) . Co-infection of certain cells in bovine or swine with BToV Breda1 type virus and PToV is a pre-requisite for the generation of new recombinant viruses by double-template switch. The generated new recombinant virus, which has tropism to bovine cells, might adapt to bovine hosts and spread among cattle.

Although the four Japanese BToVs evaluated in this study were collected from distinct regions, the sequences of genomic regions apart from the S, M, and HE coding regions were closely related to each other and showed high identity (97.3-99.4%). Based on the sequences of the S, M, and HE coding regions, these BToVs were divided into at least two groups (Fig. 3D-G) . The four Japanese BToVs showed genetic diversity, they all possessed a recombinant genome, suggesting that Japanese BToVs evolved after recombination events after or before being introduced in Japan. The Japanese and foreign BToV strains retrieved from the database clustered closely together in the phylogenetic tree, suggesting that all BToVs investigated in this study, except BToV Breda1, underwent recombination at the same genomic regions. Strains genetically closely related to Japanese BToVs are distributed worldwide.

In conclusion, we identified a natural recombination event between BToV Breda1 and PToV of Japanese BToVs. The four Japanese BToVs possessed two recombination breakpoints mapped to the 3′ end of ORF1b coding region and to the 3′ end of HE coding region and nearly 76% of the genome is similar to the PToV genome. Homologous recombination events of the virus genome are crucial for their evolution, and are significant factors for changes in host range and virulence of viruses. Our current data may provide important evidence to evaluate the epidemiological basis of ToV in cattle and swine population and to understand the mechanisms underlying the evolution of BToV. 

Characterization of epidemic diarrhea outbreaks associated with bovine torovirus in adult cows

Genome sequence of torovirus identified from a pig with porcine epidemic diarrhea virus from the United States

Genus torovirus assigned to the coronaviridae

Evolution and homologous recombination of the hemagglutinin-esterase gene sequences from porcine torovirus

The complete sequence of the bovine torovirus genome

Characterization of torovirus from human fecal specimens

Detection of bovine torovirus in fecal specimens of calves with diarrhea from Ontario farms

Recombination, reservoirs, and the modular spike: mechanisms of coronavirus cross-species transmission

Detection of bovine torovirus in neonatal calf diarrhoea in Lower Austria and Styria (Austria)

Feline coronavirus type II strains 79-1683 and 79-1146 originate from a double recombination between feline coronavirus type I and canine coronavirus

Evaluation of a novel strain of infectious bronchitis virus emerged as a result of spike gene recombination between two highly diverged parent strains

Enteric and nasal shedding of bovine torovirus (Breda virus) in feedlot cattle

Detection of bovine torovirus and other enteric pathogens in feces from diarrhea cases in cattle

Epidemiological analysis of bovine torovirus in Japan

Detection and characterization of bovine torovirus from the respiratory tract in Japanese cattle

Genetic and antigenic characterization of newly isolated bovine toroviruses from Japanese cattle

Outbreak of bovine torovirus diarrhea in Ishikawa Prefecture

Detection of bovine torovirus in fecal specimens of calves with diarrhea in Japan

Association of diarrhea in cattle with torovirus infections on farms

Identification and characterization of a porcine torovirus

First isolation of cytopathogenic bovine torovirus in cell culture from a calf with diarrhea

Severe acute respiratory syndrome (SARS) coronavirus ORF8 protein is acquired from SARS-related coronavirus from greater horseshoe bats through recombination

Fulllength human immunodeficiency virus type 1 genomes from subtype C-infected seroconverters in India, with evidence of intersubtype recombination

RDP: detection of recombination amongst aligned sequences

A modified bootscan algorithm for automated identification of recombinant sequences and recombination breakpoints

Torovirus detection in faecal specimens of calves and pigs in Hungary: short communication

Full genome analysis of bovine astrovirus from fecal samples of cattle in Japan: identification of possible interspecies transmission of bovine astrovirus

First detection and molecular diversity of Brazilian bovine torovirus (BToV) strains from young and adult cattle

Molecular epidemiology of bovine toroviruses circulating in South Korea

Breda virus-like particles in pigs in South Africa

Infectious agents associated with diarrhoea of calves in the canton of

Detection and molecular characterization of porcine toroviruses in Korea

Phylogenetic and evolutionary relationships among torovirus field variants: evidence for multiple intertypic recombination events

Nidovirus sialate-O-acetylesterases: evolution and substrate specificity of coronaviral and toroviral receptor-destroying enzymes

Sequencing, phylogenetic analysis, and potential recombination events of infectious bronchitis viruses isolated in Korea

Molecular characterization and phylogenetic analysis of the genome of porcine torovirus

MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods

Emergence of pathogenic coronaviruses in cats by homologous recombination between feline and canine coronaviruses

Evidence of recombinant strains of porcine epidemic diarrhea virus, United States

Torovirus gastroenteritis presenting as acute abdomen

Breda virus associated with respiratory disease in calves

Origin and possible genetic recombination of the Middle East respiratory syndrome coronavirus from the first imported case in China: phylogenetics and coalescence analysis

Purification and partial characterization of a new enveloped RNA virus (Berne virus)

Studies with an unclassified virus isolated from diarrheic calves

Serotype shift of a 793/B genotype infectious bronchitis coronavirus by natural recombination

Genotype shift in human coronavirus OC43 and emergence of a novel genotype by natural recombination

Analysis of a QX-like avian infectious bronchitis virus genome identified recombination in the region containing the ORF 5a, ORF 5b, and nucleocapsid protein gene sequences

This work was supported by the Grants from the Ministry of Health, Labour and Welfare of Japan and JSPS KAKENHI Grant Number 15K07718.