Diet restriction field experiment
Construction and build-up of beehives
The experiment consisted of 18 full-size hives with 8 frames each: 2 drawn combs; 5 foundation combs (templates for the construction of the drawn combs) and a honeycomb (frame full of honey). Each hive was made up of one queen bee (genetically different from each other) and 2 kg of worker bees. A raiser and a mat were placed under the lid. Bees and equipment for this experiment were provided by the apiary of The University of Western Australia in Perth.
After the construction of the hives, we immediately placed them in a Red Gum (Corymbia calophylla) forest at Gidgegannup, WA (31°47'59.6"S 116°10'13.8"E), allowing free foraging for 4 weeks. Red Gum eucalypt was chosen since it is an excellent source of nectar and pollen for honey bees40. We placed 10 colonies on January 31st, 2019 and a further 8 colonies on February 5th, 2019.
After 4 weeks of freely foraging, the hives were shifted to a burnt unproductive site at Lancelin WA (30°47'11.8" S, 115°16'48.7"E), to stop foraging activity and to control the feeding. We randomly assigned the colonies to three different nutritional treatments, with six replicate colonies each: 1) Starving - no food provided, 2) Sugar - colonies fed each week with 1 L of sugar syrup 150% (w/v), and 3) Sugar + pollen - colonies fed each week with a pollen patty of 250 g of Red Gum pollen mixed with 50 mL of 150% (w/v) sugar syrup and 1 L of sugar syrup 150% (w/v). The pollen patties were placed under the mat of each hive, and the sugar syrup was placed on a 1.10 L bucket container under the lid of the hives and over the mat. Each week 27 newly hatched worker bees were collected from each hive, and immediately placed on dry ice and stored at -80°C for biochemical analyses. As health indicators of the hives, we monitored the weight, temperature, humidity, and quantity of brood.
Data and sampling collection lasted for 4 weeks at the unproductive site when we detected that the health condition of the hives from the starving diet treatment was deficient. The hives were moved to a honey flow in a White Gum (Eucalyptus wandoo) forest (excellent source for honey production according to Coleman 196240) in Calingiri, WA (31°09'05.0" S, 116°18'05.1" E) to allow the free foraging of all the colonies. The recovery or maintenance of health was monitored with the same variables for a further 4 more weeks, and the sampling collection continued. From week 8 to 25, we shifted the hives into areas of good nectar sources, and we tracked their productivity by weighing the hives once a month (at weeks 12, 16, 21 and 25), and assessing the capped brood quantity once at week 21. Biochemical analyses were carried out for the bees collected from the diet treatment period and the recovery period (weeks 1 to 8, Fig. 1).
Nutritional disorders in the colony affect first-generation, but also the subsequent ones 41. Nurse bees feed the worker larvae from the first day of hatching from the egg until before the cell is sealed. This is from day 4 to 9 of the worker bees' development in the cell, which hatch 12 days after being capped 12. Because newly hatches bees were used for the study, they might represent the nutritional condition of the nurses that fed them on their larval stage, showing the nutritional condition of the bees 12 days ahead of the hive health indicators. Hence, the lipidomic analysis baseline was 12 days after week 0 of the experiment, which is between weeks 1 and 2 (Fig. 1). Subsequently, we used week 1 as the baseline for the lipidomic analyses, which represented bees that were still fed at the end of the build-up period, but not fed during the diet treatment period yet (Fig. 1).
Hive health monitoring
Weight, temperature, humidity, and quantity of the brood were measured each week. The hives were weighed using a digital weighing floor scale and hive weight was used as a measure of general hive strength and honey storage. We measured temperature and humidity by placing a digital monitor (remote relative humidity/temperature monitor-800027, Sper Scientific) in each hive. To measure the quantity of capped brood in each hive we took photos of both sides of each frame (Canon PowerShot A3300 IS). The photos were processed in the software “Beestly” (https://cyency.com/products/beestly/index.html) which allows the easy identification and quantification of capped brood cells. Finally, at the end of each sampling day the hives were fed according to the diet treatments.
Samples preparation for mass spectrometry analysis
The bees’ fat bodies from weeks 1 to 8 were dissected by removing the gut and the rest of the organs enclosed in the abdomen and conserving only the abdominal segment with the fat body attached at the dorsal region.
Each sample consisted of 20 abdomens (approx. 180 mg). The samples were weighed and processed for analysis in two separate batches that were created by randomization of hives. A solvent of 80% aqueous methanol (MeOH) was added at a concentration of 4:1 aqueous MeOH/bee material. The solvent contained a labelled internal standard at approximately 0.8 ppm of sphingosine-d9, tryptophan-d5, and taurodeoxycholic acid-d4. Samples were ground in 2 mL tubes containing ceramic beads (Precellys CKMix for tissue homogenizing), using a tissue homogenizer (Precellys evolution, Bertin instruments) at 6,000 g, 2x20 sec, at 0°C, with 15 sec of rest between cycles.
Untargeted lipidomic analysis by liquid chromatography high resolution mass spectrometry (LC-HR-MS)
For lipidomic analysis 100% chloroform (CHCl3) was added to the extracts, reaching a final concentration of 1:1 aqueous MeOH/CHCl3. Samples were vortexed for 10 sec and subsequently incubated for 30 min at room temperature and 750 rpm in a Thermomixer. The samples were dried under nitrogen and resuspended in 100% isopropanol (IPA) according to weight, obtaining a concentration of 4:1 IPA/bee material. The samples were centrifuged at 9,000 g, 4°C, for 10 min and the supernatant was split into 2 separate tubes, where 200 µL were destined for lipids analysis and 200 µL for fatty acid analysis. All the samples were dried under nitrogen and the samples destined for lipid analysis using liquid chromatography high resolution mass spectrometry (LC-HR-MS) were resuspended in 150 µL of IPA containing an internal standard of 1 ppm of deuterated cholesterol. Samples were vortexed and the supernatant isolated and transferred to a glass vial ready for analysis. Quality control (QC) samples were prepared by pooling 90 µL of each sample from batch 1.
Samples were analysed by high-pressure liquid chromatography (Dionex UltiMate 3000 RS) coupled to an Orbitrap Q-Exactive mass spectrometer (Thermo Fisher Scientific) with a heated electrospray ionisation probe (HESI) as the interface. Separation was achieved on a reversed phase Supelco Analytical Titan C18 column (2.1 × 75 mm, 1.9 µm particle size). The mobile phase was (Solvent A) 60:40 acetonitrile:water containing 10 mM ammonium formate and 0.1% formic acid and (Solvent B) 85.5:9.5:5 2-propanol (IPA):acetonitrile:water also containing 10 mM ammonium formate + 0.1% formic acid. The elution gradient was as follows: isocratic with 20% solvent B for 0.5 min, followed by an increase to 100% solvent B (0.5-8.5 min) and maintained at 100% B for 1 min. Then the system reverted to initial conditions (20% B) over 2 minutes and was equilibrated for 2.5 min before the next injection. The flow rate was 0.4 mL/min, the injection volume was 0.1 µL and column oven temperature was 55°C. Full scans with data-dependent tandem mass spectrometry were acquired on the Orbitrap mass analyzer. Full scans were acquired at a resolution of 70,000 at mass-to-charge ratio (m/z) 200 over the m/z range 150–2000 with the ESI conditions as follows: capillary temperature = 350°C, sheath gas = 48 (arbitrary units), auxiliary gas = 15 (arbitrary units), ion spray voltage = + 3.2 kV S-lens 60%. Tandem mass spectrometry analyses were performed at a resolution of 17,500 at m/z 200 on each sample with the collisional dissociation energy set at 40 eV. Data acquisition was carried out using Xcalibur software (Thermo Fisher Scientific).
Lipids data was analysed in MS-DIAL (version 4.18), where peak detection, identification and alignment were done using the LipidBlast adjusted Fiehn lipidomics dataset for the peak identification42. Data normalization and exportation was implemented as part of the MS-Dial software as recommended by Perez de Souza et al.43. To decrease the possibility of duplicated features and to increase the confidence of the identified compounds in the lipid analyses, only analytes that generated a MS fragmentation pattern were considered. We detected repeated diacylglycerols and triacylglycerols identities by sodium and ammonium adducts presenting the same results independently of the adduct. Both adducts display a structurally informative fragmentation pattern44. Hence, to avoid duplications, we conserved the ammonium adducts, which consistently presented a higher area abundance.
Fatty acids methyl ester (FAME) analysis
Previously dried samples for fatty acid analysis were resuspended in 305 µL of extraction buffer [MeOH containing 2% H2SO4 (v/v)]. The tubes were sealed with microtube cap locks, then mixed (using a vortex) for 10 sec and incubated for 2h at 80°C and 750 rpm in the Thermomixer. Samples were cooled to room temperature and 300 µL of NaCl [at 0.9%] and 300 µL of hexane were added. Samples were vortexed for 10 s and centrifuged at 15,700 g, for 3 min at room temperature. Next, 100 µL of the upper layer were transferred into a glass vial and used for FAME analysis by GC-MS. A QC sample was pooled from 90 µL of each sample from batch 1. The analysis was performed on a Gas Chromatograph (Agilent 6890N), equipped with a Gas Mass Selective Detector, 5975 series (Agilent Technologies). Separation was completed on a Agilent VF-5MS column (30m x 0.25mm, 0.25µm film thickness) + 10 m EZ-Guard column. The sample was injected at an initial temperature of 280°C and held for 5 min with a solvent delay of 5 min. Helium was used as a carrier gas with a flow rate of 1 mL/min. The purge flow was set to start at 2 min at 20 mL/min. The oven temperature was programmed as follows: initial temperature 70°C, hold for 2 min, then increase to 350°C at a rate of 7°C/min, and a post-run hold time of 5 min at 280°C. The transfer line temperature was kept at 250°C, and the ion source temperature was 230°C. The detector operated in scan mode from 40 to 600 Da. The MS Single Quad operated at 150°C and the collisional dissociation energy was set to 70 eV.
Fatty acids’ data was first analysed with Agilent MassHunter Qualitative Analysis Navigator (version B.08.00) to discard background peaks from blank samples. The fatty acids were identified by deconvoluting peaks with the Automatic Mass Spectral Deconvolution Identification System (AMDIS) (version 2.64, 2006) and accessing the Mass Spectral Search Program (NIST/EPA/NIH/Mass Spectral Library, Version 2.0.d, 2005) to identify each compound. After the main compound features were identified and named, we used the program Agilent MassHunter Quantitative Analysis (version B.08.00) to align and quantify peak intensities.
Desaturase genes primers, RNA isolation and Quantitative Real-Time PCR
Desaturase genes of honey bee were found in the National Center for Biotechnology Information (NCBI). We used primers previously reported for five acyl-CoA Delta (11) desaturases from Falcon et al.33 and one more designed from Vernier et al.45 Using the software SnapGene 5.2 (www.snapgene.com) we designed primers for three more acyl-CoA Delta (11) desaturase or desaturase-like, and one sphingolipid Delta (4)-desaturase (Supplementary Table 1). AmACT was used as housekeeping gene since has been reported as a stable housekeeping gene for honey bees46. Primers sequences used for each gene and specific information is listed on Supplementary data Table 1.
Bees from weeks one, three and four were dissected as previous and ground in liquid nitrogen using mortar and pistil (using 7 bees per sample and samples consisting of 3 to 4 hives per treatment). RNA isolation was done with 15 mg of ground tissue from each sample using SV Total RNA Isolation System (Promega, Madison, WI), according to the manufacturer’s specifications. 1.2 µg of RNA were used for cDNA synthesis using iScript cDNA synthesis kit (Bio-rad, 1708890). cDNA was further diluted and used for quantitative real-time PCR using QuantiNova SYBR green PCR kit (Qiagen, 208056). Values were normalized against the housekeeping gene and the relative abundance of each gene by sample was calculated by the comparative Ct method47.
Statistics and reproducibility
For hive health data, we fitted generalized linear models with quasi-poisson distribution for capped brood cell counts and gamma distribution for the rest predicted variables (hive weight, temperature, and humidity). We considered treatment (6 hives per diet treatment), time (sampling week) and the second order interaction of both as response variables. Akaike information criterion (AIC) value was used to find better-fitted model and distribution. To determine significant differences between specific time points we did post-hoc comparisons, adjusting the P-value with the false discovery rate (FDR).
Precision post-analysis tests were done to clean the MS data for both FAME and lipidomic results, determining the relative standard deviation (RSD) and the dispersion ratio (D-ratio) for untargeted analysis, discarding the compounds that were over the acceptance criterion (RSD > 20% and D-ratio > 50%), (for more details refer to Broadhurst et al.48. Lipidomics data was analysed with generalized linear mixed effect models with gamma distribution, where the logarithm of each compound peak area was the predicted variable. The response variables were treatment, time, and their interaction. The batch of laboratory analyses (two levels) was included as a random effect variable since we detected an effect in the raw data of the samples (Supplementary Fig. 1). Significance was tested with Least Squares Mean (LS mean) comparisons using FDR for adjusting the P-value. Fatty acid results did not present a batch effect; however, the machine presented a performance decrease during the running of the samples. Subsequently, we performed a QC based random forest signal correction (QC-RFSC) using statTarget package in R language49,50. The values were normalized to time 1 by treatment and Mann-Whitney-Wilcoxon test was used to evaluate significant differences. In a closer examination of the lipidomics’ data the software Perseus 220.127.116.11 was used to generate hierarchical clusters with the estimates of the models considering only the LS mean significant compounds in each treatment (P < 0.05).
The relative expression of desaturase genes was analysed using generalized linear models with gamma distribution or generalized linear mixed effect models when considering hive as random effect resulted in a better-fitted model (according to the AIC value). We ran LS mean comparisons using FDR to test significance. In every analysis assumptions of homogeneity of variance and dispersion were checked with Levene test and residuals distribution plots. Apart from Perseus software for the clustering analysis, the main statistical analyses and plots were performed in R50.