Transcriptional responses to fluctuating thermal regimes underpinning differences in survival in the solitary bee Megachile rotundata

Insects have evolved the ability to cope with harsh, long-term environmental conditions by using an array of physiological and behavioral responses. One of the most common defenses against unfavorable environmental conditions is the physiological state known as diapause, a life-history strategy characterized by diminished metabolic activity and a hiatus in development. While diapause offers a level of protection against seasonal stress, such as cold winter temperatures, exposure can still take a physiological toll when temperatures are abnormally low or when winter lasts for a long time (Hayward et al., 2014; Renault et al., 2004; Teets and Denlinger, 2013).

Exposure to low temperatures, also known as chilling, can have a wide range of physiological effects (Teets and Denlinger, 2013). The physiological effects associated with long-term chilling, such as exposure to winter conditions, are known as indirect chill injuries and have been associated with gradual failure of homeostatic processes (Lee, 2010). The downstream consequences of this phenomenon are likely to be complex, but disruption of ion homeostasis and metabolic imbalance are probable outcomes (Kostal, 2004; Koštál et al., 2006). Additionally, chill injury has been associated with the occurrence of oxidative stress (Lalouette et al., 2011; Rojas and Leopold, 1996).

While the accumulation of chill injuries can have a significant effect on the well being of the insect, periodically increasing temperatures during long-term chilling can allow for an increase in survival (Colinet et al., 2006; Coulson and Bale, 1996; Rinehart et al., 2011, 2013). These periodic increases in temperature, commonly known as fluctuating thermal regimes (FTRs), have been implicated in the repair of, or protection against, chill injury during exposure to low temperatures and have been observed across insect life stages and taxa (Koštál et al., 2007; Renault et al., 2004). Despite the common physiological mechanism associated with chill injury and its repair, few studies have measured the chill injury response at the level of gene expression (Hayward et al., 2014).

The alfalfa leafcutting bee Megachile rotundata is widely used in alfalfa seed production agroecosystems and provides an excellent model for assessing changes in gene expression during chilling. Adult M. rotundata generally emerge during early summer. Soon after emergence, the females mate once and begin provisioning brood cells (Pitts-Singer and Cane, 2011). Larvae develop within their brood cells until the fifth instar. At this point, they defecate and enter diapause as prepupae and remain dormant through the winter (Michener, 2007). Development resumes as temperatures warm in late spring or early summer; bees pupate and later emerge as adults (Pitts-Singer and Cane, 2011). Diapause in M. rotundata often terminates before harsh environmental conditions end and is followed by a physiological state known as post-diapause quiescence (Yocum et al., 2005, 2006). Although quiescent individuals are responsive to warmer temperatures, they maintain the cold tolerance that is characteristic of diapause (Hayward et al., 2005). Quiescent individuals may be roused to continue development, but will stay quiescent until more favorable conditions arise (Koštál, 2006).

The harsh winter temperatures experienced by M. rotundata during diapause and post-diapause quiescence have the potential to cause chill injury in both natural populations and commercial-use bees. In managed populations, diapausing individuals are typically overwintered at a constant temperature of 4–6°C (Pitts-Singer and Cane, 2011). Individuals reared under constant, low temperatures over long periods of time likely experience an accumulation of chill injuries, eventually leading to an increase in mortality (Kostal, 2004; Koštál et al., 2007). When exposed to daily fluctuations to a warmer temperature (20°C) during chilling, M. rotundata shows a dramatic decrease in mortality, and this decrease is maintained for many months (Rinehart et al., 2013). After extended exposure to fluctuating temperatures, emerging adults have no significant differences in quality compared with control (Bennett et al., 2013; Rinehart et al., 2013). This result suggests that fluctuating temperatures provide a protective effect during chilling, alleviating sub-lethal effects associated with long-term low constant temperature exposure. We hypothesized that: (1) long-term exposure to low temperatures causes chill injury and that subsequent exposure to fluctuating temperatures provides a protective effect by repairing and/or diminishing the effect of chill injury; and that (2) this protective effect is a result of a decrease in oxidative stress and a return to more favorable ion concentrations and metabolic rates.

In this study, we screened for the transcripts responsible for the differences in survival between individuals reared under either constant or fluctuating temperatures. RNA samples extracted from post-diapause quiescent prepupae were harvested from individuals reared under both temperature treatments before and after mortality rates began to diverge. Treatment- and time-specific transcriptome profiles were assessed using RNA-seq on the Illumina platform. We predicted: (1) that a larger quantity of transcripts would be differentially expressed in late-sampled individuals, after mortality diverged; and (2) that transcripts functioning in a stress response (i.e. oxidative stress pathways) and those functioning in diminishing other sub-lethal effects associated with chill injury would be up-regulated in individuals exposed to fluctuating temperatures.

Individuals were harvested for RNA-seq at two time points: (1) after 12 months in storage, before mortality rates diverge between temperature treatments and (2) 2 months later, after mortality had significantly decreased under constant temperatures (Fig. 1A; modified from Rinehart et al., 2013). For simplicity, individuals selected before mortality diverged will be referred to as ‘early’ and those selected after mortality diverged will be referred to as ‘late’ (Fig. 1B). Exposure to constant temperatures will be referred to as STR (static thermal regime), and exposure to fluctuating temperatures as FTR (fluctuating thermal regime) in all subsequent figures and tables. Twelve individual RNA-seq libraries (accession number: SRP047335) consisting of an early and late time point for each treatment with three biological replicates for each were generated, averaging 62,774,914 pair-end reads per library (Table 1). Both the percentage of raw reads mapped and the percentages of reads properly paired with its mate were used as initial quality metrics of the assembly. On average, approximately 80% of the raw reads generated mapped to the M. rotundata genome (accession number: PRJNA66515); of these, 88%, on average, could be properly paired with its mate (Table 1).

Four pair-wise comparisons were deemed biologically relevant because they either directly compare expression profiles of one treatment at two time points or compare the differences between two treatments in the same time point (Fig. 1B). Differential expression analysis was conducted using the Tuxedo protocol (Trapnell et al., 2012). Transcripts with an absolute value greater than or equal to log₂2 (α=0.05) were selected for downstream analysis. With this conservative significance threshold, 287 differentially expressed genes were identified among the four biologically relevant comparisons. We confirmed 79 comparisons (both differentially expressed and not) from the expression analysis using quantitative real-time PCR (qPCR; supplementary material Table S1). RNA samples assessed for qPCR validation originated from different individuals than the ones subjected to RNA-seq, providing an independent confirmation of expression differences.

In differential expression analysis, multiple transcripts may represent a single gene locus. The Tuxedo protocol averages expression values for each transcript from a gene to give the fold change value for the gene (i.e. each gene is an average of all transcripts at that locus). To include the identities of each individual transcript in downstream analyses (319 among all four biologically relevant comparisons), each transcript was treated as a unique entity for analysis and is included in all of the following iterations and numerical representations of the data.

The expression profiles show a clear asymmetry in the distribution of differentially expressed transcripts among comparisons (Fig. 2). The largest number of differentially expressed transcripts came from the comparison of individuals exposed to either constant or fluctuating temperatures before mortality diverged (Fig. 2A). We had predicted that the greatest number of differentially expressed transcripts would occur after mortality diverged because of the marked difference in mortality at that time. However, only 63 transcripts were differentially expressed at the later time point. The larger number of differentially expressed transcripts in the early comparison (217) indicates that expression profiles are changing months before observable changes in bee quality.

Few differentially expressed transcripts were identified for the within-treatment comparisons (i.e. early STR versus late STR) for both constant and fluctuating temperatures (12 and 22 transcripts, respectively; Fig. 2B). The low number of differentially expressed transcripts indicates that the profiles of samples sequenced in each treatment over the two dates are relatively static between the time points sampled. This stasis could be due to either the emergence of large-scale differences before the time-points sampled, or relatively few transcripts regulating the physiological response to the temperature treatments.

Venn diagrams were used to identify differentially expressed transcripts among the four biologically relevant comparisons that were unique either temperature treatments or time points sampled (Fig. 3). Differentially expressed transcripts in early STR versus early FTR, the comparison with highest distribution of differentially expressed transcripts (67.8%), shared 23 of its 215 differentially expressed transcripts with the late STR versus late FTR comparison (Fig. 3A). These 23 transcripts (22 up-regulated, 1 down-regulated) are specific to individuals reared under fluctuating temperatures before and after mortality diverged (Fig. 3A). Because the direction of regulation is relative to the second treatment in each comparison, an FTR-specific up-regulated transcript is also a down-regulated STR-specific transcript.

In additional to identifying transcripts differentially expressed in animals exposed to the FTR treatment across time, the Venn diagrams identified up- and down-regulated transcripts that are unique to each comparison. For example, 43 down-regulated transcripts were unique to the early STR versus early FTR comparison (Fig. 1B) and 36 transcripts were up-regulated in late STR versus late FTR that were not shared with any other comparison (Fig. 1C). Distinctive transcriptional differences exist between treatments within time points.

A gene ontology (GO) analysis, conducted using Blast2GO (Conesa et al., 2005), revealed an abundance of transcripts functioning in oxidative stress, various metabolic functions, ion homeostasis, neurogenesis immune response and several functioning in growth and development (Table 2; supplementary material Table S2). These transcripts are up-regulated before mortality diverged between temperature treatments. An ontological analysis of the 63 differentially expressed transcripts after mortality diverged (late STR versus late FTR) shows an up-regulation of transcripts functioning in similar GO classes to the earlier comparison (early STR versus early FTR). Twenty-three of the 63 differentially expressed transcripts (36.5%) in late-sampled individuals were also differentially expressed in early-sampled individuals (supplementary material Table S3).

Chilling can be deadly to insects, but periodic warm pulses during chilling can partially alleviate the damage of chill injury. The gene expression changes associated with the protective effects of these warm pulses are unknown. Assessing gene expression changes during chilling will provide validation of the already-established physiological responses and may elucidate additional mechanisms that have not been previously identified. Here, we present the first transcriptome-level analysis of increased survival under fluctuating temperatures during chilling. Megachile rotundata were exposed to either a constant, low temperature (6°C) or the same low temperature with a daily fluctuation in temperature to 20°C. When reared under these conditions, individuals exposed to daily temperature fluctuations show a marked increase in survival, suggesting some type of physiological benefit to periodic deviations away from constant, low temperature exposure. We hypothesized that prolonged exposure to a constant temperature leads to indirect chill injury in M. rotundata, and further, that the increase in survival in bees that were exposed to fluctuations in temperature during post-diapause quiescence would be the result of a decrease in the physiological stressors associated with chill injury.

First, we predicted that more transcripts would be differentially expressed in late-sampled individuals, after mortality diverged. However, differential expression analysis reveals an opposite trend; the greatest quantity of differentially expressed transcripts occurred in early-sampled individuals, before any significant differences in mortality were observed. This result suggests that sub-lethal effects accumulate prior to significant mortality in the constant-temperature-reared population. The relatively small within-treatment changes in expression profiles between the two time points, in both FTR and STR, also support the hypothesis that the physiological mechanisms responsible for counteracting the accumulation of sub-lethal effects are set into motion earlier in extended post-diapause quiescence. Additional RNA-seq sampling or extensive quantitative real-time PCR (qPCR), over a longer developmental interval, starting further upstream from this analysis, would be necessary to validate this hypothesis.

Our second prediction was that transcripts functioning in a stress response and those functioning in diminishing other sub-lethal effects associated with chill injury would be up-regulated as a result of exposure to fluctuating temperatures. We discovered transcripts functioning in ion and metabolic imbalance, which have been hypothesized to result in sub-lethal effects, and several with antioxidant properties; results that are in agreement with our initial predictions. Interestingly, several other classes of transcripts were up-regulated in individuals reared under constant temperatures, suggesting additional mechanisms for coping with stresses related to long-term low temperature exposure.

Chill injury and subsequent mortality have been associated with aberrations in ion homeostasis and metabolic imbalance in several insect species (Koštál, 2006; Koštál et al., 2007; Lalouette et al., 2011; Macmillan et al., 2012). Disturbances in ion homeostasis could help to explain some of the large-scale physiological effects seen in chill-injured insects such as abnormal muscle contraction (Yocum et al., 1994) and atypical neuronal function (Hosler et al., 2000).

Our GO analysis revealed that an abundance of transcripts involved in counteracting disruptions of ion homeostasis and energy metabolism (including lipid metabolism) (Tables 2; supplementary material Table S2), both of which have been implicated as molecular components of chilling injury physiology (Hayward et al., 2014). The return of more favorable concentrations of ions within the organisms exposed to fluctuations in temperature may provide protection against failure of neuromuscular coordination in chill-injured individuals. These results suggest that the up-regulation of these transcripts under exposure to fluctuating temperatures may allow for protection against, or repair of, chill injury and, thus, lead to an increase in survival.

Traditionally, the most obvious symptom of chill injury has been loss of coordination such as defects in crawling behaviors (Hazell and Bale, 2011). This physiological impact suggests a neurological consequence of chilling. In individuals exposed to fluctuating temperatures during quiescence, a total of 14 transcripts (Table 2) functioning in neurological patterning and development were up-regulated in early sampled individuals and either were maintained or up-regulated further after mortality rates began to diverge. The maintenance of these expression profiles through the two time points suggests, again, that these transcripts are playing a role over time, even when mortality has already diverged.

The up-regulation of these transcripts in fluctuating temperature-reared individuals suggest one of two mechanisms: either these transcripts are functioning in the repair of damage caused by chill injury or they are functioning in normal neurological development and exposure to constant chilling causes a phenocopy defect, impeding normal developmental processes during this life stage. Phenocopy defects were originally used to describe morphological deformities in Drosophila exposed to stressful environmental conditions during development (Goldschmidt, 1935). Phenocopy defects are common under high temperature stress (Mitchell and Lipps, 1978; Mitchell et al., 1979). The mechanisms responsible for the physiology associated with high- and low-temperature stress show some overlap (Yocum et al., 1994), indicating that decreased survival in individuals injured by chilling may be the result of phenocopy defects. A phenocopy-defect-centric hypothesis would suggest that the suspension of normal developmental processes might be a factor to increased mortality between the two treatments. The blocking of development during the STR treatment may lead to additional stresses not previously associated with chill injury. The initiation of transcription associated with developmental processes in the FTR treatment may mitigate this phenocopy defect.

Current estimates of the developmental threshold in M. rotundata range from 15.7–19.0°C (Kemp and Bosch, 2000; Whitfield and Richards, 1992). While these estimates, even at the upper bound, are still below the 20°C pulse experienced during the FTR protocol, although only by as little as one degree, it has been shown that it is not possible for M. rotundata to complete development at 18°C and show greatly delayed development and increased mortality at 22°C (Kemp and Bosch, 2000). While 20°C is indeed above M. rotundata's developmental threshold, the individuals harvested for the RNA-seq assay had not yet moved into the next developmental stage (pupation) and any developmental processes occurring during the exposure to the FTR treatment are not sufficient to lead to gross morphological changes in the animal. Although large-scale morphological changes were not witnessed in M. rotundata as a result of the FTR treatment, we recognize that smaller-scale developmental processes may still be occurring.

In addition to those functioning in neurogenesis, the up-regulation of transcripts functioning in other aspects of development, growth, or cell differentiation (Table 2) provide additional evidence that fluctuating temperatures may allow for normal development in preparation for pupation and adult development. Conventionally, diapause and post-diapause quiescence have been thought of as developmentally inert, but the idea that these life stages are more developmentally dynamic, have started to gain momentum (Koštál, 2006); our data provide support for this paradigmatic shift. The up-regulation of these developmental transcripts also argues for a role in the prevention of phenocopy defects caused by sub-optimal environmental conditions.

The low temperatures experienced during chilling may have the capacity to reduce the effectiveness of antioxidant enzymes (Joanisse and Storey, 1996), and thus contribute to chill injury through oxidative stress, causing damage to DNA, protein and lipid molecules (Monaghan et al., 2009). In adult Alphitobius diaperinus, exposure to fluctuating temperatures during chilling has been associated with decreased levels of reactive oxygen species (ROS) when compared with individuals reared under constant low temperatures (Lalouette et al., 2011).

In this study, we found a number of transcripts functioning in oxidation stress response, including Peroxidase, and Glutathione S-transferase, and several from the Cytochrome P450 (CPY) superfamily (Table 2). These transcripts were up-regulated in early-sampled individuals reared under fluctuating temperatures and maintained in late-sampled individuals. Exposure to environmental stressors (i.e. low temperature stress), has been associated with malfunctioning in mitochondrial respiration and the production of peroxide (Prasad et al., 1994). The up-regulation of Peroxidase witnessed in M. rotundata may serve to diminish the resulting oxidative stress associated with the aforementioned metabolic malfunction. During exposure to low temperatures, the function of antioxidant enzymes, such as those belonging to the CYP superfamily, can be diminished; leading to an increase in ROS present in the cell (Baek and Skinner, 2003). The increased expression of CPY transcripts suggests that exposure to fluctuating temperatures may help to counteract this decreased efficiency.

Furthermore, other transcripts with antioxidant properties and associations with insect survival were up-regulated when M. rotundata were exposed to fluctuating temperatures. Vitellogenin (Vg), which has been implicated in caste-specific differences in longevity in honeybees and has antioxidant properties (Corona et al., 2007) and Apolipoprotein D-like, another transcript associated with oxidative stress and survival in Drosophila (Walker et al., 2006), were both up-regulated. The up-regulation of these transcripts again supports our hypothesis that a reduction of oxidative stress plays a role in increased survival caused by fluctuating temperatures during chilling.

Increased cellular membrane damage via an accumulation of ROS to stressful levels may be associated with damage to mitochondrial membranes and may exacerbate ion and metabolic imbalances that have previously been correlated with chilling injury in other insects. Our results suggest that individuals reared under constant, low temperatures may be experiencing elevated levels of ROS and that a periodic increase in temperature may help alleviate the resulting oxidative stress via up-regulation of antioxidants.

Several other functional classes of transcripts that have not previously been linked to chill injury, including those functioning in growth and development and immune response, were up-regulated in early-sampled individuals exposed to fluctuating temperatures and were maintained in late-sampled individuals. These additional functional classes provide mechanisms leading to increased survival that have not been previously associated with a chill injury response.

Environmental stressors, such as temperature, have been associated with a decrease in immune function across taxa through trade-offs in life history strategies (Lochmiller and Deerenberg, 2000). For example, temperature stress in M. rotundata decreases the immune response to chalkbrood disease caused by Ascosphaera aggregata (Xu and James, 2012). Our study revealed a similar suite of transcripts expressed at increased levels when individuals were exposed to fluctuating temperatures, suggesting that chilling may decrease immune function and may be associated with a decrease in survival. Transcripts such as Serine protease snake-like, Protein Toll and Melanization-related protein, which were differentially expressed in our study (Table 2), have been previously implicated in temperature-specific immune responses of M. rotundata (Xu and James, 2012). Additionally, increased levels of ROS have been associated with immune function because of free radial production during phagocytosis (Nappi and Vass, 1993), suggesting another cause for the up-regulation of the antioxidants highlighted above. Antioxidants have also been shown to play a role in limiting microbial growth (Broderick et al., 2009), again suggesting another link to increased expression of antioxidants.

Our results suggest that rearing individuals under long-term overwintering conditions can lead to an accumulation of chill injury, damage to neural structures and/or defects in development, and a decrease in immune system function. When individuals are exposed to fluctuating temperatures after an extended period of chilling, insect survival increases dramatically. We have developed a working hypothesis, based on evidence from both the literature and this study, to explain the increased survival associated with fluctuating temperature exposure during chilling in M. rotundata (Fig. 4). In this model, exposure to fluctuating temperatures activates pathways functioning in oxidative stress response, ion homeostasis, metabolism and neurogenesis, leading to increased survival. This model not only supports current chill injury hypotheses, such as oxidative stress, ion imbalance, and impaired metabolic function by providing evidence for transcriptional responses, but also adds additional mechanisms such as neurological development/repair and increased immune function. While it is clear that exposure to fluctuating temperatures provides an enormous benefit to the organism and the mechanisms behind it are beginning to be understood, it is still ambiguous whether these mechanisms act in a protective fashion, by slowing or stopping the accumulation of chill injury entirely, or by repairing chill injury caused by previous long-term exposure to low temperatures.

All leafcutting bees (Megachile rotundata Fabricius 1787) for this project, derived from the 2009 field season, were purchased from JWM Leafcutter, Inc. (Nampa, ID) as loose cell bees and were of Canadian origin. Bees were stored at 6°C under constant darkness upon arrival.

Cells containing bees were placed individually into the wells of 24-well culture plates, and housed in Percival model I-30BLL reach-in incubators. In April 2010, after 7 months of storage at constant temperatures, post-diapause quiescent prepupae were placed into a fluctuating temperature treatment or left at standard storage temperatures. The constant temperature regime consisted of a 6°C±0.5°C constant temperature with a 15 h:9 h (L:D) photoperiod. Bees reared under temperature fluctuations were exposed to 6°C with a daily warm pulse of 20°C. The warm pulse consisted of a 1 h ramp up to 20°C (0.23°C min⁻¹), a 1 h incubation at 20°C and a 1 h ramp down, back to 6°C. Peak temperature occurred during the photophase of the 15 h:9 h light cycle. At monthly intervals, three of the 24-well culture plates were removed from the temperature treatments and place at 29°C to initiate pupal development. Adult emergence rates were used to measure survival (as described in Rinehart et al., 2013).

Total RNA was collected from three postdiapause quiescent prepupae reared under either constant or fluctuating temperature treatments at an early and late timepoint in 2010 (12 bees total) using the Invitrogen TRIzol protocol (Carlsbad, CA, USA). To ensure that we were assessing response to low temperatures in both treatments, individuals for both treatments were harvested during the cold phase (6°C). RNA was stored as an ethanol precipitate at −80°C until needed. Samples were then dissolved in RNAase-free H₂O and shipped at 1–10 µg total RNA at a concentration of no less than 20 ng µl⁻¹ on dry ice overnight to University Georgia.

RNA-seq libraries from both constant and fluctuating temperature treatments at both early and late time points in 2010 were prepared using Illumina TruSeq mRNA standard protocol at the University of Georgia Genomics Facility (Athens, GA, USA) and sequencing was outsourced to the University of Missouri Columbia DNA Sequencing Core facility (Columbia, MO, USA). The RNA samples were sequenced on a HiSeq2000, running HiSeq Control Software (HCS) v1.4.8. Samples were sequenced on two lanes of a HiSeq Flowcell v1.5 and libraries were clustered on a cBot v1.4.36.0 using Illumina's Truseq PE Cluster Kit v2.0 and sequenced using a 200 cycle TruSeq SBS HS v2 kit. The clustered flowcell was sequenced for 206 cycles, broken down into three separate reads. The first read was 100 cycles in length, followed by a six-cycle index read. Following the index read, paired-end resynthesis was performed also using Truseq PE Cluster Kit v2.0, which was then followed by another 100 cycles. Image analysis and base calling were performed using the standard Illumina Pipeline consisting of real-time analysis (RTA) version v1.12.4.2 and Casava v1.8 using the default settings.

Raw sequence data generated from the Illumina HiSeq2000 (Accession: SRP047335) were quality checked using FastQC (Version 0.10.1; Babraham Bioinformatics). The raw reads were aligned to the M. rotundata reference genome (accession number: PRJNA66515) using TopHat (v2.0.5). Mapped reads were assembled, with the aid of a GFF (v3) annotation file generated using the genome annotation pipeline MAKER (Cantarel et al., 2008), using Cufflinks (v2.0.2). The differential expression profiles were analyzed using Cuffdiff (v2.0.2) via the iPlant collaborative discovery environment (Stanzione, 2011). A threshold value of log₂2, (α≤0.05) was used to determine significance in the differential expression analysis. The R package cummerbund (Trapnell et al., 2012) was used for downstream analysis and the generation of differential expression figures. The Venn diagram depicted in Fig. 3 was generated using the R package ‘venn diagram’ (Chen and Boutros, 2011) and the Java-based gene ontology enrichment and functional annotation program Blast2GO was used for GO and KEGG analysis (Conesa et al., 2005; Kanehisa et al., 2012).

RNA samples for qPCR were collected from post-diapause quiescent bees during RNA-seq library preparation in one-month intervals between May 2010 and 2011 in the fashion described above. RNA samples were diluted to a concentration of 0.5 μg μl⁻¹ and then subjected to treatment with DNAase I (Invitrogen, Carlsbad, CA, USA) followed by first-strand cDNA synthesis using Super Script III first strand synthesis system for RT-PCR (Invitrogen). Three controls were used for this step: a ‘no-template control’, consisting of all enzymatic components except the RNA template, a ‘negative RT control’, lacking reverse transcriptase and a ‘positive RT control’, using a control RNA template.

qPCR was conducted on a LightCycler 480 (Roche, Indianapolis, IN, USA). SYBR Green I Master Mix with ROX (FastStart Universal, Indianapolis, IN, USA) protocol and reagents were used for the qPCR reactions. Primer design for all targets and reference genes was conducted using the IDR program from Integrated DNA Technologies (Coralville, IA, USA). qBasePLUS (Biogazelle, Ghent, Belgium) was used for analysis of the qPCR data.

The authors would like to thank Marnie Larson of the USDA-ARS, Fargo, ND for her technical assistance. We also thank Gene Robinson, Karen Kapheim and the B10 consortium for access to the unpublished Megachile rotundata genome and for permission for early submission.