# Virus Prevalence in Egg Samples Collected from Naturally Selected and Traditionally Managed Honey Bee Colonies across Europe

David Claeys Bouuaert, Lina de Smet, Marleen Brunain, Bjørn Dahle, Tjeerd Blacquière, Anne Dalmon, Daniel Dezmirean, Dylan Elen, Janja Filipi, Alexandru Giurgiu, et al.

<figure><img src="image_1.png" /></figure>

## To cite this version:

David Claeys Bouuaert, Lina de Smet, Marleen Brunain, Bjørn Dahle, Tjeerd Blacquière, et al. Virus Prevalence in Egg Samples Collected from Naturally Selected and Traditionally Managed Honey Bee Colonies across Europe. Viruses, 2022, 14 (11), pp.2442. 10.3390/v14112442. [https://doi.org/10.3390/v14112442](https://doi.org/10.3390/v14112442)

HAL Id: hal-04099644

https://hal.inrae.fr/hal-04099644

Submitted on 17 May 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

Article

# Virus Prevalence in Egg Samples Collected from Naturally Selected and Traditionally Managed Honey Bee Colonies across Europe

David Claeys Bouuaert <fnref n="1" />,* , Lina De Smet <fnref n="1" /> , Marleen Brunain <fnref n="1" />, Bjørn Dahle <fnref n="2" />,<fnref n="3" />, Tjeerd Blacquière <fnref n="4" />, Anne Dalmon <fnref n="5" /> , Daniel Dezmirean <fnref n="6" />, Dylan Elen <fnref n="7" />,<fnref n="8" /> , Janja Filipi <fnref n="9" />, Alexandru Giurgiu <fnref n="6" />, Aleš Gregorc <fnref n="10" />, John Kefuss <fnref n="11" /> , Barbara Locke <fnref n="12" />, Joachim R. de Miranda <fnref n="12" />, Melissa Oddie <fnref n="3" />, Delphine Panziera <fnref n="4" />, Melanie Parejo <fnref n="13" />, Maria Alice Pinto <fnref n="14" />,<fnref n="15" /> and Dirk C. de Graaf <fnref n="1" />

<figure><img src="image_4.png" /></figure>

check for
updates

**Citation:** Claeys Bouuaert, D.; De Smet, L.; Brunain, M.; Dahle, B.; Blacquière, T.; Dalmon, A.; Dezmirean, D.; Elen, D.; Filipi, J.; Giurgiu, A.; et al. Virus Prevalence in Egg Samples Collected from Naturally Selected and Traditionally Managed Honey Bee Colonies across Europe. *Viruses* **2022**, *14*, 2442. https://doi.org/10.3390/v14112442

Academic Editor: Huoqing Zheng

Received: 12 October 2022

Accepted: 2 November 2022

Published: 3 November 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_5.png" /></figure>

**Copyright:** © 2022 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article
distributed under the terms and
conditions of the Creative Commons
Attribution (CC BY) license (https://
creativecommons.org/licenses/by/
4.0/).

**Abstract:** Monitoring virus infections can be an important selection tool in honey bee breeding. A recent study pointed towards an association between the virus-free status of eggs and an increased virus resistance to deformed wing virus (DWV) at the colony level. In this study, eggs from both naturally surviving and traditionally managed colonies from across Europe were screened for the prevalence of different viruses. Screenings were performed using the phenotyping protocol of the 'suppressed in ovo virus infection' trait but with qPCR instead of end-point PCR and a primer set that covers all DWV genotypes. Of the 213 screened samples, 109 were infected with DWV, 54 were infected with black queen cell virus (BQCV), 3 were infected with the sacbrood virus, and 2 were infected with the acute bee paralyses virus. It was demonstrated that incidences of the vertical transmission of DWV were more frequent in naturally surviving than in traditionally managed colonies, although the virus loads in the eggs remained the same. When comparing virus infections with queen age, older queens showed significantly lower infection loads of DWV in both traditionally managed and naturally surviving colonies, as well as reduced DWV infection frequencies in traditionally managed colonies. We determined that the detection frequencies of DWV and BQCV in honey bee eggs were lower in samples obtained in the spring than in those collected in the summer, indicating that vertical transmission may be lower in spring. Together, these patterns in vertical transmission show that honey bee queens have the potential to reduce the degree of vertical transmission over time.

**Keywords:** honey bee; suppressed in ovo virus infection; vertical transmission; virus resistance

<footnote n="1">¹ Department of Biochemistry and Microbiology, Ghent University, 9000 Ghent, Belgium</footnote>

<footnote n="2">² Faculty of Environmental Sciences and Natural Resource Management, Norwegian University of Life Sciences, 1430 Ås, Norway</footnote>

<footnote n="3">³ Norwegian Beekeepers Association, 2040 Kløfta, Norway</footnote>

<footnote n="4">⁴ Bees@wur, Wageningen University & Research, 6708 Wageningen, The Netherlands</footnote>

<footnote n="5">⁵ Abeilles et Environnement, INRAE, 84914 Avignon, France</footnote>

<footnote n="6">⁶ Department of Apiculture and Sericulture, University of Agricultural Sciences and Veterinary Medicine, 400372 Cluj-Napoca, Romania</footnote>

<footnote n="7">⁷ Department of Molecular Ecology & Evolution, School of Natural Sciences, Bangor University, Bangor LL57 2DG, UK</footnote>

<footnote n="8">⁸ Taskforce Research, ZwarteBij.org VZW, 9890 Gavere, Belgium</footnote>

<footnote n="9">⁹ Department of Ecology, Agronomy and Aquaculture, University of Zadar, 23000 Zadar, Croatia</footnote>

<footnote n="10">¹⁰ Faculty of Agriculture and Life Sciences, University of Maribor, 2311 Pivola, Slovenia</footnote>

<footnote n="11">¹¹ Le Rucher D'Oc, 31200 Toulouse, France</footnote>

<footnote n="12">¹² Department of Ecology, Swedish University of Agricultural Sciences, 75007 Uppsala, Sweden</footnote>

<footnote n="13">¹³ Applied Genomics and Bioinformatics, University of the Basque Country, 48940 Leioa, Spain</footnote>

<footnote n="14">¹⁴ Centro de Investigação de Montanha, Instituto Politécnico de Bragança, Campus de Santa Apolónia, 5300 Bragança, Portugal</footnote>

<footnote n="15">¹⁵ Laboratório Associado Para a Sustentabilidade e Tecnologia em Regiões de Montanha (SusTEC), Instituto Politécnico de Bragança, Campus de Santa Apolónia, 5300 Bragança, Portugal</footnote>

<footnote>* Correspondence: david.claeysbouuaert@ugent.be; Tel.: +32-(0)-9-264-51-51</footnote>
# 1. Introduction

Disease pressure is an inherent driver of the evolution of eusociality [1,2] or social task division [3], and it forms an important component in the evolution of western honey bees (*Apis mellifera*). With the arrival of the Varroa mite (*Varroa destructor*), the virus landscape in honey bee colonies was considerably changed by the introduction of a new transmission pathway, thereby influencing virus virulence and evolution [4–13]. Through Varroa-mediated transmission, virus diseases have become one of the most important proximate causes of colony mortality and honey bee decline [14–21].

Of the 72 virus species that have been identified in honey bees [22], the most commonly occurring belong to the families Iflaviridae and Dicistroviridae [23], particularly the sacbrood virus (SBV), black queen cell virus (BQCV), acute bee paralysis virus (ABPV) and deformed wing virus (DWV). Both ABPV and DWV consist of a complex of closely related, co-circulating master variants capable of forming viable recombinants [22,24,25]. The DWV complex is best described as a group of functionally and genetically compatible minor and major variants and their recombinants based on four master strains [26], of which DWV-A and DWV-B are currently the most common [27–32]. Dynamics in the presence and abundance of honey bee viruses show strong seasonal and geographical variation [33–35]. This variation is driven by the local adaptations of the virus, host and vector species, as well as by the specific characteristics of each virus [36–40]. Together, they form a geographic mosaic of coevolution [41].

In the first years after managed colonies are left untreated against the Varroa mite, colony mortality increases considerably [42,43]. This results in strong selective pressure forcing bees, mites and the viruses to adapt to each other. Most naturally surviving populations consist of unmanaged or feral colonies [40]. In managed colonies, two approaches have been described to transition from treated colonies to naturally surviving colonies. The first consists of leaving a large number of colonies unmanaged with respect to swarming, re-queening and Varroa control [39] and allowing natural selection to take place. This is described as the ‘Bond’ test: ‘live and let die’ [39]. A second approach, named ‘Darwinian black box’, builds further on this by adding selection for strong spring development [43]. One of the best studied naturally surviving populations with regard to virus-host co-evolution is an isolated, closed honey bee population located at the tip of the Näsudden peninsula in the south of Gotland, a Swedish island in the Baltic sea. After implementing the Bond test in 1999, these honey bee colonies evolved an increased tolerance for DWV infections [38,44]. In addition, BQCV and SBV infections were less abundant in the autumn and early spring, possibly due to the reduced colony size of the Gotland colonies in these seasons [45].

Each response to a parasite influences transmission dynamics within and between honey bee colonies. As honey bees live in large groups, the transmission of viruses through trophallaxis, feeding or body contact occurs frequently. This form of transmission between individuals of the same generation is defined as horizontal transmission. However, transmission between generations by either eggs or semen is defined as vertical transmission [27]. Virus infections of queens, or their eggs, have been shown to interfere with normal egg development, to elicit a stress response in eggs [46] and to cause important health risks for the queen herself [47–50]. The importance of the honey bee queen in the viral dynamics of the colony was recently highlighted with the discovery of the ‘suppressed in ovo virus infection’ trait (SOV) [49]. This trait is described by the virus status of a sample containing 10 pooled drone eggs, and it reflects the degree of vertical transmission of viruses at the time of sampling. Colonies headed by a queen laying virus-free eggs have been found to show fewer and less severe DWV infections in almost all developmental stages of both drones and workers [49]. In addition, this potential to suppress viral infections is heritable [49] and alters the tissue specificity of DWV [50]. Drone eggs are unfertilized and therefore only reflect the vertical virus transmission of the queen to her offspring [51]. As drone eggs are only produced in the spring and summer, the choice of drone eggs compared to worker eggs includes the limitation of when samples can be collected.
The aim of this study was to compare virus infections in eggs collected from naturally surviving colonies (NSCs) and traditionally managed colonies (TMCs) across Europe, along with an analysis of the effects that queen age and sampling season have on infection patterns in these eggs. The previously described SOV phenotyping protocol [49] was used to screen for the presence of viruses after implementing two improvements: First, the detection of viral pathogens using qPCR instead of *end-point* PCR allowed for the quantification of the viral load of the eggs and lowered the detection threshold. Second, as the SOV trait is associated with increased virus resistance across DWV genotypes [50], a shift was made from screening for DWV-A only to a generic screening for the DWV complex. Overall, this research improved our understanding of how the patterns of the vertical transmission of viruses differ across Europe and in different breeding program selection strategies.

**2. Materials and Methods**

**2.1. Virus Screening in Eggs**

The 187 samples collected as part of the Flemish bee-breeding program in 2020 were used to compare the transition from an *end-point* PCR to a qPCR approach and to compare the quantification of DWV infections using primers specific to DWV-A or DWV-B and a generic DWV primer (DWV-Fam). All samples were collected following the phenotyping protocol for the SOV trait, as described by de Graaf et al. (2020), and they were screened using qPCR for DWV-A, DWV-B, DWV-Fam, SBV, ABPV and BQCV. To compare the performance of the *end-point* PCR to that of the qPCR, positive samples covering a 10¹–10⁸/10 egg range were selected and analyzed using *end-point* PCR for SBV, ABPV, DWV-A and BQCV.

**2.2. Egg Sample Collection across Europe**

Egg samples were collected in 9 countries across Europe during either the spring (the beginning of March to the end of May) or summer (the beginning of June to the end of July) of 2020, depending on the presence of drone brood in each country and in each season (Table 1). Each country sampled between 6 and 12 TMCs (colonies managed following local standard practices, including treatment against the Varroa mite) and, if present, between 4 and 13 NSCs (colonies from populations that survive without treatment against the Varroa mite). In total, 53 samples were collected from NSCs, and 160 were collected from TMCs (including 72 from Slovenia). The samples from Slovenia were collected in the scope of a different project, hence the larger sample size. The colonies were managed following local standard practices, and the queens descended from locally adapted or native stock. The treatment of the TMCs against the Varroa mite was performed with registered products in each country. The management of the NSCs was conducted according to the 'Darwinian black box' selection method [43] or the 'Bond' test [39]. From each of the 213 sampled colonies, a pooled sample of 10 drone eggs was collected following the phenotyping protocol of the SOV trait, as described by de Graaf et al. (2020). If drone eggs were not present in the spring or summer and if attempts to induce drone laying did not succeed, worker eggs were collected instead (as was the case for 34 colonies). All samples were immediately stored at -20 °C and kept in a cold chain during transport to Belgium, where they were analyzed for DWV-Fam, BQCV, SBV and ABPV using RT-qPCR. For each sampled colony, information was gathered on the sampling season (spring or summer), subspecies, queen age, beekeeping method (for both TMCs and NSCs) and the presence of clinical signs at the time of sampling. This information was used to explain possible outliers and to look for correlations between multiple factors. Additional samples were collected if the apiary was composed of more than 10 colonies with drone eggs present during sampling. Supplementary Table S2 provides an overview of the sampled populations and countries, the number of worker egg samples collected in each country and the presence of clinical signs at the time of sampling, and it lists the location of the sampled populations, including the year of establishment of the NSC populations.
<table><caption>**Table 1.** Virus prevalence in pooled egg samples from each participating country for both naturally surviving colonies (NSCs) and traditionally managed colonies (TMCs).</caption><thead><tr><th rowspan="2">Country</th><th rowspan="2">Sampling Season</th><th rowspan="2">NSC/TMC</th><th rowspan="2">No. of Sampled Colonies</th><th rowspan="2">No. of Virus-Free Samples</th><th colspan="5">No. of Samples Positive For:</th><th colspan="3">Mean Infection Load (Log10/Egg)</th></tr><tr><th>DWV</th><th>BQCV</th><th>SBV</th><th>ABPV</th><th>DWV</th><th>BQCV</th><th>SBV</th><th>ABPV</th></tr></thead><tbody><tr><td rowspan="2">Belgium</td><td rowspan="2">Spring</td><td>NSC</td><td>10</td><td>1</td><td>9 (90%)</td><td>0</td><td>0</td><td>0</td><td>4.3</td><td></td><td></td><td></td></tr><tr><td>TMC</td><td>11</td><td>4</td><td>7 (64%)</td><td>2 (18%)</td><td>0</td><td>0</td><td>6.1</td><td>4.2</td><td></td><td></td></tr><tr><td rowspan="2">Croatia</td><td rowspan="2">Summer</td><td>TMC</td><td>10</td><td>3</td><td>4 (40%)</td><td>3 (30%)</td><td>1 (10%)</td><td>0</td><td>4.7</td><td>5.3</td><td>3.3</td><td></td></tr><tr><td>NSC</td><td>13</td><td>2</td><td>7 (54%)</td><td>11 (85%)</td><td>0</td><td>0</td><td>5.6</td><td>5.5</td><td></td><td></td></tr><tr><td rowspan="2">France</td><td rowspan="2">Spring</td><td>TMC</td><td>10</td><td>2</td><td>1 (10%)</td><td>8 (80%)</td><td>0</td><td>0</td><td>5.8</td><td>6.4</td><td></td><td></td></tr><tr><td>NSC</td><td>10</td><td>2</td><td>8 (80%)</td><td>3 (30%)</td><td>0</td><td>0</td><td>6.2</td><td>5.4</td><td></td><td></td></tr><tr><td rowspan="6">the Netherlands</td><td rowspan="6">Partly in spring and summer</td><td>TMC</td><td>6</td><td>1</td><td>5 (83%)</td><td>3 (50%)</td><td>0</td><td>0</td><td>5.4</td><td>4.9</td><td></td><td></td></tr><tr><td>NSC</td><td>10</td><td>1</td><td>9 (90%)</td><td>1 (10%)</td><td>0</td><td>0</td><td>5.1</td><td>4.9</td><td></td><td></td></tr><tr><td>TMC</td><td>10</td><td>5</td><td>4 (40%)</td><td>3 (30%)</td><td>0</td><td>0</td><td>4.3</td><td>4.7</td><td></td><td></td></tr><tr><td>TMC</td><td>10</td><td>1</td><td>8 (80%)</td><td>1 (10%)</td><td>0</td><td>0</td><td>5.1</td><td>3.3</td><td></td><td></td></tr><tr><td>NSC</td><td>4</td><td>1</td><td>0</td><td>3 (75%)</td><td>0</td><td>0</td><td>4.0</td><td>4.5</td><td></td><td></td></tr><tr><td>TMC</td><td>9</td><td>4</td><td>2 (22%)</td><td>4 (44%)</td><td>0</td><td>0</td><td>4.0</td><td>4.0</td><td></td><td></td></tr><tr><td rowspan="2">Norway</td><td rowspan="2">Summer</td><td>TMC</td><td>72</td><td>27</td><td>38 (53%)</td><td>11 (15%)</td><td>2 (2%)</td><td>2 (2%)</td><td>5.3</td><td>5.2</td><td>3.2</td><td>3.8</td></tr><tr><td>NSC</td><td>10</td><td>4</td><td>5 (50%)</td><td>1 (10%)</td><td>0</td><td>0</td><td>5.5</td><td>4.8</td><td></td><td></td></tr><tr><td rowspan="2">Portugal</td><td rowspan="2">Spring</td><td>TMC</td><td>6</td><td>5</td><td>1 (16%)</td><td>0</td><td>0</td><td>0</td><td>5.0</td><td></td><td></td><td></td></tr><tr><td>NSC</td><td>12</td><td>11</td><td>1 (8%)</td><td>0</td><td>0</td><td>0</td><td>4.3</td><td></td><td></td><td></td></tr><tr><td rowspan="2">Romania</td><td rowspan="2">Spring</td><td>TMC</td><td>10</td><td>4</td><td>5 (50%)</td><td>1 (10%)</td><td>0</td><td>0</td><td>5.5</td><td>4.8</td><td></td><td></td></tr><tr><td>NSC</td><td>6</td><td>5</td><td>1 (16%)</td><td>0</td><td>0</td><td>0</td><td>5.0</td><td></td><td></td><td></td></tr><tr><td rowspan="2">Slovenia</td><td rowspan="2">Spring</td><td>TMC</td><td>10</td><td>4</td><td>5 (50%)</td><td>1 (10%)</td><td>0</td><td>0</td><td>5.5</td><td>4.8</td><td></td><td></td></tr><tr><td>NSC</td><td>6</td><td>5</td><td>1 (16%)</td><td>0</td><td>0</td><td>0</td><td>5.0</td><td></td><td></td><td></td></tr><tr><td rowspan="2">Spain</td><td rowspan="2">Spring</td><td>TMC</td><td>10</td><td>4</td><td>5 (50%)</td><td>1 (10%)</td><td>0</td><td>0</td><td>5.5</td><td>4.8</td><td></td><td></td></tr><tr><td>NSC</td><td>6</td><td>5</td><td>1 (16%)</td><td>0</td><td>0</td><td>0</td><td>5.0</td><td></td><td></td><td></td></tr><tr><td rowspan="2">Sweden</td><td rowspan="2">Summer</td><td>TMC</td><td>12</td><td>11</td><td>1 (8%)</td><td>0</td><td>0</td><td>0</td><td>4.3</td><td></td><td></td><td></td></tr><tr><td>NSC</td><td>12</td><td>11</td><td>1 (8%)</td><td>0</td><td>0</td><td>0</td><td>4.3</td><td></td><td></td><td></td></tr></tbody></table>
## 2.3. RNA Extraction and cDNA Synthesis

All samples were first homogenized in the presence of zirconium beads in 0.5 mL QIAzol lysis reagent (Qiagen, Hilden, Germany). RNA was extracted using an RNeasy Lipid Tissue Mini Kit (Qiagen) according to the manufacturer's instructions, including a DNAse step, and it was finally eluted in 30 µL elution buffer. The concentration of the total RNA was measured using Nanodrop (Isogen, De Meern, The Netherlands). Using random hexamer primers, 200 ng of RNA was retro-transcribed with a RevertAid H Minus First Strand cDNA Synthesis Kit (Thermo Scientific, Waltham, MA, USA). Honey bee β-actin was used to control RNA integrity.

## 2.4. End-Point PCR

All *end-point* PCR reaction mixtures contained 2 µM of each primer (see Supplementary Table S1), 1 mM MgCl₂, 0.2 mM dNTPs each, 1.2 U HotStarTaq Plus DNA polymerase (Qiagen) and 2 µL cDNA product. The *end-point* PCR assays were performed using the following cycling conditions: 95 °C for 5 min; 94 °C for 30 s, 55 °C for 30 s, 72 °C for 1 min, 35 cycles; final elongation 72 °C for 10 min, hold 4 °C. The *end-point* PCR amplicons were analyzed by electrophoresis using 1.5% agarose gels stained with ethidium bromide and visualized under UV light. Positive and negative controls were included in each run.

## 2.5. qPCR

The virus load qPCR determination was performed using Platinum™ SYBR™ Green qPCR SuperMix-UDG (Thermo Scientific). Each reaction consisted of 0.4 µM of each primer (sequences provided in Supplementary Table S1), 11.45 µL RNase-free water, 12.5 µL SYBR Green and 1 µL of cDNA template. All samples were run in duplicate in a three-step RT-qPCR. The thermal cycling conditions started with an initial activation stage at 95 °C for 2 min, followed by 35 cycles of a denaturation stage at 95 °C for 15 s, an annealing stage at 58 °C for 20 s and an extension stage at 72 °C for 30 s. This procedure was followed by a melting curve analysis to confirm the specificity of the product (55–95 °C with increments of 0.5 °C s⁻¹). Each plate included a no-template control and a positive control. A standard curve obtained through an 8-fold 5× serial dilution of a known amount of viral plasmid loads (range of 10⁴–10¹⁰ copies/µL) was used for absolute quantification. All data were analyzed using CFX Manager™ 3.1 software (Bio-Rad). Baseline correction and threshold setting were performed using the automatic calculation offered by the same software. The maximum accepted quantification cycle (Ct) difference between replicates was set to two Ct. The successful amplification of the β-actin internal reference gene was used to confirm RNA integrity throughout the entire procedure [52]. For each sample, the virus load for the 200 ng RNA included in the cDNA reaction was multiplied to account for the total volume of RNA per sample, and it was subsequently divided by 10 to represent data as total viral load per individual egg. The linear standard equations for the plasmid standards and primers specific for each virus were as follows: $Ct = -3.519 \times +47.762, R^2 = 0.997$ for DWV; $Ct = -4.458 \times +50.946, R^2 = 0.865$ for ABPV; $Ct = -4.260 \times +55.151, R^2 = 0.943$ for BQCV; and $Ct = -4.571 \times +45.651, R^2 = 0.938$ for SBV.

## 2.6. Statistics

The viral loads for each sample were Log10-transformed to improve data visualization. Detection thresholds for all pathogens were set at 30 Ct (corresponding to 103 copies for DWV, BQCV, SBV and APBV). Below this threshold, samples cannot be reliably quantified using qPCR [53]. RStudio version 3.6.1 was used for data analyses and visualization. Analyses of the differences in the number of infections between groups were conducted using chi-squared tests. For comparisons between the infection loads, T-tests were used. All tests were checked for and complied with the required assumptions.
## 3. Results

### 3.1. SOV Phenotyping Method

The virus detection thresholds on the end-point PCR (based on the end-point PCR of the samples positive on the qPCR along a 10¹–10⁸ copy/reaction range of the starting template) were 10² for DWV-A, 10⁷ for BQCV, 10⁶ for ABPV and 10⁸ for SBV (see Supplementary Figure S1). A total of 187 samples were screened for DWV-A, DWV-B and DWV-Fam (generic DWV primer). Of these, 153 (82%) showed amplification when screening with DWV-Fam. One sample was only amplified with the DWV-A assay, and two samples were only amplified with the DWV-B assay. Of the 153 samples that were amplified with the DWV-Fam assay, 98 (64%) were also amplified with the DWV-B assay, 2 (1%) were amplified with the DWV-A assay, and 6 (4%) were amplified with both the DWV-A and DWV-B assays. The remaining 47 samples (25%) were only amplified with the DWV-Fam assay and not with either the DWV-A or DWV-B assay. The median infection load for DWV-Fam was 5.8 Log₁₀ virus copy number/egg and was on average 1.69 Log₁₀ higher than the sum of DWV-A and DWV-B.

### 3.2. Virus Prevalence

Table 1 shows an overview of the different virus prevalence across all samples, and it presents the number of collected samples per country, the number of infections for each virus, and the mean infection load for each country and for each selection strategy (TMCs or NSCs). In total, 53 samples were collected from NSCs, and 160 were collected from TMCs (including 72 from Slovenia). Of the 213 pooled egg samples screened, most infections were with DWV (51%), followed by BQCV (25%). Only three samples were infected with SBV, and two samples were infected with ABPV. Multiple virus infections in the same sample occurred in only 14% of the samples (28 samples were infected with two viruses, and one was infected with three viruses). No virus infections were found in 35% (74/213) of the samples. The differences in infection frequencies varied considerably between and within countries. The worker egg samples from Norway and Sweden had lower infection frequencies (11/34) than the drone egg samples from other locations (61/107). The Slovenian samples were collected in the scope of a different project, hence the larger number of samples. To avoid an uneven distribution of the sample size across groups, the Slovenian samples were not included for further analyses in this study. The differences between subspecies could not be analyzed due to the high variability between countries and the hybridization between subspecies. Due to bad weather conditions or the inaccessibility of some locations, some countries (The Netherlands, Romania and Sweden) were not able to sample the requested number of colonies.

Figure 1 shows the percentage of virus infections occurring in the two sampling seasons (spring and summer). The samples collected in the spring had significantly higher infection frequencies than the samples collected in the summer for both DWV (X² (1, N = 141) = 9.4, *p* < 0.05) and BQCV (X² (1, N = 141) = 12.3, *p* < 0.05) but not for SBV.

### 3.3. Natural Survivors vs. Traditionally Managed Colonies

Figure 2 shows the percentage of virus infections (A) and the infection loads (B) for both NSCs and TMCs and for each virus. The infection frequencies were significantly higher in NSCs than in TMCs for DWV (X² (1, N = 111) = 8.6, *p* < 0.05) but not for BQCV (X² (1, N = 111) = 0.1, *p* = 0.75). No significant differences were found between the infection loads of TMCs and NSCs for DWV (*t* (68) = -0.6, *p* = 0.52; TMCs = 5.2 ± 0.4, NSCs = 4.9 ± 1.6) and BQCV (*t* (33) = 1.6, *p* = 0.11; TMCs = 4.4 ± 0.9, NSCs = 4.8 ± 0.1). On average, most infections hovered around 10⁵ for both DWV and BQCV in this sample cohort.
<figure><img src="image_6.png" /><figcaption>**Figure 1.** Percentage of virus infections in each sampling season. Significant differences are indicated with *. The numbers on top of each bar represent the number of positive samples and the total sample for each group.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 2.** Percentage of viral infections (A) and infection loads (B) for naturally selected and traditionally managed colonies. Data are provided for each virus. Significant differences are indicated with *. The numbers on top of each bar represent the number of positive samples and the total sample for each group.</figcaption></figure>

### 3.4. *Queen Age*

Figure 3A shows the frequency of infection with DWV or BQCV for each queen age and for both TMC and NSC groups. For DWV, a significant decrease in the percentage of infected samples was found in TMCs between queens aged 0 and 1 year ($X^2$ (1, N = 59) = 3.9, *p* < 0.05). This trend continued with a lower infection frequency in queens aged 2 years (1/10) than in queens aged 1 year (14/42), albeit not significant. No differences in infection frequencies were found between queen ages in NSCs. Figure 3B shows the infection loads of DWV and BQCV for each queen age for both TMC and NSC groups. As previously shown, the infection load did not differ between the two groups. Comparing infection loads between queen ages showed significantly higher infection loads in queens aged 0 years (M = 5.7, SD = 0.2) than in queens aged 1 year for DWV (M = 4.8, SD = 1.6; t (29) = 2.6, *p* < 0.05) for both TMCs and NSCs. Albeit not significant, the mean DWV infection load for queens aged 2 years (M = 5.3) was lower than the mean of queens aged 1 year (M = 6.2). No significant differences in infection frequencies or infection load were found for BQCV. The infection load did however show a similar general trend, with mean infection loads decreasing with age (5.0 in queens aged 0 years, 4.6 in queens aged 1 year and 4.4 in queens
aged 2 years). Interestingly, the spread in the infection loads of queens aged 0 years largely lacks infection loads lower than $10^7$ for DWV and BQCV. Seasonal differences in sample collection did not influence the infection frequency of DWV in queens aged 0 years ($X^2$ (1, N = 25) = 1.69, *p* = 0.16) and older queens ($X^2$ (1, N = 97) = 2.2, *p*-value = 0.14). Queen age was unknown for 19 out of 141 samples (13%).

<figure><img src="image_8.png" /><figcaption>**Figure 3.** Overview of the infection frequency (A) and the infection load (B) for each queen age from both naturally selected and traditionally managed colonies. Significant differences are indicated with *. The numbers on top of each bar represent the number of positive samples and the total sample for each group.</figcaption></figure>

**4. Discussion**

In comparison with the previously described SOV phenotyping protocol [49], this study implemented two improvements for virus screening in honey bee eggs. Measuring virus prevalence with qPCR showed, as expected, a higher detection sensitivity than when measuring virus prevalence with *end-point* PCR. The detection threshold of the qPCR was around $10^3$ (30 Ct) for all viruses, while the end-point PCR showed significantly higher detection thresholds for SBV, BQCV and ABPV ($10^8$ for SBV, $10^7$ for BQCV and $10^6$ for ABPV). This implies that phenotyping using end-point PCR underestimated the number of virus infections for these viruses. There was no difference between both detection thresholds for DWV. It should be noted that samples negative on qPCR can still be infected below the detection threshold and that positive samples might be infected with viruses in a dormant state. An important advantage of qPCR is that breeding programs can manually set threshold values based on the breeding goal and the virulence of the virus. Each breeding program can thus determine the degree of positive or negative selection desired. Infection loads in eggs are linked with the infection status of the queen [53,54] and have been shown to reduce virus infections in the colony as a whole [49]. Nevertheless, the impact of different virus infection loads in eggs on subsequent developmental stages is currently unknown. Further research, where eggs with different virus infections are reared in vitro, could improve our understanding of the impact that vertical transmission has on antiviral responses and honey bee health.

The comparison between individual DWV genotypes and the generic DWV showed a large underestimation of DWV infections in this sample cohort when screening for either one of the genotypes or the sum of DWV-A and DWV-B. This can be seen in terms of the underestimation of the number of infected samples (25% of the samples) and the lower
infection loads (on average 1.69 Log₁₀ lower). In comparison, a previous study in the UK found 40% higher DWV titers when screening with a universal DWV-complex assay than when pooling the results of screening with specific DWV-A and DWV-B assays [55]. Possible explanations for this could be the presence of different genotypes, DWV-C or DWV-D, neither of which has to date been detected in Belgium [26,27], or mutations in the primer region that hamper correct primer hybridization [56]. A study on the genetic diversity within a DWV population in a colony showed that 82% of the genome had >1 sequence variant present in the frequency of >1%, and 39% of the genome had >1 sequence variant present in a frequency of >10% [8]. In addition, shifts in the sequence space of the DWV-A quasispecies have been shown after injection in honey bee pupae [57]. The rapid shifts in the DWV quasispecies are consistent with the punctuated evolution theory, whereby the infection of a new host causes a selective sweep, followed by diversification towards an increased genetic heterogeneity that has potentially adapted to the host-specific antiviral defenses [28]. This implies that primer regions, although being in conserved regions, may evolve over time and reduce the primer amplification efficiency.

With regard to the virus prevalence in drone eggs collected across Europe, virus-free samples were found in all countries and in both TMC and NSC groups in the sample cohort studied herein. According to the SOV protocol, queens laying virus-free eggs at the time of sampling were phenotyped as SOV-positive (SOV+) [49]. The low number of infections in the Gotland population (located in Näsudden) and in the TMC (located in Sigarve) populations from Sweden was remarkable, as both groups only had one sample infected with DWV despite multiple studies recording high viral loads in the worker bees of Gotland throughout the years [35,38,44,45,58–60]. The presence of SOV+ queens across Europe serves as a possible starting point for local breeding programs to perform selection within the variation in virus resistance present within honey bee populations [49,61]. Including subsequent generations descending from SOV+ queens in breeding programs is crucial to maximize selection on the heritable genetic contribution [49] behind virus resistance. This is because the SOV trait does not differentiate between samples that are free of viruses due to the virus resistance of the queen or due to other non-genetic circumstances. In this study, virus abundance and prevalence were based on one sampling time point. However, for accurate SOV phenotyping, multiple samples could be collected at different time points to differentiate between queens that temporarily lay virus-free eggs and queens that do so consistently. This is an important consideration, as virus dynamics have been shown to change across the seasons [62].

In addition to the high prevalence of DWV, the second most common virus found in this study was BQCV. This virus is the most common cause of queen larval death [63,64], but it has not been found to cause overt symptoms in queens despite the detection of high infection loads [65]. Viruses of the ABPV complex and SBV have been found in eggs [54,66–68], but they were rarely detected in this study. The higher virulence of BQCV, SBV and ABPV [4,48,69–72] compared to the lower virulence of DWV [69] could explain why they are less likely to be transmitted vertically without causing queen supersedure or colony health issues [73].

Differences in infection patterns between countries can be caused by climatic conditions, by the seasonality of honey bee viruses [34,56,74–80] or by the low number of samples per country. This is reflected in the significant differences in infection frequencies between the spring and summer sampling seasons in this study. Interestingly, infection frequencies with DWV were higher in spring, contrasting with the generally higher infection frequencies reported for adult bees during the summer and autumn [62,75,77,80]. These findings indicate that conclusions based on SOV phenotyping should always take the time of sampling into account when interpreting results and that the SOV status can change between seasons. Other factors affecting virus abundance are nutritional quality and availability [81,82], connectivity between colonies [83], colony demography [84,85], population heterogeneity [86–90], colony management [91], the degree of local adaptation [92], individual and colony-level immune responses [93] and other stress factors (such as exposure
to neonicotinoids) [94]. In this study, worker eggs did not have higher infection loads than drone eggs despite the possible occurrence of trans-spermal virus transmission [25].

By comparing naturally surviving with traditionally managed populations in the same local context, insights can be gained into which evolutionary adaptations are needed for honey bees to survive without treatment against the Varroa mite. Typical for naturally surviving populations is that they harbor higher mite numbers that serve as an important vector for viruses [39,95]. Honey bee colonies react to these high disease pressures with adaptations in their antiviral responses or by forms of social immunity [96]. This study shows that, in this sample cohort, infection frequencies were significantly higher in NSCs than in TMCs for DWV but that infection loads did not differ between the two groups. Honey bee queens appear to avoid increased vertical transmission loads despite increased infection frequencies. Virus loads in worker bees were not studied here. Therefore, it remains uncertain if the higher infection frequencies are a result of increased virus circulation in the naturally surviving populations. The high variability between both groups in each country could be caused by previously mentioned factors affecting virus abundance, the time since colonies were left untreated [36,37] or the degree of genetic divergence between TMCs and NSCs within a country [52,62,90].

Honey bee queens accumulate viral infections and infection loads during queen rearing [97,98], during mating flights [47,99,100] or as they become older [47,51,68], despite increased immune responses [101]. In contrast to what was expected, both the infection frequency of DWV and the infection load in the eggs of queens from TMCs decreased with increasing queen age. This difference is not caused by the mortality of queens with high infection loads, as the distribution of the egg infection loads does not overlap between queens aged 0 years and 1 year. Queens from NSCs showed the same trends in infection loads across all queen ages but did not show differences in infection frequencies. For BQCV, only the infection load was significantly lower in queens aged 0 years than in queens aged 1 year from TMCs. In beekeeping practices, queens are often renewed yearly, as young queens are associated with lower winter mortality [66,102]. This study suggests that older queens from colonies that are treated against the Varroa mite might be able to adapt their antiviral responses to DWV and thereby reduce the infection loads transmitted via their eggs. If so, the frequent renewal of queens could limit this potential as opposed to selecting towards increased queen longevity.

By focusing on the role of honey bee queens, this research adds to the growing literature on the relationship between viral infections and honey bee health. Evolutionary patterns of resistance and tolerance can form the theoretical foundation to incorporate virus resilience in breeding programs. This is a promising perspective, as shown by the variability of vertical transmission over time, across queen ages and under different evolutionary conditions.

**Supplementary Materials:** The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/v14112442/s1, Table S1: Primer sequences; Figure S1: Comparison of PCR and qRT-PCR; Table S2: Additional information regarding European Screening.

**Author Contributions:** Conceptualization and methodology, D.C.B. and D.C.d.G.; formal analysis, D.C.B. and M.B.; resources, D.C.B., B.D., T.B., A.D., D.D., D.E., J.F., A.G. (Alexandru Giurgiuand), A.G. (Aleš Gregorc), J.K., B.L., J.R.d.M., M.O., D.P., M.P. and M.A.P.; writing—original draft preparation, D.C.B.; writing—review and editing, D.C.B., B.D., T.B., A.D., D.D., D.E., J.F., A.G. (Alexandru Giurgiuand), A.G. (Aleš Gregorc), J.K., B.L., J.R.d.M., M.O., D.P., M.P. and M.A.P.; visualization, D.C.B.; supervision, D.C.d.G. and L.D.S.; project administration, D.C.B.; funding acquisition, D.C.d.G. and A.G. (Aleš Gregorc). All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by Research Foundation-Flanders (FWO, research Grant G003021N), by the Slovenian ARRS research grant project N4-0192 and by the European Union's Horizon 2020 research and innovation program under grant agreement No 817622 (B-GOOD).

**Institutional Review Board Statement:** Not applicable.
**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** The data presented in this study is openly available in FigShare at https://doi.org/10.6084/m9 figshare.21485865.

**Acknowledgments:** We would like to thank Maria Bolt, Victor Kohut and Olivier Raujol for their participation in this study, as well as Viktor Deleu for her assistance in the lab. The SLU group would like to thank Åke Lyberg and Meghan Millbrath for help with the management of the colonies and the collection of the samples. UPV/EHU thanks Egoitz Galartza and June Gorrochategui for help with sampling. CIMO/IPB thanks Ana Rita Lopes for help with sampling.

**Conflicts of Interest:** The authors declare no conflict of interest.

## References

1. Nowak, M.A.; Tarnita, C.E.; Wilson, E.O. The evolution of eusociality. *Nature* **2010**, *466*, 1057–1062. [CrossRef] [PubMed]

2. Liu, L.; Zhao, X.Y.; Tang, Q.B.; Lei, C.L.; Huang, Q.Y. The Mechanisms of Social Immunity Against Fungal Infections in Eusocial Insects. *Toxins* **2019**, *11*, 244. [CrossRef] [PubMed]

3. Udiani, O.; Fefferman, N.H. How disease constrains the evolution of social systems. *Proc. R. Soc. B* **2020**, *287*, 20201284. [CrossRef] [PubMed]

4. Sumpter, D.; Martin, S.J. The dynamics of virus epidemics in Varroa -infested honey bee colonies. *J. Anim. Ecol.* **2004**, *73*, 51–63. [CrossRef]

5. Dainat, B.; Neumann, P. Clinical signs of deformed wing virus infection are predictive markers for honey bee colony losses. *J. Invertebr. Pathol.* **2013**, *112*, 278–280. [CrossRef]

6. Wegener, J.; Ruhnke, H.; Scheller, K.; Mispagel, S.; Knollmann, U.; Kamp, G.; Bienefeld, K. Pathogenesis of varroosis at the level of the honey bee (*Apis mellifera*) colony. *Insect Physiol.* **2016**, *91–92*, 1–9. [CrossRef]

7. Gisder, S.; Möckel, N.; Eisenhardt, D.; Genersch, E. In vivo evolution of viral virulence: Switching of deformed wing virus between hosts results in virulence changes and sequence shifts. *Environ. Microbiol.* **2018**, *20*, 4612–4628. [CrossRef]

8. Ryabov, E.V.; Childers, A.K.; Lopez, D.; Grubbs, K.; Posada-Florez, F.; Weaver, D.; Girten, W.; VanEngelsdorp, D.; Chen, Y.; Evans, J.D. Dynamic evolution in the key honey bee pathogen deformed wing virus: Novel insights into virulence and competition using reverse genetics. *PLoS Biol.* **2019**, *17*, e3000502. [CrossRef]

9. Norton, A.M.; Remnant, E.J.; Tom, J.; Buchmann, G.; Blacquiere, T.; Beekman, M. Adaptation to vector-based transmission in a honeybee virus. *J. Anim. Ecol.* **2021**, *90*, 2254–2267. [CrossRef]

10. Neumann, P.; Yañez, O.; Fries, I.; De Miranda, J.R.; de Miranda, J. Varroa invasion and virus adaptation. *Trends Parasitol.* **2012**, *28*, 353–354. [CrossRef]

11. Wilfert, L.; Long, G.; Leggett, H.; Schmid-Hempel, P.; Butlin, R.; Martin, S.J.; Boots, M. Honeybee disease: Deformed wing virus is a recent global epidemic in honeybees driven by Varroa mites. *Science* **2016**, *351*, 594–597. [CrossRef] [PubMed]

12. Traynor, K.S.; Mondet, F.; de Miranda, J.R.; Techer, M.; Kowallik, V.; Oddie, M.A.Y.; Chantawannakul, P.; McAfee, A. Varroa destructor: A Complex Parasite, Crippling Honey Bees Worldwide. *Trends Parasitol.* **2020**, *36*, 592–606. [CrossRef] [PubMed]

13. Mondet, F.; Beaurepaire, A.; McAfee, A.; Locke, B.; Alaux, C.; Blanchard, S.; Danka, B.; Le Conte, Y. Honey bee survival mechanisms against the parasite Varroa destructor: A systematic review of phenotypic and genomic research efforts. *Int. J. Parasitol.* **2020**, *50*, 433–447. [CrossRef] [PubMed]

14. Highfield, A.C.; Nagar, A.E.; Mackinder, L.; Noël, L.; Hall, M.J.; Martin, S.J.; El Nagar, A.; Mackinder, L. Deformed wing virus implicated in overwintering honeybee colony losses. *Appl. Environ. Microbiol.* **2009**, *75*, 7212–7220. [CrossRef]

15. Berthoud, H.; Imdorf, A.; Haueter, M.; Radloff, S.; Neumann, P. Virus infections and winter losses of honey bee colonies (*Apis mellifera*). *J. Apic. Res.* **2010**, *49*, 60–65. [CrossRef]

16. Iqbal, J.; Mueller, U. Virus infection causes specific learning deficits in honeybee foragers. *Proc. R. Soc. B Biol. Sci.* **2007**, *274*, 1517–1521. [CrossRef]

17. Martin, S.J. The role of varroa and viral pathogens in the collapse of honeybee colonies: A modelling approach. *J. Appl. Ecol.* **2001**, *38*, 1082–1093. [CrossRef]

18. Peck, D.T.; Seeley, T.D. Mite bombs or robber lures? The roles of drifting and robbing in Varroa destructor transmission from collapsing honey bee colonies to their neighbors. *PLoS ONE* **2019**, *14*, e0218392. [CrossRef]

19. Barroso-Arévalo, S.; Fernández-Carrión, E.; Goyache, J.; Molero, F.; Puerta, F.; Sánchez-Vizcaíno, J.M. High load of deformed wing virus and Varroa destructor infestation are related to weakness of honey bee colonies in Southern Spain. *Front. Microbiol.* **2019**, *10*, 1331. [CrossRef]

20. Dainat, B.; Evans, J.D.; Chen, Y.P.; Gauthier, L.; Neumann, P. Predictive markers of honey bee colony collapse. *PLoS ONE* **2012**, *7*, e32151. [CrossRef]

21. McMenamin, A.J.; Genersch, E. Honey bee colony losses and associated viruses. *Curr. Opin. Insect Sci.* **2015**, *8*, 121–129. [CrossRef] [PubMed]
22. Beaurepaire, A.; Piot, N.; Doublet, V.; Antunez, K.; Campbell, E.; Chantawannakul, P.; Chejanovsky, N.; Gajda, A.; Heerman, M.; Panziera, D.; et al. Diversity and global distribution of viruses of the western honey bee, *Apis mellifera*. *Insects* **2020**, *11*, 239. [CrossRef] [PubMed]

23. Valles, S.M.; Chen, Y.; Firth, A.E.; Guérin, D.M.A.; Hashimoto, Y.; Herrero, S.; de Miranda, J.R.; Ryabov, E. ICTV Virus Taxonomy Profile: Dicistroviridae. *J. Gen. Virol.* **2017**, *98*, 355. [CrossRef] [PubMed]

24. McMahon, D.P.; Wilfert, L.; Paxton, R.J.; Brown, M.J.F. Emerging Viruses in Bees: From Molecules to Ecology. *Adv. Virus Res.* **2018**, *101*, 251–291. [CrossRef] [PubMed]

25. Yañez, O.; Piot, N.; Dalmon, A.; de Miranda, J.R.; Chantawannakul, P.; Panziera, D.; Amiri, E.; Smagghe, G.; Schroeder, D.; Chejanovsky, N. Bee Viruses: Routes of Infection in Hymenoptera. *Front. Microbiol.* **2020**, *11*, 943. [CrossRef]

26. de Miranda, J.R.; Brettell, L.E.; Chejanovsky, N.; Childers, A.K.; Dalmon, A.; Deboutte, W.; de Graaf, D.C.; Doublet, V.; Gebremedhn, H.; Genersch, E.; et al. Cold case: The disappearance of Egypt bee virus, a fourth distinct master strain of deformed wing virus linked to honeybee mortality in 1970's Egypt. *Virol. J.* **2022**, *19*, 12. [CrossRef]

27. Mordecai, G.J.; Wilfert, L.; Martin, S.J.; Jones, I.M.; Schroeder, D.C. Diversity in a honey bee pathogen: First report of a third master variant of the Deformed Wing Virus quasispecies. *ISME J.* **2016**, *10*, 1264–1273. [CrossRef]

28. Ryabov, E.V.; Wood, G.R.; Fannon, J.M.; Moore, J.D.; Bull, J.C.; Chandler, D.; Mead, A.; Burroughs, N.; Evans, D.J. A Virulent Strain of Deformed Wing Virus (DWV) of Honeybees (*Apis mellifera*) Prevails after Varroa destructor-Mediated, or In Vitro, Transmission. *PLoS Pathog.* **2014**, *10*, e1004230. [CrossRef]

29. Dalmon, A.; Desbiez, C.; Coulon, M.; Thomasson, M.; Le Conte, Y.; Alaux, C.; Vallon, J.; Moury, B. Evidence for positive selection and recombination hotspots in Deformed wing virus (DWV). *Sci. Rep.* **2017**, *7*, srep41045. [CrossRef]

30. Moore, J.; Jironkin, A.; Chandler, D.; Burroughs, N.; Evans, D.J.; Ryabov, E.V. Recombinants between Deformed wing virus and Varroa destructor virus-1 may prevail in Varroa destructor-infested honeybee colonies. *J. Gen. Virol.* **2011**, *92*, 156–161. [CrossRef]

31. Wang, H.; Xie, J.; Shreeve, T.G.; Ma, J.; Pallett, D.W.; King, L.A.; Possee, R.D. Sequence Recombination and Conservation of Varroa destructor Virus-1 and Deformed Wing Virus in Field Collected Honey Bees (*Apis mellifera*). *PLoS ONE* **2013**, *8*, e74508. [CrossRef] [PubMed]

32. Zioni, N.; Soroker, V.; Chejanovsky, N. Replication of varroa destructor virus 1 (VDV-1) and a varroa destructor virus 1-deformed wing virus recombinant (VDV-1-DWV) in the head of the honey bee. *Virology* **2011**, *417*, 106–112. [CrossRef] [PubMed]

33. di Prisco, G.; Zhang, X.; Pennacchio, F.; Caprio, E.; Li, J.; Evans, J.D.; DeGrandi-Hoffman, G.; Hamilton, M.; Chen, Y.; Prisco, G.; et al. Dynamics of Persistent and Acute Deformed Wing Virus Infections in Honey Bees, *Apis mellifera*. *Viruses* **2011**, *3*, 2425–2441. [CrossRef]

34. Glenny, W.; Cavigli, I.; Daughenbaugh, K.F.; Radford, R.; Kegley, S.E.; Flenniken, M.L. Honey bee (*Apis mellifera*) colony health and pathogen composition in migratory beekeeping operations involved in California almond pollination. *PLoS ONE* **2017**, *12*, e0128214. [CrossRef] [PubMed]

35. Thaduri, S.; Marupakula, S.; Terenius, O.; Onorati, P.; Tellgren-Roth, C.; Locke, B.; de Miranda, J.R. Global similarity, and some key differences, in the metagenomes of Swedish varroa-surviving and varroa-susceptible honeybees. *Sci. Rep.* **2021**, *11*, 23214. [CrossRef] [PubMed]

36. Moro, A.; Blacquière, T.; Panziera, D.; Dietemann, V.; Neumann, P. Host-Parasite Co-Evolution in Real-Time: Changes in Honey Bee Resistance Mechanisms and Mite Reproductive Strategies. *Insects* **2021**, *12*, 120. [CrossRef]

37. Moro, A.; Blacquière, T.; Dahle, B.; Dietemann, V.; Le Conte, Y.; Locke, B.; Neumann, P.; Beaurepaire, A. Adaptive population structure shifts in invasive parasitic mites, Varroa destructor. *Ecol. Evol.* **2021**, *11*, 5937–5949. [CrossRef]

38. Thaduri, S.; Stephan, J.G.; de Miranda, J.R.; Locke, B. Disentangling host-parasite-pathogen interactions in a varroa-resistant honeybee population reveals virus tolerance as an independent, naturally adapted survival mechanism. *Sci. Rep.* **2019**, *9*, 6221. [CrossRef]

39. Locke, B. Natural Varroa mite-surviving *Apis mellifera* honeybee populations. *Apidologie* **2016**, *47*, 467–482. [CrossRef]

40. Le Conte, Y.; Meixner, M.D.; Brandt, A.; Carreck, N.L.; Costa, C.; Mondet, F.; Büchler, R. Geographical distribution and selection of european honey bees resistant to varroa destructor. *Insects* **2020**, *11*, 873. [CrossRef]

41. Thompson, J.N. *The Geographic Mosaic of Coevolution*; University of Chicago Press: Chicago, IL, USA, 2005.

42. Locke, B.; Fries, I. Characteristics of honey bee colonies (*Apis mellifera*) in Sweden surviving Varroa destructor infestation. *Apidologie* **2011**, *42*, 533–542. [CrossRef]

43. Blacquière, T.; Boot, W.; Calis, J.; Moro, A.; Neumann, P.; Panziera, D. Darwinian black box selection for resistance to settled invasive Varroa destructor parasites in honey bees. *Biol. Invasions* **2019**, *21*, 2519–2528. [CrossRef]

44. Locke, B.; Thaduri, S.; Stephan, J.G.; Low, M.; Blacquière, T.; Dahle, B.; Le Conte, Y.; Neumann, P.; de Miranda, J.R. Adapted tolerance to virus infections in four geographically distinct Varroa destructor-resistant honeybee populations. *Sci. Rep.* **2021**, *11*, 12359. [CrossRef] [PubMed]

45. Locke, B.; Forsgren, E.; De Miranda, J.R. Increased tolerance and resistance to virus infections: A possible factor in the survival of Varroa destructor-resistant honey bees (*Apis mellifera*). *PLoS ONE* **2014**, *9*, e99998. [CrossRef] [PubMed]

46. Amiri, E.; Strand, M.K.; Rueppell, O.; Tarpy, D.R. Queen quality and the impact of honey bee diseases on queen health: Potential for interactions between two major threats to colony health. *Insects* **2017**, *8*, 48. [CrossRef]

47. Gauthier, L.; Ravallec, M.; Tournaire, M.; Cousserans, F.F.; Bergoin, M.; Dainat, B.; de Miranda, J.R. Viruses Associated with Ovarian Degeneration in *Apis mellifera* L. Queens. *PLoS ONE* **2011**, *6*, e16217. [CrossRef] [PubMed]
48. Martin, S.J.; Highfield, A.; Brettell, L.; Villalobos, E.; Budge, G.E.; Powell, M.; Nikaido, S.; Schroeder, D.C. Global honey bee viral landscape altered by a parasitic mite. *Science* **2012**, *336*, 1304–1306. [CrossRef]

49. de Graaf, D.C.; Laget, D.; De Smet, L.; Claeys Bouúuaert, D.; Brunain, M.; Veerkamp, R.F.; Brascamp, E.W. Heritability estimates of the novel trait 'suppressed in ovo virus infection' in honey bees (*Apis mellifera*). *Sci. Rep.* **2020**, *10*, 14310. [CrossRef]

50. Claeys Bouuaert, D.; De Smet, L.; de Graaf, D.C. Breeding for virus resistance and its effects on deformed wing virus infection patterns in honey bee queens. *Viruses* **2021**, *13*, 1074. [CrossRef]

51. Woyke, J. Drone Larvae from Fertilized Eggs of the Honeybee. *J. Apic. Res.* **2015**, *2*, 19–24. [CrossRef]

52. Chen, Y.; Pettis, J.S.; Corona, M.; Chen, W.P.; Li, C.J.; Spivak, M.; Visscher, P.K.; DeGrandi-Hoffman, G.; Boncristiani, H.; Zhao, Y.; et al. Israeli Acute Paralysis Virus: Epidemiology, Pathogenesis and Implications for Honey Bee Health. *PLoS Pathog.* **2014**, *10*, e1004261. [CrossRef] [PubMed]

53. Shen, M.; Cui, L.; Ostiguy, N.; Cox-Foster, D. Intricate transmission routes and interactions between picorna-like viruses (Kashmir bee virus and sacbrood virus) with the honeybee host and the parasitic varroa mite. *J. Gen. Virol.* **2005**, *86*, 2281–2289. [CrossRef] [PubMed]

54. Amiri, E.; Kryger, P.; Meixner, M.D.; Strand, M.K.; Tarpy, D.R.; Rueppell, O. Quantitative patterns of vertical transmission of deformed wing virus in honey bees. *PLoS ONE* **2018**, *13*, e0195283. [CrossRef]

55. Bradford, E.L.; Christie, C.R.; Campbell, E.M.; Bowman, A.S. A real-time PCR method for quantification of the total and major variant strains of the deformed wing virus. *PLoS ONE* **2017**, *12*, e0190017. [CrossRef]

56. De Miranda, J.R.; Bailey, L.; Ball, B.V.; Blanchard, P.; Budge, G.E.; Chejanovsky, N.; Chen, Y.; Gauthier, L.; Genersch, E.; De Graaf, D.C.; et al. Standard methods for virus research in *Apis mellifera*. *J. Apic. Res.* **2013**, *52*, 1–56. [CrossRef]

57. Yañez, O.; Chávez-Galarza, J.; Tellgren-Roth, C.; Pinto, M.A.; Neumann, P.; de Miranda, J.R. The honeybee (*Apis mellifera*) developmental state shapes the genetic composition of the deformed wing virus-A quasispecies during serial transmission. *Sci. Rep.* **2020**, *10*, 5956. [CrossRef] [PubMed]

58. Locke, B.; Forsgren, E.; Fries, I.; De Miranda, J.R. Acaricide Treatment Affects Viral Dynamics in Varroa destructor-Infested Honey Bee Colonies via both Host Physiology and Mite Control. *Appl. Env. Microbiol.* **2012**, *78*, 227–235. [CrossRef]

59. Thaduri, S.; Locke, B.; Granberg, F.; De Miranda, J.R. Temporal changes in the viromes of swedish varroa-resistant and varroa-susceptible honeybee populations. *PLoS ONE* **2018**, *13*, e0206938. [CrossRef]

60. Locke, B.; Semberg, E.; Forsgren, E.; De Miranda, J.R. Persistence of subclinical deformed wing virus infections in honeybees following Varroa mite removal and a bee population turnover. *PLoS ONE* **2017**, *12*, e0180910. [CrossRef]

61. Bhatia, S.; Baral, S.; Vega Melendez, C.; Amiri, E.; Rueppell, O. Comparing survival of israeli acute paralysis virus infection among stocks of U.S. honey bees. *Insects* **2021**, *12*, 60. [CrossRef]

62. Faurot-Daniels, C.; Glenny, W.; Daughenbaugh, K.F.; McMenamin, A.J.; Burkle, L.A.; Flenniken, M.L. Longitudinal monitoring of honey bee colonies reveals dynamic nature of virus abundance and indicates a negative impact of Lake Sinai virus 2 on colony health. *PLoS ONE* **2020**, *15*, e0237544. [CrossRef] [PubMed]

63. Harizanis, P. Infestation of queen cells by the mite Varroa-Jacobsoni. *Apidologie* **1991**, *22*, 533–538. [CrossRef]

64. Anderson, D.; Gibbs, A.J. Inapparent Virus Infections and their Interactions in Pupae of the Honey Bee (*Apis mellifera* Linnaeus) in Australia. *J. Gen. Virol.* **1988**, *69*, 1617–1625. [CrossRef]

65. Chen, Y.; Pettis, J.S.; Collins, A.; Feldlaufer, M.F. Prevalence and transmission of honeybee viruses. *Appl. Environ. Microbiol.* **2006**, *72*, 606–611. [CrossRef]

66. Gray, A.; Adjlane, N.; Arab, A.; Ballis, A.; Brusbardis, V.; Charrière, J.D.; Chlebo, R.; Coffey, M.F.; Cornelissen, B.; Amaro da Costa, C.; et al. Honey bee colony winter loss rates for 35 countries participating in the COLOSS survey for winter 2018–2019, and the effects of a new queen on the risk of colony winter loss. *Apidologie* **2020**, *59*, 744–751. [CrossRef]

67. Ravoet, J.; Maharramov, J.; Meeus, I.; De Smet, L.; Wenseleers, T.; Smagghe, G.; de Graaf, D.C. Comprehensive Bee Pathogen Screening in Belgium Reveals Crithidia mellificae as a New Contributory Factor to Winter Mortality. *PLoS ONE* **2013**, *8*, e72443. [CrossRef]

68. Francis, R.M.; Nielsen, S.L.; Kryger, P. Patterns of viral infection in honey bee queens. *J. Gen. Virol.* **2013**, *94*, 668–676. [CrossRef]

69. de Miranda, J.R.; Genersch, E. Deformed wing virus. *J. Invertebr. Pathol.* **2010**, *103*, S48–S61. [CrossRef]

70. Genersch, E.; Aubert, M. Emerging and re-emerging viruses of the honey bee ( *Apis mellifera* L.). *Vet. Res.* **2010**, *41*, 54. [CrossRef]

71. de Miranda, J.R.; Fries, I. Venereal and vertical transmission of deformed wing virus in honeybees (*Apis mellifera* L.). *J. Invertebr. Pathol.* **2008**, *98*, 184–189. [CrossRef]

72. Li, J.; Wang, T.; Evans, J.D.; Rose, R.; Zhao, Y.; Li, Z.; Li, J.; Huang, S.; Heerman, M.; Rodríguez-García, C.; et al. The phylogeny and pathogenesis of sacbrood virus (SBV) infection in European honey bees, *Apis mellifera*. *Viruses* **2019**, *11*, 61. [CrossRef] [PubMed]

73. Sandrock, C.; Tanadini, M.; Tanadini, L.G.; Fauser-Misslin, A.; Potts, S.G.; Neumann, P. Impact of chronic neonicotinoid exposure on honeybee colony performance and queen supersedure. *PLoS ONE* **2014**, *9*, e103592. [CrossRef]

74. Di Prisco, G.; Pennacchio, F.; Caprio, E.; Boncristiani, H.F.; Evans, J.D.; Chen, Y.; Prisco, G.D.; Pennacchio, F.; Caprio, E.; Boncristiani, H.F.; et al. Varroa destructor is an effective vector of Israeli acute paralysis virus in the honeybee, *Apis mellifera*. *J. Gen. Virol.* **2011**, *92*, 151–155. [CrossRef] [PubMed]

75. Traynor, K.S.; Rennich, K.; Forsgren, E.; Rose, R.; Pettis, J.; Kunkel, G.; Madella, S.; Evans, J.; Lopez, D.; VanEngelsdorp, D. Multiyear survey targeting disease incidence in US honey bees. *Apidologie* **2016**, *47*, 325–347. [CrossRef]

76. Allsopp, M.H. Analysis of Varroa Destructor Infestation of Southern African Honeybee Populations. Ph.D. Thesis, University of Pretoria, Pretoria, South Africa, 2006.
77. Tentcheva, D.; Gauthier, L.; Zappulla, N.; Dainat, B.; Cousserans, F.; Colin, M.E.; Bergoin, M. Prevalence and seasonal variations of six bee viruses in *Apis mellifera* L. and Varroa destructor mite populations in France. *Appl. Environ. Microbiol.* **2004**, *70*, 7185–7191. [CrossRef] [PubMed]

78. Hou, C.; Rivkin, H.; Slabezki, Y.; Chejanovsky, N. Dynamics of the presence of israeli acute paralysis virus in honey bee colonies with colony collapse disorder. *Viruses* **2014**, *6*, 2012–2027. [CrossRef] [PubMed]

79. Blažytė-Čereškienė, L.; Skrodenytė-Arbačiauskienė, V.; Radžiutė, S.; Čepulytė-Rakauskienė, R.; Apšegaitė, V.; Būda, V. A three-year survey of honey bee viruses in Lithuania. *J. Apic. Res.* **2016**, *55*, 176–184. [CrossRef]

80. Porrini, C.; Mutinelli, F.; Bortolotti, L.; Granato, A.; Laurenson, L.; Roberts, K.; Gallina, A.; Silvester, N.; Medrzycki, P.; Renzi, T.; et al. The status of honey bee health in Italy: Results from the nationwide bee monitoring network. *PLoS ONE* **2016**, *11*, e0155411. [CrossRef]

81. Walton, A.; Toth, A.L.; Dolezal, A.G. Developmental environment shapes honeybee worker response to virus infection. *Sci. Rep.* **2021**, *11*, 13961. [CrossRef]

82. Dolezal, A.G.; Carrillo-Tripp, J.; Judd, T.M.; Allen Miller, W.; Bonning, B.C. Interacting stressors matter: Diet quality and virus infection in honeybee health. *R. Soc. OPEN Sci.* **2019**, *6*, 181803. [CrossRef]

83. Seeley, T.D.; Smith, M.L. Crowding honeybee colonies in apiaries can increase their vulnerability to the deadly ectoparasite Varroa destructor. *Apidologie* **2015**, *46*, 716–727. [CrossRef]

84. Brosi, B.J.; Delaplane, K.S.; Boots, M.; De Roode, J.C. Ecological and evolutionary approaches to managing honeybee disease. *Nat. Ecol. Evol.* **2017**, *1*, 1250–1262. [CrossRef]

85. Lloyd-Smith, J.O.; Cross, P.C.; Briggs, C.J.; Daugherty, M.; Getz, W.M.; Latto, J.; Sanchez, M.S.; Smith, A.B.; Swei, A. Should we expect population thresholds for wildlife disease? *Trends Ecol. Evol.* **2005**, *20*, 511–519. [CrossRef] [PubMed]

86. Keesing, F.; Belden, L.K.; Daszak, P.; Dobson, A.; Harvell, C.D.; Holt, R.D.; Hudson, P.; Jolles, A.; Jones, K.E.; Mitchell, C.E.; et al. Impacts of biodiversity on the emergence and transmission of infectious diseases. *Nature* **2010**, *468*, 647–652. [CrossRef] [PubMed]

87. Civitello, D.J.; Cohen, J.; Fatima, H.; Halstead, N.T.; Liriano, J.; McMahon, T.A.; Ortega, C.N.; Sauer, E.L.; Sehgal, T.; Young, S.; et al. Biodiversity inhibits parasites: Broad evidence for the dilution effect. *Proc. Natl. Acad. Sci. USA* **2015**, *112*, 8667–8671. [CrossRef]

88. Zhu, Y.; Chen, H.; Fan, J.X.; Wang, Y.; Li, Y.; Chen, J.; Fan, J.X.; Yang, S.; Hu, L.; Leung, H.; et al. Genetic diversity and disease control in rice. *Nature* **2000**, *406*, 718–722. [CrossRef]

89. Mundt, C.C. Use of multiline cultivars and cultivar mixtures for disease management. *Annu. Rev. Entomol.* **2002**, *40*, 381–410. [CrossRef]

90. Beaurepaire, A.; Sann, C.; Arredondo, D.; Mondet, F.; Le Conte, Y. Behavioral Genetics of the Interactions between *Apis mellifera* and Varroa destructor. *Insects* **2019**, *10*, 299. [CrossRef]

91. Bartlett, L.J.; Boots, M.; Brosi, B.J.; De Roode, J.C.; Delaplane, K.S.; Hernandez, C.A.; Wilfert, L. Persistent effects of management history on honeybee colony virus abundances. *J. Invertebr. Pathol.* **2021**, *179*, 107520. [CrossRef]

92. Francis, R.M.; Amiri, E.; Meixner, M.D.; Kryger, P.; Gajda, A.; Andonov, S.; Uzunov, A.; Topolska, G.; Charistos, L.; Costa, C.; et al. Effect of genotype and environment on parasite and pathogen levels in one apiary—A case study. *J. Apic. Res.* **2015**, *53*, 230–232. [CrossRef]

93. Brutscher, L.M.; Daughenbaugh, K.F.; Flenniken, M.L. Antiviral defense mechanisms in honey bees. *Curr. Opin. Insect Sci.* **2015**, *10*, 71–82. [CrossRef]

94. Di Prisco, G.; Cavaliere, V.; Annoscia, D.; Varricchio, P.; Caprio, E.; Nazzi, F.; Gargiulo, G.; Pennacchio, F. Neonicotinoid clothianidin adversely affects insect immunity and promotes replication of a viral pathogen in honey bees. *Proc. Natl. Acad. Sci. USA* **2013**, *110*, 18466–18471. [CrossRef] [PubMed]

95. Coulon, M.; Dalmon, A.; Di Prisco, G.; Prado, A.; Arban, F.; Dubois, E.; Ribière-Chabert, M.; Alaux, C.; Thiéry, R.; Le Conte, Y. Interactions Between Thiamethoxam and Deformed Wing Virus Can Drastically Impair Flight Behavior of Honey Bees. *Front. Microbiol.* **2020**, *11*, 766. [CrossRef] [PubMed]

96. Kurze, C.; Routtu, J.; Moritz, R.F.A. Parasite resistance and tolerance in honeybees at the individual and social level. *Zoology* **2016**, *119*, 290–297. [CrossRef] [PubMed]

97. Gregorc, A.; Bakonyi, T. Viral infections in queen bees (*Apis mellifera* carnaica) from rearing apiaries. *Acta Vet. Brno* **2012**, *81*, 15–19. [CrossRef]

98. Žvokelj, L.; Bakonyi, T.; Korošec, T.; Gregorc, A. Appearance of acute bee paralysis virus, black queen cell virus and deformed wing virus in Carnolian honey bee (*Apis mellifera* carnaica) queen rearing. *J. Apic. Res.* **2020**, *59*, 53–58. [CrossRef]

99. Amiri, E.; Meixner, M.D.; Kryger, P. Deformed wing virus can be transmitted during natural mating in honey bees and infect the queens. *Sci. Rep.* **2016**, *6*, 33065. [CrossRef]

100. Gregorc, A.; Smodiš Škerl, M.I. Characteristics of honey bee (*Apis mellifera* carnaica, Pollman 1879) queens reared in Slovenian commercial breeding stations. *J. Apic. Sci.* **2015**, *59*, 5–12. [CrossRef]

101. Chaimanee, V.; Chantawannakul, P.; Chen, Y.; Evans, J.D.; Pettis, J.S. Effects of host age on susceptibility to infection and immune gene expression in honey bee queens (*Apis mellifera*) inoculated with Nosema ceranae. *Apidologie* **2014**, *45*, 451–463. [CrossRef]

102. Van Der Zee, R.; Gray, A.; Holzmann, C.; Pisa, L.; Brodschneider, R.; Chlebo, R.; Coffey, M.F.; Kence, A.; Kristiansen, P.; Mutinelli, F.; et al. Standard survey methods for estimating colony losses and explanatory risk factors in *Apis mellifera*. *J. Apic. Res.* **2013**, *52*, 1–36. [CrossRef]