Occurrence and Molecular Phylogeny of Economically Relevant Viruses Infecting Honey Bees (Apis mellifera L.) of Bingöl Province, Turkey

From March to April 2020, a total of 128 apiaries were visited in eight different areas of the Bingöl province: 44 apiaries in Center (1), 10 in Karlıova (2), 19 in Solhan (3), 17 in Genç (4), 11 in Yedisu (5), 12 in Adaklı (6), 4 in Yayladere (7) and 11 in Kiğı (8) (Fig. 1). During the surveys, a total of 384 worker bees from three different hives of each apiary were randomly sampled from outside the beehives and stored at −20°C for subsequent use.

Fig. 1

To determine whether honey bees are infected with a certain virus, total RNA isolation and cDNA synthesis were performed as reported by Foissac et al. (2001) and Rüstemoğlu & Sipahioğlu (2019), respectively. All molecular processes were performed with the same device, an Eppendorf Mastercycler Personal (Hamburg, Germany). The synthesized cDNAs were used as a template for the PCR assays. Standard PCR reaction components in a 50 μl volume consisted of 2 μl cDNA, 36.6 μl nuclease-free water, 5 μl 10X PCR Buffer, 3 μl 25 mM MgCl₂, 1 μl 10 mM dNTP, 1 μl 20 μM virus-specific primer sets (Tab. 1) and 0.4 μl 5U/μl Taq DNA polymerase enzyme (Fermentas, Thermo Fisher Scientific, USA). The PCR thermocycling procedures were initial denaturation at 94°C for 2 min, followed by denaturation, 94°C/2 min, annealing: 55°C/1min (56°C/30s for SBV), extension: 68°C/2 min (30s for SBV) for 40 cycles and a final extension at 72°C for 5 min. In PCR tests, isolates confirming viral presence in preliminary studies were used as a positive control. The cDNA extracted from healthy bees was used as a negative control.

Upon completion of the amplification process, PCR products (15 μl) were run on a 1% agarose gel with 1.5 kb DNA marker (Fermentas, USA) and subsequently analyzed in a gel documentation device using UV light (Syngene™ UV Transilluminator 2020LM).

The virus-positive bands were gel-purified (Thermofisher Gel Extraction kit, USA) and ligated directly to the pGEM T-Easy bacterial plasmid using T4 DNA ligase enzyme (Promega, USA). The resulting recombinant plasmids (5 μl) were transformed into competent cell E. coli bacteria JM109 (100 μl) (Promega, USA) for blue-white colony selection. Recombinant plasmids containing viral inserts of correct sizes were purified from transformed bacteria growing in liquid LB medium with ampicillin derived from white colonies (GeneJET Plasmid Miniprep Kit, Thermo, USA). Approximately 25 μl of the recombinant plasmid was sequenced with the use of Next-Generation Sequencing (NGS) (Sentebiolab/Ankara/Turkey). DNA sequence data of honeybee viral genomes were submitted to GenBank.

To determine the phylogenetic relationships of SBV and BQCV isolates in Turkey's Bingöl province, the following identified sequences were first used: 700 nt and 342 nt sequences encoding a polyprotein gene of BQCV and SBV, respectively. Sequences were phylogenetically compared with other honeybee viral DNA sequences isolated worldwide with the use of the neighbor-joining algorithm by 1000-replicate bootstrap search. Multiple sequence alignments and molecular phylogeny inference were created by CLC Main Workbench 6.7.1 bioinformatics software (Qiagen). The phylogenetic dendrogram was constructed with the use of isolates identified in the current study, together with isolates from different hosts in different regions of the world and an outgroup.

We demonstrated the previously unevaluated infection rate and distribution pattern of BQCV, ABPV and SBV in eight different regions where beekeeping activities are intense in the province of Bingöl. In this study, we visited 128 apiaries for molecular determination of the presence of three viruses. We detected the presence of BQCV and SBV, respectively, but all apiaries were ABPV negative. SBV was found in nine apiaries with a prevalence of 7%, while BQCV was detected in fifteen apiaries with a prevalence of 11.7%. The overall infection rate was 18.75%. In molecular tests, double infection was also detected in an apiary with an infection rate of 0.78%. Detailed information is given in Tab. 2.

Amplified PCR amplicons of one isolate for BQCV and two isolates for SBV were successfully cloned and sequenced. The resulting BQCV sequence (MZ357974) was named Bingol 92, and the two SBV sequences were named Bingol 18 and Bingol 8 (MZ357972 and MZ357971) in the GenBank/EMBL database. Multiple alignments showed that the BQCV-Bingöl isolate had approximately 75.71%–96.58% nucleotide identity with the other BQCV world isolates retrieved from GenBank. It displayed a high sequence homology (96.58%) with the Brazil isolate (EU292211) from A. mellifera and the France isolate (MN565034) from Vespa velutina, and the least with the Polish (82.86%) (EF517519) and Swedish (75.71%) (KY363519) isolates from A. mellifera. We constructed a phylogenetic dendrogram based on multiple nucleotide sequence comparisons of approximately twenty isolates from different countries, including the Bingöl-BQCV isolate (MZ357974), with a bootstrap value of 1000 replicates using the neighbor-joining method. It was independently clustered with isolates from Brazil (EU292211), France (MN565034), Italy (MT416539) and Australia (KY465680) on the phylogenetic tree (Fig. 2). BLASTn analysis of the 342 bp sequences revealed that the Bingöl-SBV sequences (MZ357972 and MZ357971) had 100% nucleotide homology with each other, even though each was from apiaries in different locations, 92.98% with the other Turkish isolate (MH251274) from A. mellifera, and the least with the China isolate (85.96%) (JX679488). The neighbor-joining molecular phylogeny of Bingöl-SBV isolates (MZ357971 and MZ357972) was analyzed with the use of twenty isolates from different geographical regions, consistent with the results of the BLASTn analysis. They clustered with the other isolate from İzmir, Turkey (MH251274) due to the close similarity of the nucleotide sequence (Fig. 3).

Fig. 2

Fig. 3

The prevalence of BQCV infection of honey bees has been reported in different hosts worldwide. In previous studies, primer sets amplifying the specific viral gene regions, i.e. helicase, RNA-de-pendent RNA polymerase, protease and capsid protein gene, were used for accurate viral diagnosis (Zhang et al., 2012; Choi et al., 2015; Rüstemoğlu, 2020). Based on the PCR amplified 700 bp DNA sequence of the viral polyprotein gene, 11.7% and 7.03% of apiaries in this study were positive for BQCV and SBV as single infections, respectively. In addition, the BQCV prevalence found in the current study is lower than in such countries as France (86%), Austria (30%), and Syria (29%) (Tentcheva et al., 2004; Berényi et al., 2006; Abou Kubaa et al., 2018) but higher than in Spain (10.4%), England (1.4%), and Denmark (1%) (Higes et al., 2007; Baker & Schroeder 2008; Nielsen et al., 2008; Antunez et al., 2012). In Turkey, BQCV positivity was found to be considerably higher than that reported in our study, similar to those conducted worldwide. The presence of BQCV was first reported by Gumusova et al. (2010) at an infection rate of 21.42% in the Black Sea region (Samsun, Amasya, Rize). Oğuz et al. (2018) determined that the BQCV infection rate was 88.5% in twenty three of twenty six colonies in the Van province using the same pairs of primers. Kalaycı et al. (2020) tested seventy-six apiaries for BQCV in seven different provinces (Adana, Aydın, Bursa, İzmir, Kütahya, Muğla, and Manisa). Primers designed based on the partial coat protein gene produced 536 bp DNA amplicon in PCR tests, with a 28.9% (22/76) infection rate. In a comprehensive study conducted in Turkey's nine ecogeographic regions (Antalya, Bolu, Hatay, İstanbul, Kahramanmaraş, Mersin, Muğla, Sivas, and Yalova), Muz & Muz (2018) determined that the BQCV infection rate was 47.6%. SBV-parasitized honey bees were reported in numerous localities across the globe with infection incidences (Welch et al., 2009). In Turkey, SBV was found in nine of 128 surveyed apiaries with a 7.03% infection rate is lower than that in Greece (87%), and Denmark (81.25%) but higher than in Spain (1.1%) (Nielsen et al., 2008; Bacandritsos et al., 2010; Antunez et al., 2012). Furthermore, the prevalence of SBV in this study is lower than in studies that surveyed other regions of Turkey-12.2% in Hakkari and 22.3% in Burdur (Rüstemoğlu & Sipahioğlu, 2019; Kalaycı et al., 2020). Their high prevalence is probably because the sampling was conducted when brood production was increasing (de Miranda et al., 2011). However, sampling from within the hive can affect the frequency of certain viruses, especially SBV. It is noteworthy that BQCV and SBV rates, 11.7% and 7.03%, respectively, are lower than those detected in surveys conducted in different regions of Turkey, which indicates that the Bingöl province is a region where beekeeping can be good, both in terms of other factors affecting colony survival (low pollution, plant diversity, etc.) and reduction in the incidence of certain viruses.

Honey bees are commonly exposed to mixed infection by two or more viruses. In our tests, only one apiary (0.78%) was co-infected with both viruses. Simultaneous viral infections in apiaries have been reported in Turkey with different infection frequencies (Rüstemoğlu & Sipahioğlu, 2019; Kalaycı et al., 2020). This finding is also consistent with the results of surveys conducted by most researchers worldwide (Haddad et al., 2008; Nielsen et al., 2008).

ABPV is recognized as a prevalent viral agent in seemingly healthy apiaries, often responsible for sudden colony collapse, which beekeepers around the world find difficult to comprehend. Distribution of the infective agent was reported in various countries from almost every continent, with high frequencies of infection (Bakonyi et al., 2002; Bacandritsos et al. 2010). In the current study, we were unable to identify ABPV infection in any of the surveyed apiaries in Bingöl province, but screening of Turkish apiaries revealed that ABPV was present at an infection rate of 35.5% in seven regions from the Aegean region (Kalayci et al., 2020), 86.66% in Burdur province (Cagirgan et al., 2020), and 2.2% in Hakkari province (Rüstemoğlu & Sipahioğlu, 2016). However, ABPV was not detected in a previous study in the Black Sea Region (Gumusova et al. 2010), in line with the current study in the Eastern Anatolia Region, which hosts the province of Bingöl. Consequently, the studies conducted in the provinces of Hakkari and Bingöl in Turkey's Eastern Anatolia region indicated low ABPV incidence or no infection. This could be attributed to the fact that this region is not on the migratory beekeeping route or does not provide appropriate environmental conditions for honey bee disease, unlike the Aegean region, which has rich flora and suitable climate. Another possible reason is that mountainous regions and continental climate conditions restrict the flight of honey bees and thus prevent viral transmission among colonies.

National, documentation on viral disease of honey bees is adequate, but the phylogenetic relationships of Turkish isolates have not been adequately studied. Multiple alignments of the BQCV-Bingöl isolate revealed high sequence homology with the other BQCV world isolates recovered from the GenBank. Our isolate is independently clustered with isolates from Brazil, France, Italy and Australia on the phylogenetic dendrogram, which indicates that the polyprotein genome of our viral isolate is closely related to other isolates in evolutionary origin. This is probably due to the presence of substantially conserved regions in the relevant gene sequence, which is mainly due to the supplying of queen bees through cross-country trade.

Molecular characterization of Turkish-BQCV isolates was first revealed by Muz & Muz (2018), who had targeted the helicase and coat protein gene region. In this study, they determined nucleotide similarity between the helicase genes of BQCV isolates ranging from 37% to 85%. On the contrary, the coat protein nucleotide similarity of related isolates varied between 86% and 99% compared to other world isolates. These data indicate that the coat protein gene is more stable than the helicase. The extended nucleotide similarity (75.71–96.58%) in this study may be related to the higher genetic variation susceptibility and mutation probability occurring in honeybee viruses. In addition, the presence of recombination in viruses that are members of the Picornaviridae family is a common phenomenon (Tapaszti et al., 2009; Noh et al., 2013).

Cagirgan et al. (2020) reported that the Turkish-Burdur BQCV isolate was evolutionarily included in the same cluster with other Turkish isolates (BQCV-Aydın, -Muğla, and -İzmir), which could be attributed to the high sequence homology of partial polyprotein gene sequence (383 bp). Cagirgan & Yazici (2021) phylogenetically analyzed the 514 bp amplified gene sequence of six Turkish-BQCV isolates expressing structural polyprotein. Five viral sequences showed high sequence similarity to one another, while the other showed close association with the Polish isolates (EF517516, EF517517, and EF517518). Contrary to this report, our isolate (MZ357974) showed a distant phylogenetic relationship with the Polish isolate (EF517519) with 82.86% nucleotide sequence similarity. This inconsistency is probably due to the different sizes of the nucleotide sequences analyzed, 700 bp in our study and 514 bp in another study.

Although samples were collected from different apiaries, the Bingöl-SBV sequences in this study were completely identical to one another, and showed high sequence consensus with the other Turkish isolate from A. mellifera. Rüstemoğlu & Sipahioğlu (2019) recorded the nucleotide sequence encoding a partial coat protein gene of the Turkish-Hakkari SBV isolate and reported that this isolate showed 77–90% nucleotide identity. Cagirgan & Yazici (2021) defined the partial polyprotein gene sequence of Turkish SBV isolates (four-strain). When compared with the intercontinental records in the NCBI database, Turkish strains exhibited 89–90%, 88–90% and 87–89% nucleotide sequence similarity to European-Uruguayan, Asian and Korean strains, respectively. Kalaycı et al. (2019) determined that four Turkish genotypes showed a genetic variance of 10.4%, 11.4%, and 12.5% with European–South American, Asian and Korean genotypes, respectively. The nucleotide percentages obtained here are higher than the nucleotide sequence similarity results of other Turkish isolates in studies conducted in different parts of Turkey. This could be attributed to the recombination induced by the different selective pressure between different ecozones. On the contrary, the most specific feature of RNA viruses is their superior ability to generate and maintain extraordinary genetic variability, which underlies their survival and adaptation. The variation is optimized to be substantially below the deleterious mutation threshold and is produced through error-prone replication and extensive recombination between respective variants. Furthermore, viruses can incorporate large segments of genetic material into their genomes from other viruses, the host or any conceivable source. This situation results in the emergence of a dynamic and amorphous group of variants that are constantly changing and interconnected by genetic and functional relationships yet remain globally stable (de Miranda et al., 2011). Based on the neighbor-joining phylogenetic tree, Bingöl-SBV isolates (MZ357971 and MZ357972) were clustered with the other isolate from Turkey's İzmir (MH251274) due to the close similarity of the nucleotide sequence (Fig. 3).

Honey bee viruses have multiple transmission routes due to their strong adaptability. Viral RNA is vectored by many insects, including ants (Messor concolor), bumblebees (Bombus atratus), wax moth (Galleria mellonella), mites (V. destructor), honey bee lice (Braula schmitzi), darkling beetle (Alphitobius diaperinus), small hive beetle (Aethina tumida) and wasps (Vespa spp.) inside and outside of the colony (Eyer et al., 2009; Rosenkranz et al., 2010; Reynaldi et al., 2013; Li et al., 2016; Traiyasut et al., 2016; Avalos et al., 2019; Highfield et al., 2020; Rüstemoğlu, 2020). In addition, the epidemiology of BQCV in honey bees is directly related to microsporidian pathogens (N. apis and N. ceranae) (Oğuz et al., 2018). In the Bingöl province, the presence of Nosema spp. and Varroa was reported by Kutlu & Gül (2020), and the presence of ants, Vespa spp. and small hive beetle was observed in our surveyed areas.

A drug-based strategy has not yet been developed for the control of honey bee viruses, But the RNA interference (RNAi) approach is one of the most promising developments (Glover et al., 2011). However, it is unclear whether this costly achievement will be accepted by consumers in terms of its practicability in field beekeeping and any possible nucleic acid contamination in honey. Therefore, cultural control approaches and such beekeeping practices as hygienic working, pest control, timely replacement of the queen bee, and avoiding stress in honey bees interrupt viral transmission in terms of sustainable beekeeping.

To our knowledge, there has been no screening for the presence of honey bee viruses in eastern Turkey's Bingöl province to date. Therefore, our study is the first to provide evidence regarding BQCV and SBV viral diseases infecting honey bees. The spread of viral disease in honey bees is increasing gradually and is a cause for concern. Therefore, virus screening programs should be incorporated into the usual international trade control procedures and the health status of honey bees should be continuously monitored for beekeepers to prevent interregional virus spread.