Range and Frequency of Africanized Honey Bees in California (USA)

Africanized honey bees entered California in 1994 but few accounts of their northward expansion or their frequency relative to European honey bees have been published. We used mitochondrial markers and morphometric analyses to determine the prevalence of Africanized honeybees in San Diego County and their current northward progress in California west of the Sierra Nevada crest. The northernmost African mitotypes detected were approximately 40 km south of Sacramento in California’s central valley. In San Diego County, 65% of foraging honey bee workers carry African mitochondria and the estimated percentage of Africanized workers using morphological measurements is similar (61%). There was no correlation between mitotype and morphology in San Diego County suggesting Africanized bees result from bidirectional hybridization. Seventy percent of feral hives, but only 13% of managed hives, sampled in San Diego County carried the African mitotype indicating that a large fraction of foraging workers in both urban and rural San Diego County are feral. We also found a single nucleotide polymorphism at the DNA barcode locus COI that distinguishes European and African mitotypes. The utility of this marker was confirmed using 401 georeferenced honey bee sequences from the worldwide Barcode of Life Database. Future censuses can determine whether the current range of the Africanized form is stable, patterns of introgression at nuclear loci, and the environmental factors that may limit the northern range of the Africanized honey bee.


Introduction
The honey bee subspecies Apis mellifera scutellata, was introduced from southern Africa to Brazil in 1956 in an effort to breed honey bees better suited to the Neotropics [1]. Swarms soon escaped containment and began to hybridize with European honey bees [2,3]. The resulting hybrid offspring are known as "Africanized" honey bees. Africanized bees have since spread throughout much of the Western hemisphere, and the advance of African genes both north and south from the introduction site has become one of the most-studied cases in the hybridization literature (reviewed in [4,5]). The Africanized honey bee reached Mexico by 1985, Texas by 1990 and California by 1994 [6,7]. Since 1994, the Africanized honey bee has Finally, we examine the linkage, if any, between polymorphism at the mitochondrial CO1 locus used for DNA barcoding, and previously used mtDNA polymorphisms for identifying Africanized bees. Use of DNA barcodes for mitochondrial origin of honey bee specimens would facilitate global monitoring of the spread of Africanized honey bees because honey bees are a common component of insect biodiversity samples collected for barcoding.

Methods
No specific permissions were required for collecting activities as all sites were public roadsides and no sampling affected private property. No endangered species were sampled or harmed in this research.

Bee Collections
Worker honey bees were collected in groups of three bees per site. Sites constituted any floral resource upon which honey bees workers were foraging. In urban areas these resources were often blooming horticultural plant species. In rural sites, native or naturalized plant species were the usual floral resource. No collections were made from crop plants. All three workers were usually collected within 10 meters of one another (max. 50 meters), often from the same flowering shrub, tree, or patch of groundcover. Collection of three bees per site allowed us to assess mitotype heterogeneity within foraging sites. North of San Diego County we collected bees throughout California west of the Sierra Nevada range at haphazardly chosen sites usually separated by 30 to 80 kilometers. Ninety-one sites were sampled extending as far north as Southern Oregon. In San Diego County we collected bees from 88 sites usually separated by approximately 8 kilometers. All bees were collected using hand held insect nets, killed in soapy water (one drop of dish soap in 50ml H 2 0), immediately rinsed in fresh water and then pinned to dry or placed in 95% EtOH and pinned to dry within 24-48 hrs. Coordinates and other Information on the collection sites for all specimens can be found in S1a and S1b Table. We also obtained samples of ten workers from each of 23 beekeeper-managed hives in San Diego County. Hive bees came from 11 different beekeepers who represented a mixture of hobbyists and honey bee professionals. For these samples we asked beekeepers if the hive was started by ordering a colony from a commercial source or by capture of a feral swarm, and the time since last requeening of the hive. Commercial colonies are unlikely to be Africanized but swarms could arise from feral colonies and be Africanized. Our beekeepers only requeened with commercially reared European queens and thus requeening also increases the probability that colonies consist of European bees. Worker bees from managed hives were killed by freezing at -5°C, air dried, and shipped to our lab. One sample colony arrived with severe mold and was not used for morphological analysis although the mitotype of one bee was determined (see methods below). Using nets, we also obtained ten bees each from the nest entrances of ten feral colonies located in San Diego County.

Morphological measurements
We followed the methods outlined in [25] to differentiate Africanized from non-Africanized honey bee workers based on morphology. Measurements taken were lengths of the right forewing, right hindwing, and the lengths of the tibia and femur of the right hind leg (character set 2 of [25]). Bee parts were mounted on slides and images were captured using IScapture and measured using imageJ software [27]. For bees sampled at foraging sites, the factor loadings for individual worker bees from Table 6 of [25] were applied (2.5164 Ã forewing length + 1.2159 Ã hindwing length + 16.3439 Ã femur length-10.6356 Ã tibia length-36.4909). For the groups of ten bees sampled from each managed or feral colony, means of each measurement were taken, and the factor loadings for colony means recommended in Table 3

Molecular analysis
The left hind legs of dried bees were crushed under liquid nitrogen, and DNA was extracted using the Qiagen DNeasy blood and tissue kit according to the manufacturer's protocol. To assess if bees had mtDNA characteristic of sub-Saharan Africa, we amplified a region of the mitochondrial cytochrome b gene in 25 μl reactions using primers and conditions described in [24]. After amplification, 10 μl of the PCR product was digested with BglII at 56°C for 3 hours [24]. One band present after digestion and electrophoresis indicated the specimen had African mtDNA, while a two-banded phenotype was scored as European mtDNA [24,28].
To dissect further the origins of San Diego County honey bees, we amplified and sequenced the cytochrome I (COI)-cytochrome II (COII) intergenic spacer from a random subset (N = 48) of honey bee foraging workers using primers E2 and H2 of [29]. The PCR conditions were: denaturation at 94°C for 2 min, followed by 35 cycles of 94°C for 45s, 46°C for 45s, and 72°C for 45s, and then a final extension of 72°C for 5 min. PCR products were visualized on a 2% agar gel and sequenced by Retrogen Inc., San Diego, CA. The DNA barcode region of the mtDNA COI gene was also amplified from these same 48 bees in order to determine associations between COI variation and other markers. Amplifications used the protocol and the primers LCO1490 and HCO298 described in [30]. Sequence alignment and analyses were performed using the MEGA 6 software package [31].

Restriction-site analysis
North of San Diego County we performed the cytochrome b-BglII restriction analysis on 265 bees collected at 91 sites (Fig 1). Twenty-three bees (8.7%), carried the African mitochondria. South of 35.5°N latitude, bees with the African mitochondria were fairly common. North of that latitude, we found African mitochondria in only 6 of 233 bees. These occurred in scattered locations, one near the Pacific coast at the University of California Big Creek Reserve (36.074014, -121.59418), the others either in the foothills of the Sierra Nevada mountains, the Central Valley, or the easternmost regions of the San Francisco Bay Area. The northernmost African mitochondria was detected approximately 40 km south of Sacramento (38.1838840, -121.8433948; S1a Table).
In contrast to the statewide sample, honey bees throughout San Diego County were highly Africanized. Restriction analysis revealed that 64.8% of the 298 bees from San Diego County carried the African mitochondrion. There was no readily apparent spatial pattern to the where Africanized honey bees were collected, with bees carrying the African mitotype found in all habitats (scrub, chaparral, mountain forests, deserts) and at all levels of urbanization from the urban core in coastal regions to undeveloped areas in the eastern parts of the county. In addition, considerable heterogeneity in mitotype occurred at the level of the collecting site. Of the 95 sites from which three bees were collected, 55 contained workers with alternative cytochrome b-BglII restriction phenotypes indicating the presence of workers from different hives attending the same floral resource. The frequencies of sites with 0, 1, 2, or 3 bees with African mitochondria were not significantly different than expected by chance when drawing three bees at random from the observed mitotype frequencies (χ 2 = 7.0, DF = 3, NS).

Morphological analysis
Using the discriminant function for individual workers, 17.84% of the honey bees north of San Diego County were scored as African (Fig 2a). Honey bees with the African mitotype had a mean discriminant function score of -0.428, whereas those with the European mitotype had a    [25]. Bees with scores smaller than the midpoint are classified by the discriminant function as Africanized while bees with scores larger than the midpoint are classified European. The expected misclassification rate of this discriminant function for individual workers is 9.9% [25].  Table). Within San Diego County, 60.85% of honey bees were morphologically classified as African (Fig 2b), a percentage similar to that found based on mtDNA. However, the mean discriminant function score for honey bees with the African mitotype (-1.473), did not differ from that of bees carrying the European mitotype (-1.180; t = 0.99, DF = 232, NS; S2b Table).

Bees from San Diego County hives
Of the 24 managed hives sampled in San Diego County, 3 (12.5%) carried the African mitotype, a significantly lower frequency than found among workers collected at floral resources (χ 2 = 29.07692, DF = 1, P < 0.005). The three hives that contained bees with the African mitotype all originated from feral swarms brought into domestication. Of the ten feral colonies assessed, seven (70%) carried the African mitotype, a frequency not different than that found among foraging workers county-wide (χ 2 = 0.10989, DF = 1, NS; S3a Table).
Morphologically, 5 of 23 honey bee colonies from local beekeepers provided discriminant function scores indicating Africanization (Fig 3). Two of the three hives carrying the African mitotype scored as African by morphological analysis as well as a third hive known to have been started from a feral swarm but which carried the European mitotype. Of the remaining 20 Distribution of morphological discriminant function scores for means of ten workers per colony collected from feral and domestic hives in San Diego County. Arrow indicates the midpoint between the mean discriminant function scores for worker mean measurements from European and Africanized colonies [25]. Colonies with scores smaller than the midpoint are classified by the discriminant function as Africanized while those with scores larger than the midpoint are classified European. The expected misclassification rate of the discriminant function for worker means is 3.5% [25]. hives, the fact that two (10%) provided discriminant scores identifying them as African may reflect some further Africanization, but is not significantly different than the expected misclassification rate (3.5%; [25]) of the discriminant function for colony means (χ 2 = 2.50, DF = 1, NS). The mean discriminant function score of colony workers from managed hives in which the queen carried the African mitochondria was -1.082 while that for hives carrying the European mitochondria was larger 0.2965 (t = 1.7915, DF = 21, P < 0.0438). Workers from feral colonies were surprisingly large given the high frequency of African mitotypes found within them. Only 3 of the 10 colonies scored as African using discriminant function of worker means (Fig 3), and there was no difference in the mean discriminant scores of feral hives carrying the African ( x= -0.3530, N = 7) versus the European ( x= -0.3526, n = 3) mitotypes (t = 0.0006, DF = 8, NS; S3b Table).

Additional sequence analyses
For a random sample of 48 workers collected at floral resources in San Diego County we sequenced the COI-COII intergenic region (Genbank accessions KT279822-KT279869) previously used to distinguish within and among African and European mitotypes (see e.g. [24,32]) and the mitochondrial COI DNA barcode locus (Genbank accessions KT275960-KT276007). Using phylogenetic analysis (S1 Fig), Table) which are lineages from the northern Mediterranean, Middle East, and western Europe, respectively [32,33]. All 29 mitotypes classified as African by restriction analysis were found by sequencing the COI-COII intergenic spacer to represent the A mitotype characteristic of sub-Saharan Africa. Therefore perfect linkage between the cytochrome b restriction site and mitotypes defined using COI-COII intergenic sequences was confirmed [24].
At the DNA barcode region of the mitochondrial COI gene, we found a single nucleotide polymorphism (SNP) at position 2382 as defined in [34] that was in perfectly correlated with the cytochrome b-BglII restriction pattern (Table 1, Genbank accessions KT275960-6007). A cytosine at that position was found in all bees that gave the restriction pattern characteristic of the African mitotype, while a thymine was found at that site for all bees classified by restriction as European. To test whether this SNP is diagnostic on a wider sample of honey bees, we downloaded 401 high-quality Apis mellifera sequences (BIN AAA2326, S5 Table) with locality information from the Barcode of Life Database (http://www.boldsystems.org/) the worldwide repository of DNA barcode data [35]. The SNP phenotype diagnostic of the African mitochondrion was restricted to samples from sub-Saharan Africa, the neotropics, and the southern US Table 1. Linkage between mitotype sequence and restriction site variation in 48 randomly sampled honey bee workers from San Diego County. For cytochrome b, numbers denote one-(African) and two-(European) banded phenotypes after digestion with BglII. For CO1, numbers denote sequences that have a cytosine versus a thymine at position 2382. For COI-COII intergenic spacer sequences, numbers denote the frequency of previously described mitochondrial types (A = Africa, C = Carniola/Slovenia, M = Northern and Western Europe, and O = Oriental/Middle Eastern) determined by phylogenetic analysis (see S1 Fig).

Discussion
Since 1994, Africanized honey bees have expanded their range in California from the two counties bordering Mexico to the northern reaches of the Sacramento River delta. The northernmost site at which we recorded the presence of the African mitotype is approximately 250 kilometers northwest of the previously reported northernmost limit of Africanized honey bees in Mariposa County [7,26]. Africanized honey bees are thought to be capable of expanding their range by as much as 300-500 kilometers per year [6]. Since Africanized honey bees were recorded in Mariposa County in 2006, it would appear that range expansion has slowed considerably from maximal rates. More detailed sampling in subsequent years will be needed to determine both the exact range limit and whether it is stable. Africanized bees appear to have a reduced ability to survive cold winter temperatures and two attempts at predicting the eventual northern limit of the Africanized honey bee have been made. One study [36] used data from the southern New World limit of the range of Africanized honey bees to predict that the northern limit would be set where the average temperature was below 16°C during the coolest month (January). Harrison et al. [26] noted that by 2006 Africanized honey bees were already found overwintering north of this prediction and our study found them even further north. However, as the climate warms, the predicted range limit will also shift north so expansion beyond the predicted limit might be explained by recent warming. Our sampling occurred following the warmest winter on record for California (http://www. ncdc.noaa.gov/). If winter temperatures limit the range of Africanized honey bees, the northern range limit may covary to some degree with fluctuations in winter temperatures. Southwick et al. [37] used physiological measurements of energy balance to predict Africanized honey bees would be limited by the isocline of 120 consecutive days not exceeding 10°C during winter. Under this prediction, Africanized honey bees were expected to extend their range at least to the border with Canada along the mild west coast of the United States. This has clearly not yet occurred.
Our data on further range extension should be treated with caution because we surveyed only summer foraging workers and therefore cannot determine whether Africanized bees were successfully overwintering in all areas in which they were detected. Studies of overwintering of feral hives are needed to determine this. In addition, our data on the frequency of Africanized honey bees needs to be tempered by the fact that our state-wide survey assayed only freely foraging workers, not swarms or feral hives. Although we assayed workers after nearly all crops had finished flowering and many imported commercial honey bee hives had presumably left the state, domestic hives were still frequently observed in agricultural regions such as the Central Valley. Therefore we cannot assess the frequency with which feral bees are Africanized in our statewide sample, only that among workers encountered while foraging, Africanized bees are currently the clear minority north of approximately 35.5 degrees of latitude.
North of San Diego County honey bees carrying the African mitotype were smaller on average than those carrying European mitotypes. This difference is expected for early stages of introgression, though such associations may subsequently disappear as linkage between mitotype and nuclear gene content breaks down with further introgression [5,38]. In addition, worker bees in the statewide sample were on average larger than those collected in San Diego (Fig 2a and 2b). This may be due partly to the more recent arrival and lower frequency of Africanized honey bees in northern regions, but another factor is that feral European honey bees in California exhibited Bergmann's Rule (larger body size with higher latitude) prior to the invasion of Africanized honey bees [39].
In San Diego County, the African mitotype was found in nearly 65% of foraging workers and 61% of workers were classified as Africanized by morphological analysis, but there was no association between body size and mitotype. Assuming that genetic determinants of body size are at least partly controlled by nuclear genes, it appears that the process of Africanization in San Diego has likely followed the pattern previously documented in Texas [5] where initial linkage between mitotype and nuclear gene content broke down over a period of less than decade following invasion. What remained was a feral population of bidirectionally hybridized bees in which the frequencies of European nuclear genes and mitotypes averaged 25-37% with no linkage disequilibrium between genes from maternal and biparental genomes (Pinto et al. 2005). Our observed frequencies of European mitotypes and morphologically European bees are at, or slightly in excess of, the higher limit of the range documented in Texas [5]. However, in that study [5] only bees from feral colonies were studied while our collecting strategy focused on workers that may have come from either domestic or feral sources, potentially explaining the slightly higher frequency of European mitotypes and morphologies in our sample. Assays of nuclear gene content are needed, but it seems likely that rather than considering > 60% of honey bees in San Diego as Africanized, the actual situation is that most worker bees in the county represent a hybrid population containing > 60% African genes.
Our sample indicates that the African mitotype remains uncommon in domesticated hives in San Diego County, except for the few cases where hives were founded from feral swarms. By contrast the frequency of the African mitotype in the ten feral hives sampled (70%) was not different than that found among freely foraging workers in San Diego County. Because 21 of 24 (87.5%) domestic hives but, only 35% of foraging workers, carried the European mitotype, it appears that the great majority of foraging worker honey bees encountered in San Diego County do not come from domestic hives, but rather represent a large feral honey bee population. In addition, we observed heterogeneity of mitotypes among the three workers collected at each site as often as would be expected by randomly sampling the entire population. Therefore patches of flowers are often visited by workers from more than one hive. Honey bees can forage as much as 10 kilometers from their hive [40], but in our sampling, the floral resource from which bees were collected were common on the landscape, consisting of either frequently planted horticultural varieties in urban landscapes or common native species in rural ones. When floral resources are common, median foraging distances of honey bees may be considerably shorter than the maximum. For instance, a median foraging distance of 1.1 kilometers was estimated by analysis of waggle dances in a study area in Germany [40]. In a study of foraging distances of honey bee from hives on the UC San Diego campus, the median distance was 0.929 kilometers (Park, B. and Nieh, J. unpublished data). If foraging distances are similar throughout San Diego County, the distances between feral hives may often be less than 2 kilometers.
Our most puzzling data come from measurements of workers from feral hives. Only three of ten colonies were scored as Africanized by the discriminant function method [25], and bees from feral hives did not differ morphologically from bees from managed hives even though the managed colonies much more frequently carried European mitotypes. Our best explanation is that this result may be due simply to the small sample size (N = 10), although even with this small sample, finding that only 3 (30%) colonies that scored as Africanized rather than the expected 60.85%, as found in our large sample of San Diego County foraging workers is marginally significant (χ 2 = 3.995, P = 0.046).
Sequencing of the COI-COII spacer region confirmed that the cytochrome b-BglII restriction method reliably separates the "A" (African) mitotype from several mitochondrial lineages that are found in Europe and the Middle East. It should be noted that Franck et al. [32] found that the A mitotype occurs north of the Sahara in western North Africa, Sicily and the Iberian peninsula. Therefore, honey bees with ancestry from these sources could be detected as carrying the African mitochondrial lineage and might account for the occasional report of "African" mitochondria that pre-date the recognized arrival of Africanized honey bees from Central and South America. For instance, in a survey of 451 feral honey bee colonies collected across the southern United States (excluding California) prior to the recognized arrival of Africanized honey bees, one colony, collected in 1980 from Georgia carried the African mitotype [25,28]. This finding suggests some previous introduction, perhaps through port traffic. If such prior introductions succeeded in California, their offspring appear to be rare or absent in the northern parts of California and the pattern of occurrence of the "A" (African) mitotype strongly suggests little if any presence in California prior to the invasion of Africanized bees in 1994.
Our investigation of the DNA barcode region of the mitochondrial COI gene provides a new method to distinguish between African and European mitotypes. This region contains many other polymorphisms, and could potentially be used to distinguish further among African and European subtypes (see e.g. [8,32]). However the most important utility of this diagnostic SNP will be for further documentation of the spread and frequency of Africanized honey bees. DNA barcoding is becoming an increasingly popular method of assessing biodiversity used both in research and classroom contexts. Because honey bees are frequently collected during invertebrate biodiversity assessments that utilize DNA barcoding, the Barcode of Life Database will provide ongoing data on further spread of Africanized honey bees.
The pattern of Africanization we have documented in San Diego County and elsewhere in California appears consistent with patterns previously documented in Texas where Africanized honey bees first appeared in the United States [13]. In San Diego County, Africanized honey bees were first reported in 1994 and there now appears to be an admixture of European and African genes with the African component in the majority and no linkage between mitochondrial and morphological assessments of Africanization. Of great interest for further study would be complete assessment of nuclear genomic origins using molecular tools, [41], particularly whether there is a pattern as to which European and African nuclear genes persist after hybridization stabilizes across the geographic range of Africanized honey bee expansion. Secondly, we have not assessed aggressive behavior in Africanized honey bees in California. Doing so requires assessment of colony-level aggression as it is at this level that Africanized and European honey bees strongly differ behaviorally. Reports from the island of Puerto Rico suggest that Africanized bees there are less aggressive than elsewhere, perhaps due to limited penetration into the Africanized population of certain regions of the African nuclear genome [42,43]. Behavioral studies of Africanized honey bees at the northern limits of their continental range are lacking. It is thought that aggression or its behavioral correlates is a part of the ecological advantages that allow Africanized honey bees to displace European honey bees [4,6]. If so, we expect genes that control this behavior to persist in hybridized populations wherever aggression is advantageous. Conversely, whether certain European genes provide advantages to hybrids over pure African genomes remains to be studied.
Supporting Information S1 Table. Honey bee collection sites north of San Diego County (S1a) and within San Diego County (S1b). The number of honey bees, date collected, number carrying African mitochondria, location and elevation of collecting sites north of San Diego County. (DOCX) S2 Table. Morphometric measurements, discriminant function scores and outcomes, and cytochrome b restriction results for honey bees collected north of San Diego County (S2a) and within San Diego County (S2b): FWL, fore wing length; HWL, hind wing lenghth; TL, tibia length; FL, femur length (all in mm). (DOCX) S3 Table. Morphometric measurements, discriminant function scores and outcomes, and cytochrome b restriction results for honey bees collected from hives of 11 San Diego beekeepers (S3a) and from 10 feral hives within San Diego County (S3b): FWL, fore wing length; HWL, hind wing lenghth; TL, tibia length; FL, femur length (all in mm). (DOCX) S4 Table. Genbank accession numbers for honey bee mitochondrial COI-COII spacer region sequences used in S1