Multiple mitogenomes indicate Things Fall Apart with Out of Africa or Asia hypotheses for the phylogeographic evolution of Honey Bees (Apis mellifera)

I have just been thinking, and I have come to a very important decision. These are the wrong sort of bees.”—Winnie the Pooh.

A mitogenomic phylogeography of species of Apis and subspecies of Apis mellifera

The mitogenomic analysis confirms with high confidence previous inferences as to morpho-behavioral relationships among nine species as Dwarf, Giant, and Cavity-Nesting Honey Bees (Figs. 1a,b). Within the Apini, rooting with Bombus (Bombini) confirms Dwarf Honey Bees [A. florea + A. andreniformis] as outgroup to the Giant Honey Bees [A. dorsata + A. laboriosa] and the remaining, Cavity-Nesting Honey Bees. A. mellifera is the sister group to these [A. koschevnikovi + [A. nigrocincta + [A. nuluensis + A. cerana]]], at a greater phenetic distances than among any other Apis species pairs.

Mitogenomic sequences assigned to the taxon A. m. mellifera in GenBank occur in nine lineages collectively made paraphyletic by the placement of other subspecies. The phylogenetic analysis indicates that most of these apparent anomalies are artefacts of miss-assignment of individual sequences to the type subspecies. As noted, the Kenyan bee series when correctly referred to subspecies accords with the cladistic arrangement described here. The Arabian bee series referred to A. m. mellifera in GenBank falls phylogenetically into two distinct clades, one as part of the basal A. m. mellifera sensu lato and the other A. m. jemenitica sensu lato, along with Kenyan bees reassigned to the same subspecies. Another two Kenyan sequences assigned to A. m. mellifera are more closely related to A. m. simensis and A. m. unicolor, respectively, which suggests they are in fact members of those two subspecies. The latter extends the subspecies range beyond insular Malagasy.

A revised phylogeographic hypothesis for the evolution of A. mellifera

The received picture of the evolution of A. mellifera has been an “Out of Africa” model9, with as many as three Recent excursions, originally to Europe, more recently to Iberia, and in historic times to South America as so-called “Killer” or “Africanized” Bees, which resulted from accidental crossing of imported A. m. scutellata queen with local A. m. ligustica drones7. This has been challenged more recently by models of “Out of Asia”10 or North African/Middle Eastern12 origins. The alternative phylogeographic clade structure presented here indicates that Things Fall Apart with any of these hypotheses. The mitogenomic phylogeography (Fig. 3) indicates instead that A. mellifera evolved as a North-to-South expansion from Europe to Africa by way of Asia Minor and the Levantine. The type form A. m. mellifera sensu stricto originated in northern Europe, and has diversified as A. m. ligustica in southeastern Europe, including expansion into Asia Minor (A. m. caucasia). European bees then spread southward via the Levant (A. m. syriaca) into Nilotic East Africa (A. m. lamarckii), and across the Red Sea into the Yemeni coast of the Arabian Peninsula (A. m. jemenitica). Southward expansion has left Ethiopian (A. m. simensis) and Malagasy (A. m. unicolor) lineages as earlier offshoots. Sub-Saharan mainland forms constitute a single clade that comprises individuals referred to A. m. scutellata and A. m. capensis, as well as A. m. adansonii and A. m. monticola both closely related to sequences assigned to the former and latter subspecies, respectively. A. m. adansonii occurs throughout central Africa: the single available sequence from Niger may not be representative. The structure of the Mediterranean clade indicates a secondary return to Europe during a refugial period, and a more recent tertiary return to North Africa via the western Mediterranean islands.

Inferences from an mtDNA-based molecular clock for the evolution of A. mellifera

Divergence times within and among animal species, including Insecta, are routinely estimated from measured molecular divergences, including those made from mtDNA genomes29, including Apis30. Calibration of a molecular clock requires reliably dated external events, often geographic31, and accurate measurements of rates of substitution or divergence (2 × the former). Reliable geographic events are unavailable for A. mellifera30; substitutions per bp measured here over complete coding regions avoid gene-specific variance. Based on a fixed rate of 0.0115 substitutions/site/Myr (Brower 1994 in Ref.31) diversification of major lineages within A. mellifera can be dated to the Chibanian Age (late Pleistocene), 770–126 Kya. This coincides with and may be influenced by the European Günz glacial cycles. Separation of more southerly lineages from continental European subspecies occurred ca. 720 Kya, separation of African-endemic from Levantine/Nilotic/Arabian lineages ca. 660 Kya, establishment in Africa ca. 540 Kya, and separation of Euro-African Mediterranean from Sub-Saharan lineages ca. 250 Kya. Diversity within several subspecies dates to > 100 Kya, notably among multiple replicates of A. m. ligustica and A. m. jemenitica, 120 Ka and 150 Ka, respectively. Differences among other nominal subspecies originate only a few 10 s of Kya, notably the Asia Minor forms (A. m. meda, A. m. caucasia, and A. m. anatoliaca, 10–20 Ka) and Mediterranean species including insular forms (A. m. ruttneri and A. m. siciliana, 40 Kya).

Previous inferences from morphology

The standard classificatory and evolutionary system of Apis is based on meristic analysis of four continental groups, African (A), European (C), Mellifera (M), and Asian (O)8. In Fig. 5 (redrawn from Fig. 10.8 in Ref.8), the placement of individual bees in Principal Components space corresponds roughly to the four quadrants of the first two axes, clockwise from the upper left as A, M, C, and O, respectively. Re-tagged and grouped genetically, the Arabian/Nilotic clade overlies the Sub-Saharan clade, including local variants. However, the Asia Minor A. m. caucasia clade is distributed along the M → O axis, and is bisected by its Southeast Europe sister clade A. m. ligustica. Although the African and European/Asian clades are essentially non-overlapping, the Mediterranean clade overlies all three, as well as the basal A. m. mellifera clade. Despite their dispersion in PC space, A. m. caucasia and A. m. iberiensis are the least genetically diverse clades; individuals in the latter plot are almost entirely African. Consolidation of all “African” subspecies as a single group does not recognize the origins of the Euro-African-Mediterranean subspecies, nor does a single “Asian” group recognize the complex biogeographic connection of Eurasian subspecies in southeast Europe and Asia Minor, or of these with Africa via the Levant.

Figure 5
figure 5

Molecular clades of A. mellifera mapped onto the first two PCA axes of morphometric ACMO space, redrawn after Fig. 10.8 in Ref.8. Individual bees are identified to subspecies by the numeric color codes, and re-grouped to genetic clades as in Fig. 3, with additional color variants for Sub-Saharan taxa. Ruttner’s four groups [African (A), European (C [Continental]), Mellifera (M), and Asian (O [Oriental]] correspond roughly to the four quadrants, clockwise from the upper left as A, M, C, and O. The Levantine/Nilotic/Arabian clade overlays the Sub-Saharan clade , including local variants. The Asia Minor A. m. caucasia clade is distributed along the M → O axis, and is bisected by its Southeast European sister clade A. m. ligustica. The Afro-European Mediterranean clade overlies all three, as well as the basal A. m. mellifera clade.

The ACMO model is taken as support for an argument against an Oligocene/early Miocene dispersal from Europe to Africa. More recent morphological analysis of wing-venation patterns including fossil material32,33 suggests a European origin of A. mellifera antecedents from an A. cerana/dorsata type in the Miocene, in agreement with the molecular clock estimate here. However, the hypothesized evolution of the A. mellifera type in Sub-Saharan Africa late in the epoch is according to the clock too far back, and return to Europe in the Holocene is both too far forward and in the opposite direction. Despite European antecedents, this remains a “(Long) Out of Africa” model. The similarity dendrogram in Fig. 1 of Ref.32 shows [Bombus + [[A. dorsata + [A. florea + [A. cerana / A. mellifera]]]]] with cerana & mellifera forms co-mingled. This arrangement transposes A. florea and A. dorsata and therefore the Dwarf and Giant Honey Bee morphotypes (cf. Fig. 1a), and fails to resolve the deep phylogenetic separation of A. cerana / A. mellifera (cf. Fig. 1b and Supplementary Fig. S4).

Previous inferences from whole mitogenomes

Previous meta-analyses by Boardman and Eimanifar et al.11,16,17,18,19,20,21,22,25,26 and Tihelka et al.12 are compared with the present analysis in Fig. 6. The base phylogram is a Maximum Parsimony analysis calculated as in Fig. 2, with the addition of two problematic sequences noted in the text, KY926882 and KY926883, attributed to A. m. syriaca and A. m. intermissa, respectively. KY926882 “” falls outside the jemenitica cluster that includes A. m. syriaca KP163643, which is not used in either the “” or “” sets. KY936883 “” is more closely related to KY936882 than to A. m. intermissa KM458618, which occurs where predicted biogeographically. There are indications of faulty data in KY936882 and KY936883: inspection of the last quarter of the aligned sequences shows two unbroken runs of seven shared phylogenetically informative sites31 that unite KY936882 and KY936883 to the exclusion of their nominal subspecific sisters. In the Sub-Saharan clade, the A. m. capensis/A. m. scutellata pair (KX870183/KY614238) in the “” set occurs in A. m. capensis + A. m. scutellata clade I (but not II), and in the “” set MG552698 appears to be mis-referred to A. m. scutellata instead of capensis. A. m. adansonii (MN585109) “” is more closely related to the other, larger clade of A. m. scutellata. A. m. monticola (MF678581) “” is more closely related to A. m. capensis. Thus, reliance on single, atypical sequences in previous meta-analyses has led to faulty phylogenetic inferences.

Figure 6
figure 6

Boardman, Eimanifar et al. [11 et seq.] root their species-inclusive trees at the midpoint, rather than by an outgroup, which as drawn implies that the Dwarf and Giant Honey Bee species are sister groups, rather than that the former are the outgroup to the remaining species (cf. Fig. 1b). Their alignments include rDNA along with Coding Region sequences, which may contribute to the relatively weaker bootstrap support for some nodes: among subspecies of A. mellifera, my trial rDNA alignments were problematic at best, and these regions were excluded from the analysis here.

Tihelka et al.12 inferred their phylogeny by a BI site-heterogeneous mixture model. The root with respect to A. cerana falls between A. m. intermissa and other African (A) subspecies (notably A. m. scutellata), in support of an African origin, in contrast to Figs. 2 and 3. A. m. mellifera and A. m. iberiensis (M) are on the same side of the root as A. m. intermissa. The Asian (O) subspecies are further along the backbone, followed by the European (C) subspecies. The Asian form indicated by “?” is unidentified: the likeliest candidate seems to be A. m. caucasia. Of 11 nodes in their tree, only the C group is found in Fig. 2 here, and none corresponding to groups A, M, O, or any of the three-member clades. Exclusion of third-positions SNPs in their P12 model yields essentially the inverse of the trees in Fig. 2, with African lineages first derived, and A. m. mellifera as outgroup to the Eurasian lineages. Their model transposes the African and Eurasian lineages to either side of A. m. mellifera. At least two nearest-neighbor subspecies pairs are markedly different: East African A. m. simensis/A. m. unicolor, A. m. intermissa with respect to A. m. sahariensis/A. m. iberiensis. Relationships among subspecies lineages are similar to those in Fig. 2 (equivalent to the P123 model), but with much lower bootstrap support. However, within the A. m. simensis-inclusive clade, structure and support for the African and Mediterranean lineages collapse in the P12 with respect to the P123 model, and bootstrap support for key branches is < 50%. Especially, MN119925 A. m. unicolor is transposed as a long branch from outside to inside these clades. Of course, not all first-position SNPs are substitutions, nor all third-position SNPs silent. A Maximum Parsimony analysis of inferred amino acid substitutions at 125 sites across representatives of 22 subspecies (Supplementary Fig. S5) gives a clade structure is similar to Fig. 2 here or that from their P12 model (except as to rooting), with bootstrap support for African and Mediterranean lineages again < 50%. Bootstrap support for other clades in Fig. 2 is strong.

As shown by the molecular clock in Fig. 4, pairwise differences among nominal subspecies within the Southeast European, Asia Minor, and Mediterranean clades are less than those between individuals within other subspecies, even where within-clade relationships are well-defined. Where these local subspecies were originally defined by perceived morphological differences, review of their alpha and beta taxonomy may be indicated, with a view towards synonymizing local forms in gamma taxonomy: names matter27.

Additional phylogeographic inferences from partial mitogenome sequences

Besides the 22 subspecies of A. mellifera with complete mitogenomes compared in the main MS, GenBank includes six additional subspecies with partial sequences from the 5′ end of the ND2 region. Their names and provenance are A. m. adami (Crete), A. m. cecropia (Greece), A. m. cypria (Cyprus), A. m. macedonica (Greece into Romania), A. m. pomonella (Kazakhstan), and A. m. sicula (Sicily). Over the first 574 bp common to all these subspecies, in combination with the 22 subspecies’ mitogenomes examined in the main MS, there are 39 variable sites of which 19 are phylogenetically informative sensu Nei34 (Supplementary Fig. S6).

Sequences from subspecies A. m. adami, A. m. cecropia, and A. m. macedonica are identical to those from A. m. ligustica and carnica in the Southeast Europe clade, as is the sequence from A. m. cypria except for a single autapomorphic Y/G SNP. The sequence from A. m. pomonella is identical to that of A. m. meda in the Asia Minor clade. The sequence from A. m. sicula is identical with that of the reference A. m. mellifera, except for a T SNP otherwise found only in the African superclade.

Inclusion of these short sequences tests the phylogeographic hypothesis. The similarity of the Central Asian A. m. pomonella and Middle Eastern A. m. meda extends the Asia Minor clade further eastward, separately from A. m. sinisxinyuan. A. m. sicula is cladisitically distinct from A. m. siciliana in the Mediterranean clade, and the co-occurrence of the two forms on Sicily suggests a secondary colonization from the westerly islands after in initial occupation from the European mainland. Likewise, the continental association of A. m. adami and A. m. cypria from Crete and Cyprus, respectively, is consistent with a dispersal corridor from southeastern Europe through the eastern Mediterranean and into the Levant, and suggests an eastern limit of the re-colonization of the insular Mediterranean at Sicily. The phylogeography proposed here accommodates these data well: further insight may be expected from complete mitogenomes.

Previous inferences from nucDNA genome sequences

The “Thrice Out of Africa” hypothesis9 follows a phylogenetic scheme that retrieved the same four groups (ACMO) as the then current continental geographic scheme based on meristics and morphology8. Their molecular network transposes the alphabetical order to an MAOC backbone, rooted with respect to a composite outgroup within the (A)frican subspecies, which along with Mellifera (M) was separated from the European (Continental) and Asian (Oriental) subspecies. Subsequent analyses of various nucDNA markers and combinations of A. mellifera subspecies broadly identify the MAOC groups9,10,23,24,32,33,35, however placement of the root on various branches according to internal and (or) external criteria varies, so as to suggest alternative geographic origins of A. mellifera. These include alternative interpretations of the Thrice Out of Africa data29 [root indeterminate between M + A and O + C], as well as Thrice Out of Asia9 and North African/Middle Eastern origin models10,32,33,35, for example with A. m. jemenitica in the Y lineage as outgroup to other subspecies of A.

There are fundamental differences between the MAOC backbone and mtDNA genome phylogenies (Supplementary Fig. S7), both with respect to placement of the root and allocation of subspecies to clades. No MAOC model places A. m. mellifera as outgroup to other subspecies: the mtDNA data place A. m. scutellata distal rather than proximal to the root of A. mellifera evolution. Closely related subspecies and geographically contiguous subspecies within the Mediterranean mtDNA clade here, including A. m. iberiensis and A. m. intermissa9,32, and A. m. ruttneri24, are dispersed over the M, A, and C clusters, respectively. The A. m. jemenitica/A. m. lamarckii/A. m. syriaca clade here is also dispersed over the A(Y), A(L), and O clusters. On the other hand, pairings of A. m. ligustica and A. m. carnica (C) and A. m. anatoliaca and A. m. caucasia (O) are consistent between nuc- and mtDNA data. I note that allocation of subspecies among rooted MAOC clusters conforms closely to the original continental model8, which was not cladistic in approach or form. The mtDNA-based model explains a more complex continental distribution, where European, Asian, and African subspecies are of multiple origin.

Classification of nucDNA SNP markers from eight of the subspecies examined here24 inter alia overlays A. m. anatoliaca and A. m. caucasia in the O group, and separates A. m. iberiensis from A. m. mellifera (three clusters) in the M lineage. A. m. ruttneri, the only representative of the African A group, is placed midway between A. m. mellifera and the C group, rather than close to A. m. intermissa and A. m. iberiensis in the Mediterranean clade as here. That is, lineages paired in the two subspecies of the Eurasian clade and variants within the Mediterranean clade here are separated by SNP data in contrast to the clade hierarchy here.

Contrasts between nucDNA, mtDNA, and even morphological phylogenies are not unknown29. Maternally-inherited mtDNA phylogenies have the virtue of tracing maternal evolutionary lineages36, and may thus be particularly reliable for phylogeographic inferences about the origins and spread of “queen”-dispersed eusocial insect. A phylogeographic origin Out of Europe inferred from prior affinity with central Asian cavity-nesting bees seems more parsimonious that disjunct Sub-Saharan or Asian origins.

Read more here: Source link