Behavioral and genetic correlates of heterogeneity in learning performance in individual honeybees, Apis mellifera

Learning an olfactory discrimination task leads to heterogeneous results in honeybees with some bees performing very well and others at low rates. Here we investigated this behavioral heterogeneity and asked whether it was associated with particular gene expression patterns in the bee’s brain. Bees were individually conditioned using a sequential conditioning protocol involving several phases of olfactory learning and retention tests. A cumulative score was used to differentiate the tested bees into high and low performers. The rate of CS+ odor learning was found to correlate most strongly with a cumulative performance score extracted from all learning and retention tests. Microarray analysis of gene expression in the mushroom body area of the brains of these bees identified a number of differentially expressed genes between high and low performers. These genes are associated with diverse biological functions, such as neurotransmission, memory formation, cargo trafficking and development.


Introduction
Some of the pivotal components of animal cognition are the adaptability to new environmental conditions, capacity to maintain successful changes in behavior caused by new learning, and the use of memories for further behavioral control.In honeybees, a classical model of learning and memory studies, a large body of data on how these cognitive functions are essential for their foraging success and survival of their colonies has been documented [1,2].Evolutionary selection for the improvement of cognitive abilities work on the individual level and are based on genetic variations that are favorable for behavioral traits related to sensory processing, learning and memory formation [3,4].The search for such cognitive adaptations can be performed both on a population level and on an individual level.However, our current understanding of the ecological and evolutionary drivers of intraspecific variability in cognitive abilities and the consequences of this variability on the individuals and populations is restricted [5][6][7][8][9] due to limited research work done on these topics.Studies across different vertebrate species recognized consistent among-individual variability in cognitive capacities [10][11][12][13][14], albeit tracing invertebrate cognition through the lens of individual variability in behavior is a much less explored territory, in spite of the increasing evidence of incredible cognitive abilities of invertebrates [1,[15][16][17].
Selection for behavioral traits usually involve many animals.Studies on learning and memory in Drosophila have led to the discovery of multiple genes related to learning and memory [18]; however, the olfactory conditioning procedures usually engage many individuals at the same time, e.g. in the T-maze paradigm [19], which makes it impossible to relate individual performances to their corresponding gene expression patterns.Olfactory conditioning protocols in Drosophila used by others also failed to decipher the heterogeneity in the expression of learned behavior among the individuals.Instead, the results were interpreted as the expression of probabilistic behavior at the individual level [20,21].Honeybees, on the other hand, can be selected for fast learning and then selectively bred via the thelytokous development of Apis mellifera capensis [22]; a highly effective procedure that can be applied more frequently to understand the contribution of genetic variability to the variation observed in cognitive ability.Tait and colleagues [23] tested individual honeybees for associative acquisition, non-associative sensitization, and short-term memory, and constructed a general heritable factor similar to the factors derived from mammals [24][25][26].The proboscis extension response (PER) training is an effective experimental procedure that allows testing individual honeybees for a whole range of associative, non-associative and rule learning related paradigms in well controlled sequence of events [27,28].However, learning and memory retention scores are expressed in all these studies as average scores of the groups of animals.Motivated by the study of Gallistel and coworkers [29], Pamir and colleagues [30] analyzed a large dataset on PER conditioning in the honeybee and found that once an animal started eliciting the conditioned response (CR), it retained the CR in the subsequent conditioning trials with a high probability.Importantly, the gradual acquisition in the groups of animals that results from the recruitment of individuals to the learner population, as predicted by Gallistel et al. [29], was also found for the honyebees.Additionally, individual animals switched abruptly to a learned state during the acquisition phase after several acquisition trials, and then continued performing as learners [30].Further, another study reported that bees switching to the state of learner early during the acquisition phase performed better in the tests of memory retention, even though they did not discriminate better between the rewarded and unrewarded conditioned stimuli during the differential conditioning [31].These results have revealed some of the important aspects of learning behaviour of individuals, which were previously unnoticed in the analyses of groupaveraged learning data.However, no clear picture has emerged on the different classes of honeybees that learn odours with different efficacies.
We therefore attempted to quantify this population heterogeneity, with the goal to understand the factor(s) that structures this among-individual heterogeneity, by a sequential conditioning procedure that allows for finer classification of individuals on the behavioral level based on their performances during two subsequent phases of differential olfactory conditioning and several retention tests; evaluating the strength of olfactory learning and sensitivity.Gene expression patterns were then determined in the brains of the honeybees that belonged to different performer classes.
Multiple studies had previously linked different behaviors in honeybees to their genetic background, such as the developmentally-regulated division of labor in the colony [32][33][34], learning and memory processes in adult bees [35][36][37], egg-laying and foraging behavior of the workers in queenless colonies [38], specialized tasks such as guarding and undertaking [39], foraging [40,41], scouting and recruiting [42] as well as hygienic behavior against the Varroa mite [43].It has also been shown that manipulation of transcript levels in the honeybee brain influences memory performances [44][45][46].Furthermore, a comparative transcriptomic study in the wasp, Polistes metricus also revealed genetic associations of behavioral tasks and substantial overlap of genes involved in the division of labor in the honeybee and foraging in the wasp [47].In the harvester ant, Pogonomyrmex californicus, large differences in gene expression, including genes involved in neuronal functions and chemosensory perception, also have been reported for those queens showing different degrees of aggressive behavior [48].We therefore selected the honeybees based on our previous work on their multiple forms of learning and memory functions and the brain parts involved in these performances [1,49], particularly the mushroom bodies for genetic analysis, which are known to be involved in higher order processing of information, learning and memory in both the bee and Drosophila [50][51][52].The aim of our study was to investigate the nature of behavioral heterogeneity that prevails in the natural populations of honeybees and to access whether the heterogeneity is associated with the expression patterns of genes in the brains of honeybees.

Honeybee lines
A backcrossed genetic line of honeybees (Apis mellifera carnica), selected for improved hygienic behavior, was used in this study (La ¨nderinstitut fu ¨r Bienenkunde (LIB), Hohen Neuendorf, Berlin [53].This line was generated by artificially inseminating a queen from a non-hygienic (+/+) colony with the sperms from the hygienic drones (H).Heterozygous (H/ +) queens were selectively raised from the F1 progeny and backcrossed with the sperms of their hygienic (H) paternal populations.These backcrossed heterozygous (H/+) queens then were used to raise four colonies, (colony 67, 73, 98 and 299) characterized by worker bees with homozygous (H/H) and heterozygous (H/+) genetic backgrounds for the hygienic trait.It is expected that these backcrossed colonies have more homogenous genetic backgrounds than the natural colonies, which is advantageous for the behavioral and genetic investigations.Indeed, previous studies have revealed that heterogeneity in olfactory learning performance of honeybees has a genetic basis [54][55][56].
The experiment took place between July and October 2010.During the course of our experiment, colony 67 became aggressive and was quarantined.From the time point of its confinement, we excluded a total of 32 bees from this colony from the data analysis.

Procedure for sequential olfactory PER conditioning
Forager honeybees were caught at the entrance of the colonies on the day before sequential conditioning.Bees were harnessed in small plastic tubes as described elsewhere [27,57,58], fed until satiation with 30% (W/V) sucrose solution and kept overnight inside a humid Styrofoam box in the dark.The next morning, they were placed in front of the experimental arena for 30-45 min for adaptation.The sequential conditioning protocol consisted of two identical phases, each composed of one round of differential conditioning followed by two sets of retention tests (Fig 1).During the differential conditioning, bees were trained to learn and discriminate an odor (conditioned stimulus, CS+) paired forwardly with the presentation of a sucrose reward (a 30% sucrose solution was used as unconditioned stimulus, US) from a non-reinforced odor (CS-) in the course of 12 consecutive conditioning trials (6 CS+ and 6 CS-).The two pure odor stimuli were presented alternately starting with the CS+ with an inter-trial interval (time gap between consecutive CS+ and CS-trials) of 8 min.During the CS+ trials, the odor was presented for 5 s and the US was applied 3 s after the odor onset by first touching the antenna with a toothpick soaked with the sucrose solution to elicit proboscis extension, followed by feeding through the proboscis for a total of 4 s.Thus, there was an overlap of 2 s between the odor stimulation and the US.During the CS-trials and the retention tests, the odors were presented alone for 5 s.During conditioning, a PER elicited during the first 3 s of odor stimulation before the onset of the US was considered as a conditioned response (CR).A 20 min break was given after the differential conditioning, which followed two sets of retention tests.In each test, bees were presented with increasing concentrations of the CS+ and CSodors on filter paper (10 −3 and 10 −2 dilutions in paraffin oil and the pure odors) starting with the CS+.Each set of retention test was terminated by stimulating the animals once with paraffin oil and once with the filter paper alone.There was no time gap between the two consecutive sets of retention tests.The second round of differential conditioning began after a pause of 30 min and a different odor pair was used (Fig 1).The ability of the bee to elicit a CR was evaluated by visualizing a PER after stimulation of the antennae with sugar water at the end of the fourth set of retention tests.Only bees that elicited PER were considered for data analysis.Bees were sacrificed immediately after this fitness-test, by placing them at -20˚C for 10 min followed by preservation at -80˚C until the gene expression study.The whole procedure took a total time of 6 h and 30 min to complete.
The protocol includes two phases, each composed of one differential olfactory conditioning and two sets of retention tests.Conditioning was performed with the pure odors (CS+, CS-) and the retention tests with increasing odor concentrations, terminating with the pure odors (10 −3 , 10 −2 , pure).Responses to the filter paper (FP) and paraffin oil (Oil) were also recorded.The time intervals between the conditioning and retention tests and between the two phases are shown.
We used two volatile odors, β-ocimene or OM (purity >90%, Sigma Aldrich, Germany) and phenethyl acetate or PEA (purity 99%, Sigma-Aldrich, Germany), for the first phase of the procedure.These odors were reported as specific for healthy brood (β-ocimene) and larvae infected by the chalkbrood pathogen, Ascosphaera apis (phenethyl acetate) and hygienic bees discriminate well between the volatile odors released from the healthy and diseased broods https://doi.org/10.1371/journal.pone.0304563.g001[59][60][61][62].We used two fatty acids, oleic acid or OA (purity >99%, TCI-Europe NV, Belgium) and linolenic acid or LA (purity >70%, TCI-Europe NV, Belgium) for the second phase of sequential learning that are reported to be involved in the process of kin-discrimination in honeybee [63].These two odors have similar chemical structures and thus were expected to be harder to differentiate.During the two rounds of differential conditioning, odor contingencies were reversed to understand whether each odor can be learned well either as a CS+ or as a CSby the bees.Odor stimuli were delivered with a 20 ml plastic syringe containing a filter paper soaked with 10 μl of the respective odor.

Analysis of behavioral data-Quantification of an individualized performance score
A total of 171 honeybees were included in the data analysis.We scored seven behavioral variable to evaluate the olfactory learning performance of the individual honeybees in our assay.Acq-1 and Acq-2: These two variables scored for the rate of CS+ odor learning respectively during the acquisition trials of the first (Acq-1) and the second (Acq-2) differential conditioning.The total number of CRs to the CS+ odor was normalized to the maximum number of CS + responses possible (maximum 5-responses) for an individual bee between the 2 nd and the 6 th CS+ trial.The PER responses to the CS+ and CS-odors during their respective first conditioning trials were not conditioned and thus were not considered.The scores ranged from 0 to +1 for these two variables.DisCond-1 and DisCond-2: These two variables scored for odor discrimination ability of the bees respectively during the acquisition trials of the first (DisCond-1) and the second (DisCond-2) differential conditioning.They were calculated as the total number of PER to the CS+ stimuli that were not followed by a response to the CS-stimuli, normalized to the maximum number of possible responses (maximum 5-responses) between the 2 nd and the 6 th acquisition trials.For variable-1 to -4, we only scored the responses of the bees starting from the 2 nd conditioning trial of the CS+ and CS-because a response was considered to be a conditioned response only after the 1 st conditioning trial.The scores ranged from 0 to +1.DisT-1,-2 and DisT-3,-4: These two variables scored for odor discrimination during the first and second sets of retention tests (DisT-1,-2) of phase-I, and the third and fourth sets of retention tests (DisT-3,-4) of phase-II of the sequential conditioning.The total number of responses to the CS+ stimuli (2 dilutions and pure odors) during the retention test that were not followed by a response to the CS-stimuli (2 dilutions and pure odors) was scored and normalized to the maximum number of possible responses (maximum 3).The scores ranged from 0 to +1.An individualized performance score (P-score) was calculated for each bee by summing the scores of these six behavioral variables.These scores ranged from 0 to +6.

Statistical analysis of behavioral data
Conditioned responses of the bees were analyzed by repeated measures ANOVA followed by the Tukey HSD posthoc test [64][65][66][67].P-value was corrected during the posthoc test following the Bonferroni procedure.Wilcoxon matched-pairs test was applied to compare the scores of the behavioral variables for the same set of bees between the two phases of sequential conditioning.The total number of PERs to the CS stimuli during the conditioning and retention tests was compared between the two phases using the G-test.Lilliefors test was performed to understand the nature of distribution of P-scores in the dataset.Mann-Whitney U test was used to compare the scores of behavioral variables between the high and low performing bees.Spearman's rank order correlation coefficients were calculated between the scores of the seven variables to understand the nature of associations between them.Hierarchical multiple regression analysis was performed to find out how much of the variability in P-score is explained by each of the six variables.Statistica version 5.0 (Dell Software, USA) and SPSS version 23 (IBM) were used for statistical analysis and MatLab (MathWorks, USA) for figure making.The analysis of significant gene overlap between the present study and [68] was performed with a χ 2 test (www.socscistatistics.com/tests/chisquare/default2.aspx).The test was applied by considering that the honeybee genome comprises 13440 genes.The differences were considered statistically significant when p < 0.05.

Tissue preparation and total RNA extraction
The head capsule of the bee was separated from the rest of the body and a window was cut in the cuticle.The manipulation was performed on a metal slide cooled on dry ice.The head was partially lyophilized in a vacuum chamber at -20˚C for 2 hours and then fixed with Tissue-Tek on a metal slide cooled on dry ice.The brain was exposed by scraping off the hypopharyngeal glands covering the brain and the dorsal part of the central brain containing the mushroom body was recovered and immediately homogenized in 400 μl Trizol (Life Technologies, Schwerte, Germany) in a Teflon glass homogenizer.After phase separation, the aqueous phase containing the total RNA was recovered, mixed in equal volume with 100% ethanol, loaded on a column (RNeasy MinElute Cleanup, Qiagen, Hilden, Germany), treated with DNAse (Qiagen, Hilden, Germany), washed and eluted in water.Extracted total RNA was stored at -80˚C until further use.

Microarray analysis
The methodology described by Gempe and colleagues was applied [43].The honeybee wholegenome oligonucleotide microarray (Design: UIUC Honey bee oligo 13 K v1, Accession: A-MEXP-755), containing 28,800 oligos that represent 13,440 genes was derived from annotations of the entire honeybee genomic sequence [69].Ten microliters of each total RNA sample was amplified prior to labelling following the manufacturer's protocol (MessageAmp II aRNA Amplification kit, Ambion).We hybridized 3 μg of each amplified labelled RNA sample to a single microarray slide.The slides were scanned (SureScan Microarry scanner, Agilent, Waldbronn, Germany) and raw hybridization signals were extracted (GenePix Pro 6.0 software, Agilent Technologies).Transcription-level data were processed and analyzed using the LIMMA 2.16 software package (https://www.bioconductor.org/packages/3.3/bioc/html/limma.html).The quality of hybridization was evaluated using the raw expression data from control probes spotted on each slide.Transcription-level data were corrected for background signal ("normexp" function) [70] and intensity-dependent bias was detected ("normalize within arrays" function with the default print-tip LOWESS normalization method) [71].Finally, the log-transformed expression ratios were calculated and duplicate spots were averaged using the "avedups" function.We used a design matrix that incorporated the P-score, the colony conditions and linear models using the Bayesian fitting option.All microarray data were MIAME-compliant, and the raw data have been deposited in a database (ArrayExpress, EMBL-EBI-https://www.ebi.ac.uk/arrayexpress/experiments/E-MTAB-7909/).Differences in gene transcription that resulted from behavior or from the type of backcross were specified as separate contrasts using linear models.P-values were adjusted for multiple testing with a 5% false discovery rate (FDR).Functional annotation of gene sets that fell into similar categories of gene ontology (GO) terms for molecular processes and biological functions were identified using DAVID (http://david.ncifcrf.gov/)[72,73] and KEGG [https://www.genome.jp/kegg/],which included an enrichment analysis of GO terms.We used the gene annotations from the UIUC Honey Bee oligo 13 Kv1 annotation file.

RT-PCR
In addition, confirmatory expression analyses were also performed for selected transcripts using RT-qPCRs according to the previously published protocol [74].Primers are listed in S1 Table and actin was chosen as reference gene (TIB Molbiol, Berlin, Germany).Fold change was calculated with the 2 -ΔΔCt method.Statistical differences were evaluated using Student's t test for independent measures.

Results
We trained two independent groups of bees in the sequential conditioning paradigm in two phases with each composed of one differential conditioning.Two different pairs of odor were used in each phase; PEA and OM were used for the first phase, LA, and OA for the second phase.To balance odor contingencies, Group-1 (N = 85 bees) bees were trained to discriminate between PEA as the rewarded odor (CS+) and OM as the unrewarded odor (CS-) in the first phase of sequential conditioning while in the second phase, LA was rewarded (CS+) and OA was unrewarded (CS-) (Table 1).Group-2 (N = 86 bees) was conditioned with OM as CS + and PEA as CS-, followed by OA as CS+ and LA as CS-.We first compared the performances of the two groups during the two differential conditioning experiments and found no significant difference between Group-1 and Group-2 when the responses, for each group, were pooled for the two phases, as revealed by a 2 × 2 × 6 (group × stimulus CS+/CS-× trial) ANOVA for repeated measures (F 5,3400 = 1.03, p > 0.05).We then pooled the data of all of the four sets of retention tests for each odor-group and compared these groups.The results of a 2 × 2 × 3 (group × stimulus CS+/CS-× trial) repeated measures ANOVA showed no significant interaction effect between group, stimulus and trial (F 2,2728 = 1.38, p > 0.05) for the retention tests.Since no significant differences were found in the overall responses of the bees of the two groups to the CS odors during the two differential conditioning and four sets of retention tests, we pooled the data of the two groups for all further analyses.

Learning and memory performance of the honeybees in the sequential conditioning
In the first differential conditioning, bees learned the discrimination task (CS+ vs. CS-) successfully.A 2 × 6 (stimulus × trial) ANOVA for repeated measures produced significant stimulus and trial effects as well as a significant interaction effect (F 5,1700 = 99.32,p > 0.01) (Table 3).Similarly, significant stimulus, trial and stimulus × trial effects were also found for the second differential conditioning (Table 2).Hence, the results of the repeated measures ANOVA test demonstrate that bees learned the discrimination task well between the CS+ and the CS-odors during the two differential conditioning (Fig 2A and 2F).The response levels to the CS+ increased across the conditioning trials as the responses to the CS-decreased.The Table 1.Honeybee groups used for conditioning in the sequential procedure.

Differential responses during the two phases of sequential conditioning
We found that honeybees, which underwent sequential conditioning, learned the odor stimuli but they showed differential responses in the two phases.To investigate the differences in conditioned responses between the first and second phase, we compared the odor discrimination scores quantified for the acquisition trials of the first (DisCond-1) and the second (DisCond-2) differential conditioning.We found that bees showed significantly higher odor discrimination in the first compared to the second differential conditioning (Wilcoxon matched-pairs test: Z = 4.11, T = 2715, p < 0.01) (Fig 3A).This was due to a significantly higher total number  ).
The box and whisker plots (displaying the median values, upper and lower quartiles, and the minimum and maximum values of the variables) are showing the distributions of odor discrimination during the first (DisCond-1, red) and the second (DisCond-2, blue) round of differential conditioning (A) and together during the first and second (DisT-1,-2; red), and the third and fourth (DisT-3,-4; blue) sets of retention tests (B).Asterisks denote statistically significant differences (p < 0.01).The number of animals is indicated (N).

Classification of the honeybees into various performance categories based on their individual P-scores
Overall performance levels of the individual honeybees in the sequential conditioning procedure are represented by a cumulative performance score (P-score), which is the summation of the six quantified variables (Fig 4).The distribution of P-scores was found to differ significantly from a Gaussian distribution (Lilliefors test: p < 0.01) (Fig 4).It varied between 0 and 5.46 with a mean value of 3.12 +/-1.34SD.Only 27 bees (*15.8%)scored � 1.78, (< 1-SD of the mean), and 28 bees (*16.4%)scored � 4.46(> 1-SD of the mean) whereas the majority of bees (116 bees, *67.8%) had scores between 1.78 and 4.46, corresponding to +/-1-SD of the mean score.We compared the performance levels of the bees with low (P-score<2.6)and high P-scores (P-score> = 2.6), which were categorized as low and high performers (Fig 4).A cutoff P-score of 2.6 was selected based on the lowest score of the high performers and the highest score of the low performers that resulted from the selection of bees for the microarray analysis (Figs 4 and 6).The acquisition and retention test performances were substantially better for the high performer bees than the low performers throughout the two phases of sequential conditioning (Fig 5, see Table 3 for the statistical analysis of between-group performance).The low performer bees showed an overall reduced number of responses to the CS odors throughout the procedure.The most notable effect was observed during the first differential conditioning, with overall very low response levels and no discrimination between the pure odors; as well as between their dilutions during the retention tests-2 and -4 (Fig 5, see S3 Table for the statistical analysis of within-group performance).We applied the Mann-Whitney U test and compared the scores of the six variables viz., Acq-1, Acq-2, DisCond-1, DisCond-2, DisT-1,-2, and DisT-3,-4 of between these two categories of bees.The results showed that high performer bees scored significantly higher than the low performers for all of the six variables (Table 3).
We also separately analyzed the responses of these two groups of bees to filter paper and paraffin oil.We found that high performers also showed significantly higher total number of responses to these stimuli than the low performer bees (G = 24.48,df = 1, p < 0.001) (Fig 5).These results demonstrate that high performer bees performed consistently better than the low performers throughout the two phases of sequential learning and retention procedures.The high performers showed early and consistent responses to the CS+ odors and discriminated well between the CS+ and CS-odors during the two differential conditioning.They also successfully discriminated between the different concentrations of the CS+ and CS-odors during the four sets of retention tests.In addition, high performer bees also presented a higher responsiveness to filter paper and paraffin oil stimulation than the low performers.
The bars represent the number of bees constituting all of the P-score categories.The dotted black line represents the threshold in P-score (2.6), used for the selection of low performer (P-score<2.6) and high performer bees (P-score = >2.6).N represents the total number of animals.
Line  A sample of 24 bees (green hexagons with bee-identification numbers) out of the 171 bees (red squares) was selected for the genetic analysis as shown in the scatter plot.Among these 24 bees, 12 were high performers with P-scores ranging between 2.6 and 5.4 and the other 12 were low performers with P-scores ranging between 0 and 2.5.

Predictor analysis of the P-score
Next, we investigated how well the scores of the six quantified variables, related to olfactory learning and memory, predict the P-score.We first looked into the distributions of the six variables and found their distributions to deviate significantly from a Gaussian distribution (Lilliefors test: p < 0.01).Thus, we performed nonparametric correlation analysis between these variables and the P-score and found a high correlation (Spearman's ρ = 0.88) between the speed of acquisition of the CS+ (Acq-1) and odor discriminability (DisCond-1) during the first phase of differential conditioning (S4 Table ).The speed of acquisition of the CS+ and odor discriminability during the first phase also showed high correlation coefficient values with the P-score (both showed a ρ-value of 0.85) (S4 Table ).The variable, DisT-1,-2 showed the next highest correlation with the P-score (ρ = 0.7) followed by Acq-2 (ρ = 0.69) (S4 Table ).
To understand the capacity of each of the six variables to predict the P-score, we performed regression analyses to find out how much of the variability in P-score is explained by each of the variables separately.Due to significant correlations between all of the six variables (S4 Table ), hierarchical multiple regression analysis was performed to understand the capacities of the individual variables to predict the P-scores.The six regression models, incorporating all of the variables, showed that Acq-1 alone explained the maximum amount of variance (73.3%) in the P-score (Table 4).Like, Acq-1, Acq-2 also explained a significant amount of variance (12.8%) in the P-score (Table 4).Thus, together the rates of acquisition of the CS+ odors during the two phases of sequential conditioning explained 86.1% of the variability in the P-score.The other four variables also contribute to the variance of the P-score: DisCond-1 (3.5%), Dis-Cond-2 (3.7%), DisT-1,-2 (4.2%) and DisT-3,-4 (2.5%) (Table 4).The six variables significantly increased the predictive capacities of the subsequent regression models (Table 4) and together they explained 100% variance in the P-score (shown by the R 2 -value of Model-6) because Pscores are unweighted summations of these six variables (Table 4).These results illustrated that the rate of learning of the CS+ odors during the first phase of differential conditioning or Acq-1 had the highest predictive capacity for the P-score.

Microarray analysis of high and low performer bees
Twenty-four bees were selected for the microarray analysis based on their performance during the behavioral tests, the colony and the date of the experiment (Fig 6 and Table 5).Bees from the colonies 67, 98, 73 and 299 were used for this analysis.The total RNAs extracted from the mushroom bodies were processed to generate the microarray probes.Two sets of probes were generated for each bee, one labelled with Cy5 and the other with Cy3.In this manner, each animal was analyzed twice by swapping colors on each array, to reduce dye-related artefacts.The last bee of the list was analyzed with the first one, thus the analysis scheme was designed to build a loop.High performer bees had P-scores ranging between 2.6 and 5.4 and low performers between 0 and 2.5 (Fig 6 and Table 5).Eight bees each from the colonies 98 and 73 and four bees each from the colonies 67 and 299 were analyzed (Table 5).The cut-off value of 2.6 is subjective between high and low performers.Indeed, the selection of bees for the gene expression analysis was performed by keeping the colony and seasonal effects, which are other sources of variation in gene expression, as low as possible.To this end, high and low performers bees selected from the same colony at the same time of the year were compared on the one array slide (Fig 6 and Table 5).To reduce the colony effects further on the identification of genes, high performer bees with P-scores ranging from 3.9 to 5.3 from the four colonies were also compared on 12 arrays (S5 Table ).In addition, each colony was also analyzed with three bees collected from July to October 2010 and dye swap loop design was applied (S5 Table ).
The microarray analysis of the behavioral design generated a list of 687 genes and the analysis of colony specific differences (genotype design) generated a list of 265 candidates.The behavioral list was sorted by removing 32 gene candidates, which were also present in the genotype list, reducing the behavioral list to 652 candidates that were differentially expressed between the high and low performer bees with adjusted p-level varying between 0.0025 and 2.9×10 −7 (S6 Table ).Most of the identified genes only showed slight differences in expression levels.The expression of 371 genes was reduced in the high performer bees, varying between 0.75× and 0.95× and the expression of 281 genes was increased between 1.05× and 2.46×.A list Selection of the pairs of bees from the four colonies (Colony), conditioned together (on the same day) in the different months of the experimental year (Month), for the behavioral analysis was based on their high and low P-scores.The column 'Bee' is showing the identification number of the bees.Processed total RNAs were labeled with Cy3 and Cy5 for microarray analysis.
https://doi.org/10.1371/journal.pone.0304563.t005 of 45 identified genes, including the candidates with the highest changes in expression levels, is presented in Table 6.They potentially encode proteins implicated in diverse biological functions, including neurotransmitter receptors for dopamine, tyramine, GABA and acetylcholine, proteins implicated in the synaptic release of neurotransmitters (VAMP, syntaxin, neuroligin) and amino-acid transporters (GABA, excitatory amino-acid transporters).Several candidates were implicated in cellular processes important for memory formation (proteasome complex, TOR pathway, protein translation) and in development (the transcription factor mblk-1, minibrain).Other interesting candidates, including calmodulin, kinases, and phosphatases were also implicated in messenger pathways important for memory formation.Interestingly, genes, which are regulated by growth hormone or juvenile hormone and malvolio, a gene implicated in the division of labor, also showed increased expressions in the high performer bees.Surprisingly, RNAs known to be absent from the mushroom body were also found in the list.The blue opsin gene, normally expressed in the visual system probably originated from a contamination with the ocellar tract.In addition, genes including, α-amylase, major royal jelly protein 5 (mrjp5) and glucose oxidase also showed the strongest differential expression with their expression levels doubled in the high performer bees.These RNAs probably originated from a contamination with the hypopharyngeal glands, , which are adjacent to the MB.Some genes of the whole genome honeybee microarray could be linked to the orthologues of the Drosophila Flybase.This way, a list of 418 Flybase entries were generated from the sorted honeybee behavioral list (S6 Table ).KEGG Automatic Annotation Server (KAAS) analysis revealed 600 enriched terms, including those related to sugar metabolism (e.g.pentose, fructose, galactose and sucrose), lipids metabolism (e.g.sphingolipid, arachidonic acid, linoleic acid), amino acid and nucleic acid metabolism (e.g.spliceosome, RNA transport, Ribosome) (S7 Table ).The results also highlighted terms directly related to classical neurotransmitter systems (GABA, acetylcholine, glutamate, dopamine serotonin), olfactory transduction, axon guidance, synaptic transmission (e.g.SNARE, Tight junction, GAP junction), the proteasome complex and several signaling pathways (e.g.Ras, MAPK, Ca 2+ , cAMP, cGMP, mTor).Many terms related to hormonal control (e.g.insulin, oxytocin, prolactin, adrenergic) and food processing (salivary-, gastric acid-, pancreatic-secretion, etc.) were also described.Several terms not directly related to the neuronal metabolic pathways underlying memory processes were also highlighted (e.g.immune response, longevity).
The behavioral gene list was further submitted to DAVID and the functional annotation chart revealed 76 statistically significant terms.Some of these were related to microtubuleassociated complex, mitochondria, neurogenesis, plasma membrane, splicing, GTP activity, oxidoreduction and nucleic acid metabolism (S8 Table ).The first eleven clusters, presenting enrichment scores between 1.45 and 0.81, highlighted terms related to small GTP-binding protein domain, mitochondrion, PDZ domain, glycoside hydrolase, α/β hydrolase, spliceosome, protein phosphatase, actin cytoskeleton, neuroblast, small GTPase superfamily and lipid metabolism (S8 Table ).

Validation of gene candidates by RT-PCR
Two low and two high performer bees from each of the four colonies (S9 Table ) were selected for the validation of nine selected genes by reverse transcription quantitative real time PCR technique (RT-qPCR) (Table 6).Low performers had scores varying between 0 and 3.33 and high performers between 3.8 and 5.4.These bees were collected between July to October 2010.None of the tested genes showed significant difference between the high and the low performer groups (Student t-test, p > 0.05) (S1 Fig).A selection of 45-identified genes from the sorted behavioral list that includes candidates with the strongest changes in expression levels between the high and low performers.The table shows the accession numbers and identity of the mRNAs deduced from the genome annotation (Honeybee mRNAs), the change in expression levels between the high and low performers (Fold Change), the averaged expression levels measured on the microarray and the adjusted p-values denoting significant changes in expression levels between the high and low performers, calculated from the microarray analysis.Selected transcripts were further validated by RT-qPCR (S1 Fig) https://doi.org/10.1371/journal.pone.0304563.t006

Discussion
A deeper understanding of the variations in individual behavior within species is attracting more attention owing to their significant contributions to the ecological traits of the species and their evolution [9,75,76].Heterogeneity in behavioral types or personality, which is related to the difference in cognitive ability of individuals, is a heritable feature in animals that controls important life-history traits, such as food intake, detecting prey, predators and mates, labor division, ability to respond better to the environmental change, and overall fitness [76][77][78][79][80]. Honeybees show individual variability in behavior, cognitive task specialization and consistency in learning proficiency across learning paradigms and stimulus modality [81][82][83][84][85]. Furthermore, an early study in honeybee showed substantial genetic contribution to the variance observed in learning phenotype, which supports the heritability of among-individual difference in cognitive ability [22].Hence, this species is an appropriate model system to study the nature of heterogeneity in cognitive abilities.In honeybee, it is already established that the group-averaged analysis of the learning data impedes the visualization of this heterogeneity in individual's learning, viz. the learning rate, final level of learning, and the strength of retained memory [30,31].Previously, response latency has been described to have the potential to underlie the heterogeneity in learning dynamics of the individual bees, which was explained by a heterogeneous Rescorla-Wagner model [30,31].However, these analyses could not capture a clear picture of the identity of a salient behavioral feature(s) that determines the individualistic variability in learning ability in the honeybee populations.Here, we addressed this issue through a complex form of olfactory learning protocol and analyzed the data to understand the contributions of different learning-related features to the overall performance levels of the individuals.Further, correlations between learning performance and the expression patterns of learning-related genes in the brains of these individual honeybees were evaluated.

Methodological consideration
We trained the honeybees in a sequential olfactory conditioning protocol in which the animals received two consecutive stages of differential training with two different pairs of odors with each differential conditioning stage followed by two sets of memory retention tests in which we applied several dilutions of the rewarded (CS+) and unrewarded (CS-) odors.The protocol was specifically designed to test the responses of the individual bees rigorously for a longer period to find out the strength in their learned responses.The tasks of discrimination between the CS+ and CS-odors during the two stages of conditioning were easier compared to the four sets of retention tests because of the use of pure odors instead of their dilutions during the conditioning; thus, the 'Aha' effect of learning led to clear separation in responses of the bees to the odors [86,87].Individual honeybees however had to perform consistently well throughout the procedure to reach a higher cumulative performance score (P-score) or they failed at different time points to score lower.This way we were able to distinguish between the low and high levels of olfactory learning performers present in the natural populations.Bees from all four colonies learned the odor stimuli during the conditioning and discriminated well in the retention tests.However, the discrimination task was more difficult when the cuticular volatile fatty acids (OA and LA) were used compared to the brood volatiles (OM and PEA).We speculate that learning of differential contingencies of these two odors is difficult in the perceptual space of honeybees owing to the structural similarities of these two compounds (both are 18-carbon, unsaturated essential fatty acids) [88].This could be a prime reason for the bees to generalize more between these two odors during the second half of the procedure in addition to the possible involvements of the non-associative components, like sensitization and hunger due to prolong testing of the animals for more than 6 hours with repeated sucrose stimulations.The associative components together with the non-associative counterparts make our procedure more stringent and specific in discriminating between the low and high levels of learning performances.We quantified the speed of learning of the CS+ odors and the levels of odor discriminability both during the two sets of conditioning and four sets of retention tests and investigated whether there is a single salient feature that persistently contributes to the overall learning performance of the individual bees.This could be considered as an important feature that controls the heterogeneity in learning dynamics among the individuals of the natural populations.There are other complex olfactory learning protocols that were previously applied on the honeybees [66,89] but that we did not use because these protocols were either not recording the bee's responses for a longer period or they incorporated extinction components in learning, both of which were avoided.It is conceivable that an extinctionlearning component prevailed in our sequential conditioning design as we performed two consecutive sets of retention tests after each of the two differential conditionings.Indeed, bees responded less to the highest dilution of the CS+ odors (10 −3 ) in the second compared to the first set of retention test (in the paired-retention tests).However, no difference in memory recall was found for the training concentrations of the CS+ odors in the two consecutive sets of retention tests.Thus, we conclude that no extinction learning took place in our experimental design in which we offered no time gap between the consecutive sets of retention tests.

Speed of learning explains the heterogeneity in learning and memory performance of the individual honeybees
Very few studies have systematically investigated the nature of heterogeneity or individualistic variation in the learning ability of honeybees.Tait and colleagues have reported individualistic variation in honeybee in the domains of olfactory and navigational learning as well as maneuverability performance; however, the authors did not conclude anything with regard to the source of variation in the learning behavior of individuals [23].Our analysis clearly demonstrates a large heterogeneity in the learning ability of the individual honeybees, which is structured with a continuous range of performance.Further, the variation in the speed of odor learning appears as the major component for this heterogeneity.We found a non-normal distribution of the P-score in our experimental population, which may be due to low sample size.Interestingly, normal distribution was reported for several cognitive capacities in human populations [90].The distribution of P-score in our study might reflect different degrees of specialization for particular learning modalities among the honeybees, which are at different developmental stages and are involved in different tasks [82,91,92].Indeed, variation among the individuals can emanate from the development of specific experience-dependent biases and thus, differences in physiologies to detect, learn and respond to particular external stimuli [93].The large heterogeneity in the learning dynamics of individual honeybees as seen in ours as well as in the previously generated learning data sets [30,31] can probably be attributed to several factors, such as stimulus sensitivity, speed of learning, and memory retention capacity.We found that high performer bees performed consistently well throughout the 56 trials of odor conditioning and retention tests compared to the low performers.They have mastered the tasks of odor discrimination during the conditioning and showed specific responses to the rewarded odors during the retention tests as well as responded more to the air puff stimulations without any odors (filter paper only or with paraffin oil) than the low performers (Fig 5).Importantly, the responses to the filter paper also varied between individual bees and thus contributed to the individual variability.These results together with higher responses of the high performers compared to the low ones to the highest CS+ dilution (10 −3 ) during the retention tests indicate possible higher olfactory sensitivity among the high performers.
The feature that correlated most strongly with the P-score was a higher acquisition rate during the first phase of conditioning.The speed of learning of the CS+ odor during the first differential conditioning (Acq-1) best predicted the cumulative performance of the individuals among all six quantified variables.Individuals with higher speed of CS+ odor learning not only discriminated the odors well during the first phase of differential conditioning, they also had higher retention scores and higher odor sensitivity.Although the honeybees also performed successfully in the second differential conditioning phase, only the speed of acquisition during the second differential conditioning (Acq-2) correlated strongly with the speed of acquisition during the first differential phase.Acq-2 was also the second strongest predictor of the cumulative score.The discrimination index during the second acquisition phase (Dis-Cond-2) and the discrimination during the second retention test (DisT-3,-4) correlated with Acq1.These variables are more influenced by the capability of an individual to adapt its behavior to new conditions and it adds an independent value to the selection of high performers.
Fast and slow learning is elaborately discussed in behavioral ecology in the context of behaviors that involve speed-accuracy trade-off.A fast learner can take more risk to gain more only for a short-term and even perhaps lacks accuracy, whereas a slow learner has better accuracy and it can perform the duty more safely with the notion for a long-term gain [75,94].The advantage of having slow learners has also been reported in a study on bumblebees where fast learning of visual information was found to be associated with lower number of days of foraging compared to the slow learning bumblebees.This suggests that superior cognitive abilities may not be beneficial to build the reserve of the colonies [95].
We were unable to test the change of conditioned responses over time shortly after learning [96].However, our results clearly show that bees, which started responding early to the CS + odor compared to the late responders during the first conditioning, showed stronger and more specific odor memories in the short-term recall tests.The findings on memory specificity are contradictory to the findings of Pamir and colleagues [31] who used an absolute training procedure and found no difference in the specificity of odor memories between the early and late responders.We assume that this difference in results could be due to the differences in our conditioning methods and odor stimuli used, because the present sampling method is more controlled.We thus surmise from our results that the mechanisms underlying the speed of learning of odor stimuli are correlated to those controlling their discrimination and consolidations of memories.

Genetic analysis in the mushroom body area of the high and low performers
An analysis of differences in gene expression in the mushroom body (MB) region of the low and high performers was performed holistically using microarray.The analysis identified 652 differentially expressed genes; most of them were characterized by small differences in expression levels between the high and low performer bees.This indicates that many genes each with small amplitude effects may contribute to the high performer phenotype.We selected nine RNAs from the transcript list originally generated by the microarray method to further verify their abundance by RT-qPCR.The presence of all nine transcripts could be confirmed by RT-qPCR, but no significant differences in abundance was detected in the high and low performer groups by this technique (S1 Fig) .These results illustrate that behavioral traits are influenced by small changes in expression levels of multiple genes [97].This particularity is also a feature described for neurological disorders [98,99].
The KEGG and DAVIDD analyses revealed similar terms.Nervous system specific terms were related to plasma membrane and synaptic transmission, illustrated by neurotransmitter receptors (GABA, dopamine, octopamine, acetylcholine, glycine and serotonin) and to synaptic components related to neurotransmitter release (e.g.SNAREs).Terms related to LTP, LTD, and neurogenesis were also revealed as well as terms related to olfaction, taste and photo-transmission.The latter might indicate contamination with tissues adjacent to the MB.Indeed, some of the terms can also be related to increased metabolic activity, like sugar metabolism, or to mitochondria and cellular respiration.The analysis of the gene list also indicates that the MB dissections were contaminated with hypopharyngeal glands as we identified higher levels of the major royal jelly protein 5 gene (mrjp5) in the high performers.Indeed, all mrjps are produced in the hypopharyngeal glands except for mrjp3 that is produced in the Kenyon cells of the MB of nurse bees [100][101][102].Equally, α-amylase and glucose oxidase, both of which were overexpressed in the high performers, are also produced in the hypopharyngeal glands [103,104].However, several intracellular signaling cascades, known to be implicated in learning and memory, were also highlighted, including several second messenger pathways (TOR signaling and the ubiquitin proteasome complex).These pathways regulate, among others, translation-and transcription-dependent processes that are important in memory processes and are highlighted in different lists.However, these processes are not restricted to neurons and the increased RNA abundances might have arisen from the activities in adjacent tissues.
The hypopharyngeal glands are important for food processing in honeybees [105][106][107][108].The genes mrjp1-4 and 7 are more active in the nurse bees compared to the foragers [101,102].In foragers, α-amylase and glucose oxidase are associated with the conversion of nectar to honey [103,104].The change in expression of these enzymes is dependent on the age and the task.Interestingly, malvolio, a protein influencing labor division in honeybees, was also overexpressed in the high performers [109,110].This manganese transporter was discovered in Drosophila and is known to influence sucrose responsiveness [111].Thus, we hypothesize that the high performer bees are in the transition from nursing to foraging activities, which is supported by the highlighted genes and terms related to hormonal pathways.It is also known that immunity changes during this transition [104,112] and terms related to immune response were also highlighted in the gene ontology analysis.
In a previous study, Whitfield [68] and colleagues analyzed behavioral maturation in the honeybee by comparing the gene expression profiles in the brains of bees in different settings; nurses and foragers of specific ages compared in crossed conditions as well as nurses treated with different chemicals, described for their roles in the onset of foraging, Methoprene (an analogue of juvenile hormone), manganese and cGMP.Both settings were reported to associate with the molecular pathways of malvolio and foraging.However, cAMP was found to play no role in the onset of foraging [68].We compared our gene list to the gene lists of them.A comparison of our list of genes with these gene sets, characterized by a significance level < 0.05, did not reveal any significant overlap.By restricting the gene sets to candidates, characterized by a significance level < 0.001, significant effects were described (Table 7).A significant correlation was found for Methoprene treatment, a condition also implicated in the transition from nursing to foraging.However, there was no significant overlap with manganese or cGMP treatment or when comparing 17 days old foragers and nurses (d17F/d17), which are conditions specific for foragers.Significant overlap was found with gene sets specific for juveniles, which are still active with hive duties but able to perform the transition from nursing to foraging (d4/d8 and d12/d17).This indicates that high performers belong to a group of young bees making the transition to foraging activities.The animals were caught at the end of the afternoon, when mostly foragers perform outbound flights.Thus, these bees might correspond to young individuals beginning their foraging activity.It is known that young foragers and precocious forager honeybees learn olfactory stimuli better than the nurse bees [82,113].However, less is known about differences between young and old foragers.
Only a few studies addressed this question [114,115].The correlation with the gene set specific for the cAMP second messenger pathway can explain the higher sensitivity of the high performer bees compared to their low performing counterparts.Indeed, this second messenger pathway is associated to modulatory neurotransmission mediated by the biogenic amines, including serotonin, octopamine, dopamine, and tyramine that are implicated in arousal, sensitivity, and learning performances [82,116,117].The cAMP second messenger pathway is linked to some of these neurotransmitter systems, which might explain why the high performer bees are responding more to olfactory stimulation with filter paper and paraffin oil.Although, the selection of high and low performers was biased for olfactory sensitivity, the differential conditioning procedure demonstrates that learning and memory are specific.The brains of the high performers present differences in several neurotransmitter systems, components involved in synaptic transmission and in signaling pathways that are important for learning and memory.Further investigations are needed to describe more precisely the differences between the young and old foragers.

S1 Fig. RT-qPCR analysis of the nine-candidate genes for the high and low performer bees.
RT-qPCR analysis of expression differences of the 9-selected gene candidates (represented by respective standard IDs): Comparison of the mRNA levels of the selected genes between the high and low performer bees.Each bar represents the relative mRNA levels +/-SE in the low (green) and high (red) performers, normalized to the levels of the low performers.The differences in expression levels measured by microarray are also given (yellow bars).(DOCX) Comparison of the gene list with adjusted p value < 0.001, described in this study (n = 326) and in the study of Whitfield and colleagues (2006) [69].The number of differentially expressed genes identified for specific conditions are presented (Total differentially expressed genes).There were treatments with Methoprene, Manganese, cGMP, cAMP, and cAMP vs cGMP.Age effects were also considered in the nurse bees of different ages (days, d).Seventeen days old nurses and foragers were also compared (d17F/d17).We evaluated if the number of overlapping genes between the sets identified in this study and those of Whitfield (Gene overlap) was random with a χ 2 test.Significant differences have a p value < 0.05. https://doi.org/10.1371/journal.pone.0304563.t007

Fig 1 .
Fig 1.The protocol of sequential conditioning.
of bees were trained along the two rounds of differential conditioning in the sequential conditioning paradigm.PEA and OM were used as CS odors in the first conditioning and OA and LA were used in the second conditioning.The contingencies (+: reinforced, -: unreinforced) of the odors were reversed for the two groups of bees during each of the two rounds of differential conditioning to keep the balance in odor contingency.https://doi.org/10.1371/journal.pone.0304563.t001differences in response levels to the CS+ and CS-increased across the successive trials of the four sets of retention tests, which correlates with the increasing concentrations of the CS odors from the first to the third trial of the retention tests (Fig 2B, 2D, 2G, and 2I, see S2

Fig 2 .Fig 3 .
Fig 2. Acquisition and memory retention performances of the bees in the sequential conditioning.https://doi.org/10.1371/journal.pone.0304563.g002 graphs, show the percent PER to the CS+ (red line) and CS-(blue line) stimuli for the high performer (N = 118) and low performer (N = 53) bees during the two phases of sequential conditioning and retention procedures.The first set of eight subfigures (A to H) represent the PERs to the odors during the first phase.Line graphs of the high and low performers for the same odor trials are given in a pairwise manner, such as Acquisition-1 (High and Low), Retention Test-1 (High and Low), Retention Test-2 (High and Low), Filter Paper and Paraffin Oil (High and Low) to visualize the differences in performance between these two groups of bees.The arrangement of the line graphs is the same for the two groups of bees for the second phase of sequential learning (I to P).Total number of responses to the filter paper soaked with paraffin oil (first data point) and to filter paper only (second data point) are given forthe first (G, H) and the second (O, P) sequential learning phases.The overall responses measured after each retention phases were summed for the high and low performer bees.Repeated measures ANOVA followed by Tukey's HSD posthoc test were performed to analyze the performances of the two groups of bees during the conditioning and retention stages (see S3 Table for

Table 3 . Results of the Mann-Whitney U test. Variable Comparison Between High and Low Performer Bees
mean scores for each of the six variables for the two selected groups of honeybees.https://doi.org/10.1371/journal.pone.0304563.t003

Table 2 . Results of the repeated measures ANOVA for the two rounds of differential conditioning for the pooled data.
Fisher statistic values (pooled data) from the repeated measures ANOVA tests are shown for the two differential conditioning.All p-values correspond to the stimulus, trial and stimulus × trial effects are significant (p < 0.01).Line graphs are, showing the percent PER to the CS+ (red line) and CS-(blue line) stimuli for the pooled population (N = 171) during the two differential conditioning (A and F) and four sets of retention tests (B, D, G, and I) of the sequential conditioning.Responses to filter paper soaked with paraffin oil (first data point) and filter paper only (second data point) are shown in C, E, H and J. Asterisks denote statistically significant differences (p < 0.01) between the PERs to the CS+ and CS-odors.See Table2 and S2Table for statistics.

Table 4 . Results of the hierarchical multiple regression analysis.
The results show the values of R 2 , changes in R 2 , the consecutive regression models with their variables, corresponding model statistics, and the β-values for the six variables.Acq-1 explains the highest, 73.3% variance in the P-score.Together the six variables explained 100% variance in the P-score (shown by the R 2 value for Model- 6).No F-statistic and p-values are shown for DisT-3,-4 by SPSS.For model-1 to -5, statistics are found significant (p < 0.01).https://doi.org/10.1371/journal.pone.0304563.t004