An apiary was established at the Baton Rouge, LA, USDA Honey Bee Breeding, Genetics and Physiology Research Laboratory (30°22′56″N, 91°10′40″W) with nine colonies (three from each of the following stocks: Italian, Pol-Line, and Russian) using packaged bees (~0.90 kg or 7,000 worker bees per package). Pol-Line and Russian queens were sourced from the same USDA lab; mite-susceptible Italian queens were sourced from a commercial queen producer. Packages were installed on 23 April 2019 and supplied with sugar water for the first month via in-hive feeders. Colonies were equalized as necessary on 13 May and 5 June with additional brood and food (combination of pollen and nectar) frames. Colonies were sampled a minimum 6 weeks after the last brood supplement to allow time for the worker populations to reflect queen genetics. Colony and pollen forager sampling was completed once per month during the weeks of 22 July (time 1), 19 August (time 2), and 9 September 2019 (time 3). Nectar foragers were only collected for time 1 to streamline later collection periods.

Colony-level pollen preferences were determined using pollen traps for 4 days, with the pollen collected daily, mixed with previous days' samples from the same colony for that sampling period, and stored at −20°C for nutritional analyses. Throughout the course of study observation, colonies were not treated with miticide in order to obtain an accurate representation of naturally occurring Varroa mite levels. To estimate Varroa mite levels per 100 bees, ~300 nurse bees were sampled from brood frames per colony (stored at −20°C), washed using detergent, then bees and mites were counted (74, 75). A pool of 50 adult bees per colony was simultaneously sampled from two brood frames then stored at −80°C for colony viral analyses.

To sample nectar and pollen foragers, colony entrances were closed between 09:00 and 11:00, and returning foragers were collected using scintillation vials (1 forager/vial). Each vial was placed on ice for 2 min to reduce movement and prevent swallowing of the honey stomach contents. Pollen loads were removed from the corbiculae of each bee bearing pollen (N = 24/colony/time) and stored at −20°C for nutritional analyses. For each nectar forager (N = 24/colony), honey stomach content was collected into a pre-weighed microcapillary tube by gently squeezing the abdomen. Nectar load was assessed via weight measurements of the microcapillary tubes after sample collection (76). Microcapillary tube contents were then exuded onto a digital refractometer (MISCO Palm Abbe Digital Fluid Refractometer, Solon, OH, USA) for analysis of sugar content (Brix) (77). Due to maximum refractometer detection limitations, nectar loads >85% sugar content were designated as >85%. All collected pollen and nectar foragers were stored individually at −80°C for subsequent forager viral analysis.

For colony-level pollen nutritional analyses, 10 g of pollen collected via pollen traps from each colony for each time point were sent to the Agricultural Chemistry Laboratory (Louisiana State University, AgCenter, Baton Rouge, LA, USA) for analyses of crude protein (5 g) and fat content (5 g). For forager-level pollen analyses, pollen pellets (from both corbiculae of the same individual bee) were weighed and then one pellet was analyzed for protein content and one for lipid content. For protein analysis, pellets were individually homogenized using a handheld pestle, vortexed, and analyzed for total soluble protein using a Bradford Assay (Bio-Rad DC protein assay kit). Protein standards (0.28, 0.24, 0.2, 0.16, 0.12, 0.08, 0.04, and 0 mg/mL) were created using BSA (Bio-Rad). Five microliters of each sample and standard in triplicate were combined with 245 μL of Bradford reagent (Bio-Rad), incubated for 5 min, and analyzed using a spectrophotometer (SpectraMax Plus Microplate Reader) at 595 λ. For analysis of lipid content, pellets were individually dried in a desiccator at room temperature for 24 h, weighed, and analyzed using a modified chloroform-extraction (78, 79). Samples were vortexed with 0.2 mL 2% sodium sulfate. One milliliter of chloroform-methanol (2:1) was added, vortexed, then centrifuged 2,180 G for 5 min at room temperature. Three hundred microliters of the layer below the supernatant was combined with 0.6 mL of deionized water, vortexed, centrifuged again, and incubated at 90°C for 20 min. Three hundred microliters of 10 N sulfuric acid (Thermo Fisher Scientific, Waltham, MA, USA) was added, then samples were incubated at 90°C for 20 min followed by a 2-min ice bath for samples to reach room temperature. Hundred microliters of each sample was read on a spectrophotometer (SpectraMax Plus Microplate Reader) at 540 λ. Samples were compared against canola oil standards (0.20, 0.18, 0.15, 0.12, 0.09, 0.06, 0.03, and 0.015 mg/mL).

For colony-level viral analyses, nurse bee pools (N = 27 samples each, representing a pool of 50 bees/colony to generally approximate colony viral loads) were placed in 50 mL homogenization vials pre-fitted with ceramic beads then homogenized using a bead mill (Omni-Inc Bead Ruptor Elite, Kennesaw, GA, USA). For forager-level viral analyses, a subset of bees (N = 378 foragers total) from each colony and time was selected based on being at the high (N = 162 bees, N = 6 bees/colony/time) and low (N = 162 bees, N = 6 individual bees/colony/time) range of individually collected pollen protein levels or the high (N = 27 bees, N = 3 individual bees/colony) and low (N = 27 bees, N = 3 individual bees/colony) range of sugar content (Brix). Individual bees were homogenized by hand at 5°C in a metal bead bucket using sterilized plastic pestles. Following homogenization, 200 μL of Promega Lysis buffer and 200 μL of Promega Homogenization solution (Promega Corporation, Madison, Wisconsin, USA) were added to each sample. All samples were briefly vortexed then centrifuged for 10 min at 4°C at 14,000 rpm. Total sample RNA was then extracted from 400 μL cleared lysate using the Maxwell RSC 48 simplyRNA cartridges per manufacturer recommended protocol (Promega Corporation, Madison, Wisconsin, USA). RNA was stored in 0.6 mL elution tubes wrapped in parafilm (Bemis NA, Neenah, Wisconsin, USA) at −80°C until cDNA synthesis.

Before cDNA synthesis, frozen RNA samples were thawed on 5°C metal beads, briefly vortexed, then centrifuged. Each RNA sample was quantified via spectrophotometry (NanoDrop One, Thermo-Fisher Scientific Inc., Waltham, Massachusetts, USA) twice using 1 μL of sample each time. The mean ng/μL NanoDrop One readings were calculated per sample then used to determine the volume of sample RNA template and nuclease-free water required to reach a sample concentration of 100 ng of RNA. cDNA was then synthesized using 0.2 mL PCR strip tubes in two steps using Qiagen Quantitect Reverse Transcriptase kits (Thermo-Fisher Scientific Inc., Waltham, Massachusetts, USA). For step one, 2 μL of gDNA wipeout was added to the mix of RNA and water for a total reaction volume of 14 μL per sample. Samples were incubated per manufacturer protocol at 42°C for 2 min in a Bio-Rad T100 Thermal Cycler (Bio-Rad, Hercules, California, USA). Samples were then briefly vortexed and centrifuged before the addition of 6 μL of Step 2 Master Mix consisting of 4 μL 5X Buffer, 1 μL of RT Primer mix, and 1 μL of RT enzyme per sample. Samples were again briefly vortexed and centrifuged then placed into the Bio-Rad T100 Thermal Cycler (42°C for 25 min then 95°C for 3 min) per manufacturer protocol. Transcribed cDNA was in strips tubes were wrapped in parafilm and stored at −80°C until RT-qPCR.

For analyses of both colony and forager viral profiles, samples were analyzed for eight viruses (primers in Table 1): Acute bee paralysis virus (ABPV), Black queen cell virus (BQCV), Chronic bee paralysis virus (CBPV), Deformed wing virus genotype A (DWV-A), Deformed wing virus genotype B (DWV-B), Israeli acute paralysis virus (IAPV), Kashmir bee virus (KBV), and Lake Sinai virus (LSV), following well-established protocols (85–87). The reference gene β-actin was used to ensure sample quality (88). Each sample was replicated three times per primer pair for RT-qPCR analyses. All RT-qPCR reactions consisted of 5 μL SsoFast Universal SYBR Green supermix (Bio-Rad, Hercules, California, USA), 3 μL nuclease-free water (Promega Corp., Madison, WI), 0.5 μL forward primer, 0.5 μL reverse primer, and 1 μL cDNA from each sample. All reactions were run in Bio-Rad CFX96 or CFX Connect Thermal Cyclers (Bio-Rad, Hercules, California, USA) with all reactions of a specific primer occurring in the same machine. The PCR thermal cycling protocol for the DWV-A and CBPV primer pairs was 95°C for 5 min followed by 40 cycles of 95°C for 5 s and 53.5°C for 10 s then 72°C for 10 s; while the protocol for ABPV, β-actin, BQCV, DWV-B, KBV, and LSV was 95°C for 5 min followed by 40 cycles of 95°C for 5 s and 52.5°C for 10 s then 72°C for 10 s. The PCR cycling protocol for the IAPV primer pairs was 95°C for 5 min followed by 40 cycles of 95°C for 5 s and 53.5°C for 10 s then 72°C for 10 s. The thermal protocols included a melt-curve dissociation analysis to confirm product size. BQCV, CBPV, DWV-A, and DWV-B results were quantified using the Standard Curve Method using linearized plasmid constructs for pooled colony samples as well as individual nectar and pollen foragers. Quantified virus levels were log transformed for analyses and graphical representation; all other viruses were counted as positive for any cycle threshold (Ct) value registered at <40 cycles.

All analyses were conducted in JMP Pro 16.0.0 Pro (SAS Institute, Cary, NC, USA). Forager viral loads were compared to those of their respective colonies with the Kruskal-Wallis test. General associations among variables were determined using Spearman's rank correlation coefficient. To determine contribution of individual and colony level virus loads to individual foraging (% protein, % lipid, protein to lipid ratio (P:L), and nectar weight], we used the generalized regression function with a gaussian distribution. Independent variables included sampling month, Varroa mites per 100 bees (mites), and virus panel data [the total number of viruses (Virus), BQCV levels, CBPV levels, DWV-A levels, and DWV-B levels] for both the foragers and the colony as a whole. The variables representing the presence/absence of ABPV, BQCV, CBPV, IAPV, LSV, and KBV were not included as variation was not large enough to be meaningful, but these data can be found in the Supplementary Material. Due to the 85% cutoff for Brix measurements, nectar sugar content was grouped into low (0–33.0%), medium (33.1–66.0%), and high (66.1–85%) categories. The contribution of virus panel data to sugar content was evaluated using an ordered logit model in the generalized regression function using the same independent variables as above. Bee stock was added to some figures for demonstration purposes but was not included in the models due to small sample size per stock.

The relative contribution of forager and colony disease load to individual forager selection was initially evaluated by using model selection criteria (-log likelihood, AICc, BIC, and R²). The full model (as described above) was compared against models including only forager and only colony virus panel data. The relative contribution of colony vs. individual forager viral infections to individual foraging decisions was further analyzed using structural equation models (SEMs). In SEMs, measured variables (manifest variables) can be combined into a larger representative latent variable that describes some larger, unmeasured factor (89–91). Based on our a priori hypothesis that individual and colony level disease can contribute differentially to foraging preferences, we contrasted two sets of models for nectar and pollen foraging. Individual forager (i) and colony (c) manifest variables (virus number and log-transformed BQCV, CBPV, DWV-A, and DWV-B levels) were combined into two respective latent factors representing viral infection of the two sources (Forager and Colony). To account for potential interactions between observed colony and individual viral infection, all models allowed Forager and Colony latent factors to covary without specifying directionality. For nectar foragers, we modeled both Forager and Colony disease latent factors with direct relationships to a Nectar latent factor (nectar weight and sugar content manifest variables). We also repeated the nectar SEMs with nectar weight and sugar content manifest variables individually. As these models were not a good fit for nectar, we only present them in the Supplementary Materials. For pollen foragers, we modeled time as directly related to both latent disease factors and modeled both latent disease factors with direct relationships to the Pollen latent factor (% protein and % lipid manifest variables). We then modeled the latent disease factors against the individual pollen manifest variables. Model fit was assessed with a combination of the comparative fit index (CFI), root mean square error of approximation (RMSEA) with priority given to models closest to optimal values of CFI (>0.9) and RMSEA (<0.1) (92, 93). Models were estimated using maximum likelihood with 1,000 iterations. All SEM results are presented as standardized estimates to directly compare the contributions of colony and individual forager viral infection to foraging decisions.

The range of bee viruses and foraged pollen macronutrient profiles from our study generally matched that of previous studies. Most colonies tested positive for a suite of viruses—ABPV, BQCV, CBPV, DWV-A, DWV-B, IAPV, and LSV—similar to previous studies (97–99). The overall range of colony-foraged pollen protein (Range = 10.5–18.4%) and lipid contents (Range = 0.7–3.14%) were within the range of prior work on US honey bees (100). The variation in individually foraged pollen protein (Range 0.03–33.95%) and lipid contents (Range = 0.095–15.384%) was greater than in colony-collected pollen. This is expected as each forager sample may represent pollen from a single plant species targeted by that individual on a single day (101, 102), while the colony sample was an amalgamation of individual samples over 4 days. Given that all colonies were based within the same apiary, colonies had an equal opportunity to forage from the same floral resources (103) and differences were due in colony-level variation such as viral infection (6, 46–48).

The first goal of this work was to determine if Varroa mites and viral loads were correlated with differences in pollen macronutrient (protein and lipid) and nectar foraging. Prior work has shown that high-protein pollen can increase honey bee survival when faced with infection while high-lipid diets may increase susceptibility to pathogens (41, 57, 72). Therefore, we expected that higher mite and viral loads would generally increase foraging on higher-protein pollen. With the Spearman's correlation we found that higher mite numbers were positively associated with pollen protein content, but no association was found using generalized regression. However, the mite loads in our colonies remained relatively low throughout the study and may not have been variable enough to obtain the most reliable data in this regard. The viruses differed in their correlation with individually foraged pollen protein and lipid contents based on the type of virus and if the individual or colony was infected. The correlations and regressions indicated that colony DWV-B levels were negatively associated with the pollen P:L ratio and positively with lipid content, whereas the total number of viruses at the colony level was negatively associated with lipid content and positively with P:L ratios.

In addition to altering pollen foraging preferences, microbial infections have the potential to change bees' ability to perceive sugar content and nectar foraging behaviors (37, 67, 104, 105). While we expected that higher viral loads would also be associated with an increase in nectar sugar content, we found that this was only the case with colony DWV-B infections. Nectar load weights also significantly increased with higher colony BQCV levels and marginally increased with higher colony DWV-B levels, another indication that different stressors may induce differential foraging (53, 56). However, individual forager DWV-B infections tended (though marginally) to decrease nectar load weight, indicating some tradeoff of individual ability and colony health (47, 67).

These foraging interactions based on virus type suggest differential immune responses and subsequent nutritional requirements in response to different viruses (106–109). For instance, honey bee transcriptional and immune responses differ greatly from bacteria and Nosema to viral infections and between viral infections such as IAPV and DWV (110, 111). Further, energetic costs to honey bees differs with parasite or pathogen identity (104, 108). Work by Annoscia et al. (38) indicated that bees infested with Varroa mites, often associated with DWV levels, and fed pollen diets had up-regulated genes related to lipid metabolism; whereas, Rutter et al. (112) indicated that IAPV interacted with diet via carbohydrate metabolism-related genes. The impact of diet on immune function in the face of viral stressors may also vary with time of year and gut microbiota community structure (113–115).

Not only does diet alter honey bee ability to overcome viral challenges, but it also alters honey bee physiological responses to another stressor—pesticides. Individual bees fed lower P:L diets exhibited longer lifespans than those fed high P:L diets when both sets were exposed to pesticides, potentially due to these diets also altering expression of genes like Vg (vitellogenin) and Defensin-1 (49). Viruses and Varroa infestation can interact with pesticide exposure to make bees more susceptible to the other stressor; though the strength of these interactions depends on the combination of pesticide and virus (116). Varroa mite infestation and DWV infection increase bee susceptibility to insecticides and vice versa (117, 118). Similarly, thiamethoxam exposure increased CBPV loads (119), and thiacloprid increased BQCV viral loads in honey bees (120). The interactions of pesticides and viral infection may be due to how both types of stressors impact bee immune response and gene expression (116, 121). This is demonstrated particularly well by how winter vs. summer bees respond to either stressor—winter bees are generally more susceptible to both stressors and have lower immune gene expression than summer bees (113, 115). Similarly, colonies in agricultural landscapes (where they are likely to encounter pesticides) with access to more floral resources in the summer, have higher Vg expression and greater winter survival rates than colonies with access to poorer food resources in the summer (122). This highlights the need to study the complex interactions of multiple stressors and the potential benefits of increased access to high quality nutrition or more targeted nutritional supplement strategies based on colony states or environmental pressures.

Bees are known to adjust gross foraging preferences based on the demands of the colony and social immunity requirements (66, 123, 124). Therefore, the second aspect of this study was to determine if individual honey bee foragers selected food items based primarily on their individual or colony's viral infections. Our data indicate that despite the tensions between individual benefits and social immunity (47, 63), the impact of colony viral infections appeared to outweigh individual infection in both the SEMs and regression comparisons. This may be due, in part, to the differences in infection levels between foragers and colonies—forager viral infections except for BQCV had lower levels and were less variable than those found in the colony-representative nurse bee samples, meaning that nurse bee infections may be more acute or there is survival bias by the time bees become foragers. So, while individual forager infection might reduce foraging efficiency (125), collective foraging decisions are directed by social immunity needs.

The importance of social immunity in collective foraging decisions is also compatible with prior observations of collective foraging of forager groups being dispatched to several food patches rather than individuals comparing food quality across patches (15, 24, 126, 127). These foragers make individual decisions about a food patch that depend on nutritional content and availability of food items that then accumulate across foragers for a collective outcome (128). However, if collective foraging is impacted by colony viral infection, individual foragers must be able to detect shifts in colony viral loads. This may occur via infection or parasitism-induced alterations to the bee cuticular hydrocarbon profile, signaling to nestmates (like foragers) changes in health status that may result in nestmate (forager) behavioral changes (129–132). For instance, in colonies with high Varroa mite loads, foragers changed dancing locations to be closer to the frame periphery relative to colonies with lower mite levels (133). Like the overall colony need for food resources, changes in chemical cues based on infection and/or parasitism may indicate to foragers that macronutrient preferences need to be adjusted.

Since several of the viruses found to influence bee foraging are vectored by Varroa mites (134, 135), we expected to potentially observe foraging differences based on the interaction of mite resistant bee stocks and viruses (52, 136, 137). Preliminary data (Figures 3, 4; Supplementary Figures 1, 2, 5–9) indicate that two mite-resistant stocks may be more responsive to viral infection in terms of P:L ratios and nectar load weights than the susceptible Italian stock. However, the mite resistant stocks did not appear to fully align with each other in how and for which viruses they changed their foraging. Given that the mechanism of mite resistance differs between Pol-Line and Russian stocks (136, 137), we might expect that the two genotypes also differ in their immune response and related nutritional requirements when infected with mite-vectored viruses (138). While preliminary due to low colony sample sizes per stock, these observations indicate that overall mite resistance may influence the host's likelihood to forage in response to pathogens (86, 103, 139, 140) but that the genetic separation of mite-resistant genotypes may further differentiate responses to viruses (52, 141). The incorporation of host genetics into future investigations of host-pathogen interactions is critical for understanding more global effects of pathogen challenges, particularly given the potential of Varroa mites becoming resistant to current miticides and in relation to mite-vectored viruses (142–146).