Queens and inseminations
We performed an experiment near Athens, Georgia consisting of four treatments in a 2x2 factorial arrangement with two levels of genetic Varroa resistance (semen from VSH-selected drones or non-selected wild type drones) and two levels of polyandry (high or low). There was one dedicated research apiary in which all treatments were represented. Each colony consisted of a 5-frame Langstroth nucleus hive and was started with one 2-lb (0.9-kg) package of worker bees from one supplier practicing a common V. destructor treatment regimen, thus normalizing incipient parasite levels across the experiment. Colonies were arranged in one large circle at the apiary to minimize effects of bee and parasite drift (Dynes et al. 2019a). The experiment was set up in Apr 2018 and ran until Nov 2018. The year before, we made up two sets of dedicated drone-source colonies for future inseminations, one set headed by wild type unselected queens and the other headed by queens from the VSH line.

In spring of 2018, we reared one wild type virgin queen and inseminated her with one wild type male. This single-drone inseminated queen was used to rear supersister virgin daughters who were subsequently inseminated and used in the experiment. This procedure minimizes maternal variation in experimental colonies (Harbo 1985). Supersister virgins were each randomly assigned one of four treatments in the 2x2 design: (1) VSH semen/low polyandry, (2) VSH semen/high polyandry, (3) wild type semen/low polyandry, or (4) wild type semen/high polyandry.

We had nine strong VSH drone-source colonies come out of winter, a number that dictated the multiples for our polyandry ranges. For our low polyandry treatment, we inseminated each supersister virgin with one male from each of the VSH drone source colonies (mo=9), and for our high polyandry treatment each queen was inseminated with 6 males from each VSH drone source colony (mo=54). Supersister virgins in the wild type treatments similarly received one drone from each of 9 wild type drone source colonies (mo=9) or 6 from each wild type source colony (mo=54).

Each virgin queen was emerged and housed in her 5-frame nucleus colony before and after the artificial insemination (AI) procedure, following methods of Cobey et al. (2013). Semen from 9 or 54 males, as per protocol, was collected into one common capillary tube. A different tube was used for each prescribed polyandry/genetic line; no contamination across treatments was possible. However, more than one queen was inseminated from any given tube. The volume of semen in each tube was increased by 15% with physiological saline (recipe 2.2.1 in Cobey et al. (2013)), thoroughly mixed in an Eppendorf tube to encourage a homogeneous drone representation per batch, re-drawn into the capillary tube, then used to inseminate 2-20 queens at an average dose of 4 μL mixed semen each, thus removing effects of semen volume (Niño et al. 2012).

Queens were narcotized three times with CO2 – once on the day of insemination and once on each of two successive days thereafter; CO2 narcosis stimulates egg laying (Mackensen and Roberts 1948). Experimental queens were housed in their home 5-frame nuclei and observed until all competent egg layers were identified. We waited six weeks thereafter to allow worker populations to turn over to progeny of the experimental queens.

Dependent variables
Colonies were sampled at each of four time points (Jul, Aug, Sep, Nov). We sampled colonies to determine brood production, brood survival, Varroa mite levels, and adult bee populations.

Brood production and worker bee population were derived by visually summing proportions of whole deep frames covered by brood or workers (Delaplane et al. 2013), converting frames of brood to cm2 by the observation that one deep Langstroth comb (both sides) = 1760 cm2, and converting frames of adult bees to bee populations with the regression model of Burgett and Burikam (1985). When necessary, we converted cm2 brood to cells of brood with the conversion of 3.8 cells per cm2 (Delaplane et al. 2013).

Brood survival was measured by removing a comb of open brood, laying over it a sheet of transparent acetate, using a felt-tip permanent marker to mark the location of 100 brood cells each containing a 1st or 2nd instar larva, then returning the comb to the hive; three days later the comb was retrieved, the same acetate laid on top of it, and the number of surviving brood cells recorded.

Relative Varroa mite numbers are derived by inserting sticky sampling sheets into bottom board hive inserts and recording the number of Varroa mites trapped after 24 hours (Dietemann et al. 2013).

Patriline determination
As a check on our success at creating two discrete classes of polyandry with AI, we were able to genotype workers from 18 colonies to determine effective realized paternity (me) in both polyandry classes.

A 50-worker sample was taken from each of 8 colonies in the mo=9 polyandry class and 10 colonies in the mo=54 polyandry class. DNA was extracted from the right hind leg for all honey bee specimens using a Qiagen DNeasy Blood & Tissue Kit. Hind legs were thoroughly pulverized with micro-scissors and digested with proteinase K in a 56°C water bath for 3 hours with sample agitation every hour. Extracted DNA was stored at -80°C until microsatellite amplification. Ten variable microsatellite loci were screened for all samples using two multiplexes: Plex 1 and Plex 2. Plex 1 included loci A107, A113, AP043, A024, and A006; Plex 2 included loci A28, A88, AP66, AP81, and B124 (Shaibi et al. 2008; Delaney et al. 2009). Both plexes were amplified for one cycle at 95°C for 7 min, 30 cycles of 95°C for 30 sec, 54°C for 30 sec, 72°C for 30 sec, and a final extension at 72°C for 60 min. A 10 μL final reaction volume containing 5 μL of PCR Master Mix (Promega, Madison, WI), 1.0-2.5 μL of fluorescent dye-labeled primer, 0.9 μL of nuclease-free water, and 2 μL of DNA extract per sample was analyzed at the University of Delaware Sequencing and Genotyping Center at the Delaware Biotechnology Institute with an Applied Biosystems 3130 XL Genetic Analyzer via capillary electrophoresis.

Microsatellite repeat sizes were scored using Geneious Prime® 2020.0.5 software (Biomatters Ltd). The number (N0) and frequency (pi) of full sibling patrilines was determined using raw microsatellite data in the software program COLONY 2.0.6.5 (Jones and Wang 2010). The observed patriline number (N0), proportion of each patriline within each sample (pi), and number of worker bees in each respective sample (n) were used to calculate me after (Nielsen et al. 2003; Tarpy et al. 2015):

m_e=((〖n-1)〗^2)/(∑_(i=1)^(N_0)▒〖p_i^2 (n+1)(n-2)+3-n〗)

Identifiers for columns of interest in filename "Delaplane polyandry asymptote data for Dryad:"

1. A, state, just GA a placekeeper
2. B, “nuc” short for “nucleus” is identifier for individual hive, one experimental unit. Nucs were sampled repeatedly at different time trials (column M)
3. E, genetic type
4. G, poly type
5. I, ERP or effective realized paternity. For 18 hives only (!) we genotyped workers “after the fact” to see what the “real” paternity was. You can see that AI is very inefficient. The poly=54 group especially was all over the map. Nevertheless, we thought this could be an informative covariate.
6. K, semen batch. The semen is drawn into little capillary tubes and used to inseminate queens. Because of genetic recombination, none of these mixtures are clones; hence, semen batch is another source of random variation. Numbers apply within state.
7. M, trial, months at which dependent variables were measured. With honey bees there is always a huge seasonal effect as the colony grows and declines. Note ahead that not all dependent variables were measure all months.
8. W, beepop, adult honey bee population, along with
9. X, sqrtbeepop, transformed via square root.
10. AA, cm2 honey bee brood, along with
11. AB, brood square root transformed
12. AD, percentage brood survival, along with
13. AC, survival arcsine transformed
14. AE, parasitic mite counts, along with
15. AG, mites square root transformed
16. AH, percentage of mites reproducing (possibly an indicator of VSH effect), along with
17. AK, reproduction arcsine transformed
18. AL, cm2 comb constructed, along with
19. AN, comb square root transformed

2-7 are structural elements. 8-19 are dependent variables of interest.

Relevant column identifiers for second dataset, "Delaplane colony survival data for Dryad:"

1. B, nuc
2. C, whether colony survived to last sampling episode
3. E, genetic type
4. H, poly type