Introduction
The black queen cell virus (BQCV) causes a disease that targets queen larvae. The main symptoms of infection include: pale yellow skin, gradual hardening of the epidermis into cysts, dark pupae body colour, sudden pupae death and queen bodies becoming grey and then black. Even the inner honeycomb wall can become brown or black (Bailey & Woods, Reference Bailey and Woods1977). BQCV affects not only the morphology of bees but also their physiology and behaviour, eventually leading to the death of the entire colony. BQCV was first isolated from queen prepupae and pupae, which were found dead in their cells. Recent investigations have shown that BQCV is a widely prevalent virus that targets multiple species of bees, including non-Apis hymenopteran pollinators. Especially in western countries and in Apis mellifera, BQCV is the most prevalent of all viruses affecting honeybees (Tentcheva et al., Reference Tentcheva, Gauthier, Zappulla, Dainat, Cousserans, Colin and Bergoin2004; Berényi et al., Reference Berényi, Bakonyi, Derakhshifar, Köglberger and Nowotny2006; Baker & Schroeder, Reference Baker and Schroeder2008; Kajobe et al., Reference Kajobe, Marris, Budge, Laurenson, Cordoni, Jones, Wilkins, Cuthbertson and Brown2010). BQCV mainly caused the death of queen larvae, yet the infected worker bees were apparently healthy (Bailey, Reference Bailey1982), and it acted as a kind of chronic infection. In recent years, it has been reported that BQCV had spread throughout Europe, Asia, Oceania, the Americas and Africa (Antúneza et al., Reference Antúnez, D'Alessandro, Corbella, Ramallo and Zunino2006; Bereényi et al., Reference Berényi, Bakonyi, Derakhshifar, Köglberger and Nowotny2006; Chantawannakul et al., Reference Chantawannakul, Ward, Boonham and Brown2006; Chen & Siede, Reference Chen and Siede2007; Higes et al., Reference Higes, Esperón and Sánchez-Vizcaíno2007; Baker & Schroeder, Reference Baker and Schroeder2008; Forgach et al., Reference Forgach, Bakonyi, Tapaszti, Nowotny and Rusvai2008; Haddad et al., Reference Haddad, Brake and Migdadi2008; Kukielka et al., Reference Kukielka, Perez, Higes, Bulboa and Sánchez-Vizcaíno2008; Nielsen et al., Reference Nielsen, Nicolaisen and Kryger2008; Teixeira et al., Reference Teixeira, Chen, Message, Pettis and Evans2008; Sanpa & Chantawannakul, Reference Sanpa and Chantawannakul2009; Welch et al., Reference Welch, Drummond, Tewari, Averill and Burand2009; Kajobe et al., Reference Kajobe, Marris, Budge, Laurenson, Cordoni, Jones, Wilkins, Cuthbertson and Brown2010; Hong et al., Reference Hong, Woo, Choi, Han, Kim, Kim, Han, Lee, Lee and Byeon2011; Kojima et al., Reference Kojima, Toki, Morimoto, Yoshiyama, Kimura and Kadowaki2011; Runckel et al., Reference Runckel, Flenniken, Engel, Ruby, Ganem, Andino and DeRisi2011). BQCV primarily infects A. mellifera, a species of bee imported and bred in most provinces in China. Observing the symptoms of colonies at the beginning of BQCV infection is difficult, as symptoms were most frequently observed after the colony had died, so BQCV is a serious virus to honeybees. The effects of this virus have largely been overlooked and China also lacks an effective preventative and controlling measure.
To reduce the large number of missing A. mellifera in China, a study assessed the genetic characteristics and variability of BQCV from the molecular level, laying the theoretical foundation. Thus, further in-depth study of the pathogenicity of BQCV was done, and preventative measures can be pursued as soon as possible. BQCV is classified into the Cripavirus genus within the Dicistroviridae family by Mayo (Reference Mayo2002). This classification has been used in the present study. Pupae were found to contain large numbers of isometric virus particles 30 nm in diameter. They had no capsule membrane and had a positively icosahedron spherical structure (Leat et al., Reference Leat, Ball, Govan and Davison2000). The genomic RNA could be directly used as a template mRNA to induce viral protein translation and consisted of two open reading frames (ORFs): the 5′-proximal ORF (ORF1) encoding putative replicase polyprotein and the 3′-proximal ORF (ORF2) encoding capsid polyprotein. ORF1 encoded the putative helicase, 3C-like cysteine protease and RNA-dependent RNA polymerase (RdRp) enzymes. The nucleotide sequence of the RdRp was analysed to determine its classification. There was a non-coding region between the two ORFs, an internal ribosome entry site (IRES), which facilitated the translation initiation of ORF2 and helped determine the precise nucleotide location of ORF2. ORF2 encoded four capsid proteins (CPs), with molecular weights determined to be 34, 32, 29 and 6 kDa, respectively.
This study was the first to isolate and report the whole genomic sequence of the BQCV China-JL1 strain and compare the sequence diversity of BQCV genotypes in different countries. It was found that the BQCV China-JL1 strain had a greater difference in dividing the nucleotide position ORF than the other six strains commonly found in the GenBank (South Africa: AF183905, South Korea: JX149531, Hungary 10: EF517515, Poland 4: EF517519, Poland 5: EF517520, Poland 6: EF517521). As for the six strains, they all contained two relatively large ORFs, while apart from these two ORFs, the BQCV China-JL1 strain also had a small ORF, temporarily called ORF3, located between the two ORFs. To determine whether ORF3 could be characterized as ORF1 or ORF2, or whether it was a separate ORF indeed, BLAST analysis was employed, as well as DNASTAR software for ORF3 analysis. The amino acid sequences corresponding to the positions of the other six strains were compared. The results consistently showed that ORF3 was contained in ORF1 and encoded the RdRp.
Results
Phylogenetic analysis of BQCV China-JL1 strain complete genome sequences
Analysis of the BQCV China-JL1 strain sequence revealed a polyadenylated genome of 8358 nt, not including the poly (A) tail (fig. 1). The genome contained a high proportion of A and U nucleotides, composed of 29.41% A, 30.51% U, 18.26% C and 21.82% G. Using DNASTAR software, the ‘Find ORF’ function was applied to the six whole genome sequences. The results showed that all six genomic sequences contained two large ORFs. The position of the first BQCV China-JL1 strain may form ORF1, which is characterized by 546–4676 nt, length 4131 nt and an ORF2 length 2454 nt (5750–8203 nt). Comparisons of the other six genome sequences revealed that the BQCV China-JL1 strain ORF1 was about 753–837 nucleotides shorter than theirs, but between the two largest ORFs of the BQCV China-JL1 strain, located at 4891–5433 nt position, there was a length of small ORF 543 nt, temporarily called ORF3. The nucleotide sequence of ORF3 was analysed by DNASTAR software using BLAST and comparative analysis. The purpose was to determine whether ORF3 was vested in ORF1 or ORF2 or indeed a single ORF.
Fig. 1. Comparison between BQCV China-JL1 strain (KP119603) and BQCV South Africa (AF183905) genome sequence. (a) A schematic diagram was created of the BQCV China-JL1 strain genome structure. It shows the nucleotide position and sizes of ORF1, ORF2 and ORF3. Regions of the genome proposed to encode a helicase, protease and RdRp are indicated. Right square indicate the nucleotide position of the strain IRES. The CUU initiation codon of the IRES in ORF2 of BQCV China-JL1 strain is shown. (b) A schematic diagram of the BQCV (SA) genome is shown.
The Chinese genotype showed a high level of similarity (93%) with the Hungary 10 genotype and a low level of similarity (86%) with the South African reference genotype. It was 92% identical to the South Korean genotype, 88% identical to the Poland 4 genotype, and 89% identical to the Poland 5 and Poland 6 genotypes, respectively. A detailed analysis of the different genomic regions, however, revealed an interesting variance in sequence identity. In the complete genome phylogenetic tree, all the BQCV genomes were separated from the deformed wing virus (DWV) complete genome sequence. Within the BQCV genomes, the China-JL1, Hungary 10 and South Korean genotypes formed a subgroup, and all the European, Poland 4–6 and South African genotypes formed a second subgroup (fig. 2). Within the first subgroup, the genotypes further divided into two branches. One branch contained the South Korean reference genotype, but the second branch divided into two further branches, Chinese BQCV and Hungary 10 genotypes.
Fig. 2. Phylogenetic tree of the complete genome sequences of Chinese BQCV-JL1 strain. The phylogenetic tree depicts Chinese, Korean, South African and central European BQCV isolates. The nucleotide sequences of the complete genomes of BQCV isolates from the investigated Chinese isolate, a South Korean strain, a South African reference strain and four central European genotypes (Poland 4–6 and Hungary 10) were aligned using the ClustalW programme. The phylogenetic tree was constructed with the DNASTAR software. The full genome sequence of DWV isolate was used as the out-group. The accession number of each viral sequence is indicated in parentheses.
Comparison of the variability between different BQCV genomic regions ORF1
BQCV particles contained a single stranded RNA genome of 8550 nt, consisting of two ORFs: the 5′-proximal ORF (ORF1) encoding a putative replicase polyprotein, and the 3′-proximal ORF (ORF2) encoding a capsid polyprotein. ORF1 encoded the putative helicase, 3C-like cysteine protease, and RdRp enzymes. ORF2 encoded four CPs(CP1–CP4) (Leat et al., Reference Leat, Ball, Govan and Davison2000). The sequence of each fragment was compared with the Chinese BQCV strain using DNASTAR to determine the location and highest identity, and then BLAST was used to analyse these sequences.
Regions encoding the helicase protein
The region located at 1981–2368 nt (388 nt) encoded the helicase protein of the South African reference strain, while the position of the Chinese strain in this study was between 1788 and 2149 nt. By BLAST-searching the GenBank database with this sequence, we obtained 22 sequences with a homology ranging from 86 to 94%. Among them, only one strain was derived from South Korea (South Korea:JX149531), and the others were all from Central European countries (Austria, Hungary and Poland). It was then found that the sequences encoding the helicase protein of European genotypes in GenBank corresponding to the position from the Chinese strain began at the 186th nt. This fragment size was 560 nt, so it was deduced that the region encoding the helicase protein of the Chinese strain is as big as that of the European genotypes, also 560 nt and similarly located between 1973 and 2532 nt.
Regions encoding the 3C-like cysteine protease
The regions encoding the 3C-like cysteine protease were analysed. The encoding sequence of 3C-like cysteine protease in China-JL1 located at nt 3171–3705. Therefore, this sequence was selected and compared with the Genbank database using BLAST. The results showed that except for the strain Poland4 (EF517519), the 3C-like cysteine protease encoding sequence of China-JL1 had a homology ranging from 84 to 93% with the other five complete genome sequences published in the Genebank. This result suggested that by now, yet no researchers have submitted the encoding sequence of the 3C-like cystein protease in BQCV to GenBank. Comparing the amino acid sequence of the China-JL1 to that of the other five strains, the China-JL1sequence (containing 178 amino acids (aa)) had a very high homology of 98.3–99.4% with other sequences, with the different aa occurring at site 39, 80, 117 and 156, indicating the high conservation of this sequence.
Regions encoding the RdRp enzymes
The ORF3 of the Chinese strain had a 543 nt length (nt 4891–5433), which was translated into 181 aa. The DNASTAR software MegAlign programme was used to compare the reciprocal 181 aa of ORF1 of each complete genome sequence of six strains to the ORF3 of the Chinese strain. The aa sequence achieved consistency between 97.3 and 98.4% (fig. 3). There were subtle differences in the sequence at the positions of 10, 77, 80, 93, 106, 117, 119, 143 and 177 nt (fig. 4). This part was located at the position of ORF1 in the other six complete genome sequences, so the segment sequence was further identified as part of ORF1.
Fig. 3. Homology of ORF3 of the Chinese BQCV-JL1 strain. The DNASTAR software MegAlign programme analysis of the ORF3 region is shown, including China and the other six strains of the complete genome sequence of the amino acid sequence identity. The ORF3 of the Chinese strain had a 543 nt length (nt 4891–5433), which was translated into 181 amino acids (aa). The DNASTAR software MegAlign programme was used to compare the reciprocal 181 aa of ORF1 of each complete genome sequence of 6 strains to the ORF3 of the Chinese strain. The amino acid sequence achieved consistency between 97.3–98.4%.
Fig. 4. Amino acid sequence analysis of ORF3 of the Chinese BQCV-JL1 strain. The DNASTAR software MegAlign programme analysis of the ORF3 region is shown, including China and the other six strains of the complete genome sequence of the amino acid sequence identity.
The South African BQCV strain had an RdRp located between 4606 and 5410 nt, while the position of the BQCV China-JL1 strain was 4413–5217 nt (805 nt). All the strains had similarities to this 805 nt fragment with 88–93% consistency. The Chinese genotype showed a high level of similarity (93%) to the Hungary 10 genotype. It was 91% identical to the Poland 5, Poland 6 and South Korea genotypes and 90% identical to the Poland 4 and the South African reference genotype. The Chinese strains could be divided into a sequence of 5′- and 3′-ends of two parts, each about 400 nt size. The six UK strains (GU825959–GU825963, DQ434991, fragment sizes of 379–427 nt) overlapped with the 5′-end partial sequence of the fragment of the Chinese strain. For the five France strains, two of them (AY669847–AY669848, fragment sizes of 350 and 365 nt, respectively) overlapped the 5′-end partial sequence of the fragment of the Chinese strain, and the 606 nt size of three strains (AY230509–AY230511) overlapped the 3′-end partial sequence of the fragment of the Chinese strain. In summary, the nucleotide sequences of the BQCV China-JL1 strain coding the RdRp were longer sequences in theory, but the sequences associated with RdRp nucleotide sequences in GenBank was shorter, with only parts of the regions encoding RdRp.
A phylogenetic tree was built to analyse the 5′-end of the section (fig. 5) and showed that the Chinese BQCV, Hungary 10 and British genotypes (GU825962) formed as a subgroup. Within it, the genotypes further divided into two branches. One branch contained the Hungary 10 genotype, but the second branch divided into two further branches, Chinese and British genotypes.
Fig. 5. Phylogenetic tree analysis of the RdRp region of the Chinese BQCV-JL1 strain. The phylogenetic tree shows the BQCV isolates from China, Britain, South Korea and Central Europe. The nucleotide sequences were aligned using the ClustalW programme. The phylogenetic tree was then constructed using the Mega5 package and the maximum likelihood method with Kimura's 2-parameter model and a bootstrap value of 1000 replicates. The accession number of each viral sequence is indicated in parentheses.
ORF2 regions encoding the CP4
As shown in fig. 6, Poland 4 and Brazil strains (EU292211) alone or in a small group were closely related. Others included Chinese genotypes in another small group. The Chinese BQCV in this study showed the closest genetic relationship with another Chinese strain (KJ599571) and was also closely related to other Korean strains (EU375537, EU770973, JX149531, JN542438).
Fig. 6. Phylogenetic tree analysis of the capsid polyprotein region of the Chinese BQCV-JL1 strain. The phylogenetic tree shows the BQCV isolates from China, Britain, South Korea and Central Europe. The nucleotide sequences were aligned using the ClustalW programme. The phylogenetic tree was then constructed using the Mega5 package and the maximum likelihood method with Kimura's 2-parameter model and a bootstrap value of 1000 replicates. The accession number of each viral sequence is indicated in parentheses
Discussion
The molecular comparison of different isolates is an exact and reliable method because it directly detects changes in genetic information during the course of evolution, and the resulting data can be evaluated statistically. There only six complete genome sequences of BQCV in GenBank, including the South African reference strain (AF183905, 8550 nt) (Leat et al., Reference Leat, Ball, Govan and Davison2000); the Central European genotypes, notably Hungary 10 (EF517515, 8335 nt), Poland 4 (EF517519, 8333 nt), Poland 5 (EF517520, 8298 nt) and Poland 6 (EF517521, 8302 nt) (Tapaszti et al., Reference Tapaszti, Forgach, Kovago, Topolska, Nowotny, Rusvai and Bakonyi2009); and the South Korean strains (JX149531, 8392 nt) (Reddy et al., Reference Reddy, Noh, Choe, Kweon, Yoo, Doan, Ramya, Yoon, Nguyen, Nguyen, Quyen, Jung, Chang and Kang2013).
In 2001, Benjeddou and Leat designed a pair of oligo primers selected for reverse transcription polymerase chain reaction (RT–PCR), specific for the capsid region of the South African reference strain (AF183905). The synthetic sequence was located at 7850–8550 nt (701 nt). It is well known that the complete genome sequence of the South African reference strain was 8550 nt in length. The last nucleotide of the fragment was just the last of the whole genome sequence, and so were the other five complete genome sequences published in GenBank. The present study sequenced the whole genome sequence of the Chinese BQCV strain and found it to be 8358 nt in length. The region encoding the capsid protein was located at 7658–8358 nt, similar to the South African reference strain. It was thus speculated that the strain attained was BQCV. After BLAST searching, we obtained 102 homologous sequences all belonging to BQCV strains, and since no whole genome sequence of BQCV Chinese strains has been reported, this study for the first time isolated a Chinese BQCV strain named BQCV-JL1. Using the description of the results in this study, it was speculated that ORF3 was part of ORF1 and part of the regions encoding RdRp. The ORF3 fragment was just at the 3′-end of the sequences encoding the RdRp protein (4891–5433 nt), further indicating that the fragment was part of ORF1. Three French genotypes (AY230509–AY230511) reported in GenBank, which were sequences that encoded RdRp, shared 90–91% of identity with the ORF3 sequence in this study. The three French genotype sizes were 606 nt, and the ORF3 size was 543 nt. Therefore, it was part of ORF1 and even part of the regions that encoded RdRp. According to the results of the comparison of ORF3 using BLAST, it was found that the initial position of the three French genotypes was 68–71 nt, which terminated at the full length sequence published. Hence, adding 63 nucleotides to the ORF3 made its length 606 nt. The key step was to locate the 63 nt to put it at the 5′- or 3′-end position. The two cases were respectively attempted. The results showed that adding the 5′-end to the Chinese strains changed the initial position of the French sequences. Adjusting the position further revealed ORF3 was located at position 4824–5427 nt (604) of the whole genomic sequence. The three French genotypes started with the first nt positon. The regions encoding the RdRp of the six complete genome sequences published in GenBank were just the penultimate nucleotides of their ORF1, respectively. It was thus deduced that the position of ORF1 of the Chinese strain was at 5428 nt, but the three encoded nucleotides (5426–5428) were ATG, not a termination codon. However, when one nucleotide was moved to the 3′-end, the three continuous nucleotides happened to be the termination codon TGA. Thus, the ORF1 of the Chinese strain was located between 546 and 5429 nt (4884 nt), only three nucleotides less than the ORF1 of European genotypes, six nucleotides less than the South Korean strain, and 84 nt less than the South African reference strains.
The nucleotides between 4674–4676 nt of the Chinese strain were TGA, a termination codon, thereby resulting in the first possible formation of the ORF1, which was shorter than the other six strains. So it was inferred that this position formed the substitutions and insertions/deletions, and then acted as a termination codon. As mentioned in the analysis above, ORF3 was actually part of ORF1, and further analysis comparing the nucleotide sequence of the BQCV China-JL1 strain to the other six strains revealed that the final position of ORF1 was 546–5429 nt.
Leat first confirmed the complete genome sequence of the South African BQCV strain in 2000. Before this, researchers isolated a strain from the original BQCV strains (AF125252, 815 nt, England, 1999). A virus with picorna-like biophysical properties was isolated from South African honey bees by Leat et al., On the basis of serology, it was identified as an isolate of BQCV. Deduced aa sequences for ORF1 and ORF2 were most similar to the non-structural and structural proteins, respectively, of the Drosophila C virus (DCV), Rhopalosiphum padivirus (RhPV), Himetobi P virus (HiPV) and Plautia stali intestine virus (PSIV). It was proposed that BQCV belonged to the group of picorna-like, insect-infecting RNA viruses constituted by DCV, RhPV, HiPV and PSIV. It has been proposed that positive-stranded RNA viruses should be classified based on a comparison of the aa sequences of their RdRp. However, their relationships with established members of the picorna-like virus lineage have only recently been examined based on comprehensive genome sequence information. Comparing the obtained BQCV South Africa strain (AF183905) to the above four viruses revealed the following commonalities: all of them were bicistronic viruses, with a 5′-end ORF encoding the replicase polyprotein, and a 3′-end ORF encoding the capsid protein. It has been reported that the translation of 3′-proximal ORF (ORF2) in PSIV was initiated from the IRES, and the comparison of the IRES sequences derived from these five strains revealed high similarity between BQVC and PSIV, indicating that BQCV may also have an IRES inducing the translation of ORF2. By searching the ORF of the South Africa strain (AF183905) using DNASTAR, we found that the ORF1 and ORF2, respectively, located between 658–5625 (4968 nt) and 5942–8395 (2454 nt); whereas Leat et al., reported the starting position of ORF2 was 5834. This may due to the fact that the 5637–5836 sequence might be the IRES of BQCV that initiates the translation of ORF2, thus the 5834–5836 bases (CCU) may act as the initiation codon of ORF2 translation, since in the whole genome sequence of PSIV, the CUU locating at the corresponding position also initiated the translation of PSIV ORF2. Therefore, we suggest that the translation ORF2 in South Africa strain was initiated from 5834 nt rather than 5942 nt. According to this, the BQCV China-JL1 strain had a 200 nt (5444–5644 nt) fragment, with high identity similarity (97%) to the South African BQCV genome. The delicate difference was at the position of 5482 nt (Chinese BQCV-JL1 strain was A, but was missing in the case of the South African BQCV genome), and the other four differences were at the beginning of that fragment (BQCV China-JL1 strain was GTAA, while South African BQCV genome was TCTC). Without the four nucleotides, it would be reasonable to suggest that this region (5448–5644 nt) of the BQCV China-JL1 strain also acted as an IRES to mediate the initiation of ORF2 translation. The CCU bases (5642–5644) of BQCV China-JL1 aligned with the CCU initiation codon of the South African BQCV and the CUU initiation codon of the PSIV, suggesting that the CCU bases (5642–5644) acts as an initiation codon for the translation of ORF2 in BQCV China-JL1. The ORFs of this BQCV China-JL1 strain were searched for using DNASTAR software, and the first potential position of ORF2 was between nucleotides at position 5750–8203 (2454 nt). Hence, the ORF2 was ultimately located between nt positions 5642 and 8203.
On the basis of the phylogenetic analysis at the regions encoding the helicase protein, the Chinese BQCV strain did not show any close genetic relationship with the others. At the regions encoding the RdRp, it showed a close genetic relationship with the British strain (GU825962), and the capsid protein showed a close genetic relationship with the other Chinese and Korean strains. A phylogenetic analysis of the complete genome sequences of BQCV showed that the Chinese BQCV strain had a high level of similarity (93%) to the Hungary 10 genotype. Britain was one of the Western European countries, Hungary was a Central European country, and South Korea and China were Asian countries. Overall, the Chinese BQCV strain investigated in this study showed the closest genetic relationship with the European and Asian countries and a high divergence from the South African reference genotype. The study material for the Caucasus bee diseases and Italian bee hybrids showed that the genetic origin of the BQCV China-JL1 strain was from Europe or Asia. In addition, based on the description of Reddy et al. (Reference Reddy, Noh, Choe, Kweon, Yoo, Doan, Ramya, Yoon, Nguyen, Nguyen, Quyen, Jung, Chang and Kang2013), they speculate that the South Korea strain (JX 149531) was originated from European countries, thus we further speculate that BQCV China-JL1 strain was also originated from Europe. The 8358 nt length of the BQCV China-JL1 strain contained two large ORFs: the 5′-end of the ORF1 encoding the putative replicase polyprotein, and the 3′-end of the ORF2 encoding capsid protein. ORF1 had its starting position at 546 nt and terminated at 5429 nt. ORF2 was positioned between 5642 and 8203 nt. These two ORFs were separated by an intergenic region with a length of 212 nt (5430–5641). The 5′-untranslated region (UTR) was located at 1–545 nt (545 nt), and the 3′-UTR was located between 8204 and 8358 nt (155 nt). The conserved areas within the 5′-UTR region were highly similar (94–99%). Between nt positions 683 and 692, and 1183 and 1251, the South African reference contained continuous sequences that were 10 and 69 nt long, respectively. This part of the sequence was missing from the investigated Chinese genotype. Interestingly, these continuous sequence sections are also missing from the Korean genotype and all European genotypes. The two ORFs in the investigated Chinese genotype were separated by a 212 nt-long intergenic region between positions 5430 and 5641 and were highly conserved.
The BQCV-JL1 strain was found to share high sequence identity (93%) with the Hungary 10 genotype at the whole-genome level and analysis of the nucleotide and amino acid sequences revealed that the BQCV-JL1 strain also shows close genetic relationships with the South Korean strain, suggesting that both the BQCV-JL1 and South Korean strains may have migrated from European countries. BQCV-JL1 strain was different from the other six strains in dividing the nucleotides positions of ORF, which was because of the gene mutation.
Materials and methods
Sample collection
Honeybee samples were collected from a private apiary, Changchun City, Jilin Province, China. The studied field did not involve endangered or protected species and was located specifically at longitude 125.23629 and latitude 43.90445. The honey bees used in this study were Italian (A. mellifera ligustica) and Caucasian (A. mellifera caucasica) hybrid honey bees, and darkening or blackening of the body or pupae can be observed with naked eyes. They were collected between 2010 and 2014 and stored at −80°C until analysis.
Extraction of viral RNA
Ten honeybees infected with BQCV samples were homogenized in 10 ml sterile phosphate-buffered saline using a mortar and pestle. The samples were centrifuged for 30 min at 14,000 rpm, and the supernatant was used for RNA extraction using a QIAamp Viral RNA Mini Kit (Qiagen, Germany), according to the manufacturer's instructions, and a matching instrument. Total RNA was eluted in 30 µl of elution buffer and stored at −80°C until use.
RT–PCR amplification of viral RNA
Purified RNA was used to synthesize cDNA using a Prime Script II 1st Strand cDNA Synthesis Kit (TaKaRa, Japan) according to the manufacturer's protocol. A pair of oligo nucleotide primers was designed based on a previous report of the partial sequence of the BQCV capsid region (Benjeddou et al., Reference Benjeddou, Leat, Allsopp and Davison2001). Ten sets of primers were designed by Generay Biotechnology Co., (Shanghai) based on the sequence of the South African (AF183905) strain (Leat et al., Reference Leat, Ball, Govan and Davison2000), the Hungary 10 (EF517515) strain (Tapaszti et al., Reference Tapaszti, Forgach, Kovago, Topolska, Nowotny, Rusvai and Bakonyi2009) and the South Korea (JX149531) strain (Reddy et al., Reference Reddy, Noh, Choe, Kweon, Yoo, Doan, Ramya, Yoon, Nguyen, Nguyen, Quyen, Jung, Chang and Kang2013) to amplify overlapping PCR products that comprised the complete genome sequence of selected BQCV genotypes. The sequences, orientations, and locations of the primers, as well as the expected product sizes, are shown in tables 1 and 2. PCR was conducted using the following thermal cycling conditions: one cycle of initial denaturation at 95°C for 5 min, followed by 30 cycles of denaturation at 95°C for 1 min, annealing at 55°C for 1 min and extension at 72°C for 1 min. Finally, samples were maintained at 72°C for 10 min and kept at 4°C until electrophoresis was performed. Negative controls (H2O) were included in each PCR amplification. The amplified products were electrophoresed in agarose gel containing ethidium bromide, visualized under UV transillumination and photographed using a Kodak Digital Science Electrophoresis Documentation and Analysis System. The product sizes were determined with reference to a 2000 bp molecular weight ladder (Enzynomics, Korea).
Table 1. The oligo primer pair was selected for reverse transcription polymerase chain reaction (RT–PCR), specific for the capsid region of the selected BQCV genotypes.
Table 2. Oligo primer pairs were designed against overlapping sequences for the complete sequencing of the selected BQCV genotypes.
Nucleotide sequencing and computational analysis
Specific amplification products were excised from the agarose gel using a QIAquick Gel Extraction Kit (Qiagen, Germany) according to the manufacturer's instructions. Purified products were sequenced by Generay Biotechnology Co. (Shanghai). The generated consensus sequences were then identified and compared with sequences in GenBank using BLAST (Basic Local Alignment Search Tool) on the National Centre for Biotechnology Information (NCBI). The nucleotide sequence of each fragment was compiled and aligned to build a continuous complete sequence using DNASTAR and CodonCode Aligner software.
Phylogenetic tree construction and analysis of BQCV
A phylogenetic tree was constructed for the complete genome of the Chinese BQCV genotype along with the previously reported South African, South Korean, Poland 4–6 and Hungary 10 BQCV genotypes (AF183905, EF517515, EF517519, EF517520, EF517521 and JX149531). Multiple sequence alignment was performed using the ClustalW algorithm. The phylogenetic trees were constructed with the Mega5 package (Tamura et al., Reference Tamura, Peterson, Peterson, Stecher, Nei and Kumar2011) using the maximum likelihood method with a bootstrap value of 1000 replicates.
Nucleotide sequence accession numbers
The BQCV sequence that encoded the capsid protein described in this study was submitted to GenBank under accession number KM206738, and the accession number of BQCV China-JL1 strain is KP119603.
Acknowledgements
This study was supported by the National Natural Science Foundation of China (31001065) and the AQSIQ Technology Program (2013IK033).
