Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Genes Suggest Ancestral Colour Polymorphisms Are Shared across Morphologically Cryptic Species in Arctic Bumblebees

Abstract

Our grasp of biodiversity is fine-tuned through the process of revisionary taxonomy. If species do exist in nature and can be discovered with available techniques, then we expect these revisions to converge on broadly shared interpretations of species. But for the primarily arctic bumblebees of the subgenus Alpinobombus of the genus Bombus, revisions by some of the most experienced specialists are unusual for bumblebees in that they have all reached different conclusions on the number of species present. Recent revisions based on skeletal morphology have concluded that there are from four to six species, while variation in colour pattern of the hair raised questions as to whether at least seven species might be present. Even more species are supported if we accept the recent move away from viewing species as morphotypes to viewing them instead as evolutionarily independent lineages (EILs) using data from genes. EILs are recognised here in practice from the gene coalescents that provide direct evidence for their evolutionary independence. We show from fitting both general mixed Yule/coalescent (GMYC) models and Poisson-tree-process (PTP) models to data for the mitochondrial COI gene that there is support for nine species in the subgenus Alpinobombus. Examination of the more slowly evolving nuclear PEPCK gene shows further support for a previously unrecognised taxon as a new species in northwestern North America. The three pairs of the most morphologically similar sister species are separated allopatrically and prevented from interbreeding by oceans. We also find that most of the species show multiple shared colour patterns, giving the appearance of mimicry among parts of the different species. However, reconstructing ancestral colour-pattern states shows that speciation is likely to have cut across widespread ancestral polymorphisms, without or largely without convergence. In the particular case of Alpinobombus, morphological, colour-pattern, and genetic groups show little agreement, which may help to explain the lack of agreement among previous taxonomic revisions.

Introduction

Recently authors have sought to reduce subjectivity in the practice of species discovery or ‘delimitation’ [1] by comparing the results from applying multiple criteria in an explicitly integrative taxonomy [2]. This is a search for congruent results from multiple criteria as corroboration, even though the different sets of characters employed (from different disciplines) may have different evolutionary constraints requiring different evolutionary explanations [2]. Indeed, not all studies using different criteria can be expected to agree, because relationships between criteria and speciation events differ [1, 3, 4]. It has even been acknowledged that those criteria that relate most closely to speciation could form a minority pattern in some cases [2]. But for some arctic bumblebees, taxonomic revisions by experienced specialists have shown unusually poor agreement on the number of species, we suggest in part because of poor agreement among results from applying different criteria.

Bumblebees (genus Bombus Latreille) have been described as morphologically homogeneous compared with other bees [5]. It is therefore unsurprising that bumblebee species were diagnosed initially in terms of the striking colour patterns of the hair on the dorsum of the body [6]. But then, more than a century ago, it was realised that the groups supported by skeletal morphology are actually more clearly circumscribed, and that within these morphological species, bumblebees are highly variable in colour patterns [710]. This is supported by recent genetic studies [1114], reaffirming the long-held view that bumblebee colour patterns do not necessarily diagnose species [15]. Detailed studies of bumblebee variation have led to the realisation that species can be cryptic in both morphology and colour pattern [16, 17], now confirmed by genetic studies [18, 19].

Particularly intriguing is the observation that bumblebee colour patterns are often impressively similar among species [8, 2026]. In many cases, a single species may show different colour patterns in different geographical regions, in each of which they may resemble closely other species as members of regional colour-pattern groups [24, 27]. Several possible explanations have been proposed (reviewed in [26, 28]), including Müllerian mimicry [22, 29]. Recent advances in understanding evolutionary relationships have begun to clarify the chronology of the evolution of some of these resemblances [30], an important pre-requisite to discerning causes [31]. In cases of regional colour-pattern resemblance, the similar species are often only distantly related, so that the similarity has been interpreted as the result of evolutionary convergence [26]. We explore another possible explanation for the evolution of colour-pattern resemblance, as an alternative to convergence: that in some situations a polymorphism within an ancestral species may have been inherited in parallel by several descendent species.

The bumblebees of the subgenus Alpinobombus Skorikov are a small group of closely related species [32]. Previous taxonomic revisions (including global revisions and broader regional faunal revisions with global referencing) of Alpinobombus by some of the most experienced bumblebee specialists have relied on both morphology (especially the more variable characters of the male penis valve, volsella and gonostylus, and of the female oculo-malar area and hind tibia) and on the colour patterns of the hair and have all differed in their conclusions on the number of species present (Table 1). Several recent faunal lists have accepted as valid just the five more clearly morphologically-diagnosable species (interpreted in a broad sense, from the shape of the female head, male genitalia, etc.): B. alpinus, B. polaris, B. balteatus, B. neoboreus, and B. hyperboreus [10, 3335]. But variation in colour pattern prompted some authors to question whether at least two more species might be added to this list: B. pyrrhopygus (but using the junior synonyms B. diabolicus Friese and B. alpiniformis Richards for bees with the Scandinavian colour pattern of this species) and B. kirbiellus, even though both lacked clearly discrete diagnostic morphological characters [3638]. The colour patterns of Alpinobombus species have been described as being both variable within species while showing close similarities among species [36, 37, 3943]. These species are unusual even among bumblebees [10, 44] for occupying extreme arctic or alpine habitats [40, 43, 45]. Therefore, at least in the far north, they co-occur with few or no other bumblebee species [36, 37, 43, 46], which simplifies the system for analysis if mimicry were involved.

thumbnail
Table 1. Lists and numbers of Species of the Subgenus Alpinobombus from previous taxonomic revisions.a

https://doi.org/10.1371/journal.pone.0144544.t001

We use sequences for parts of the COI and PEPCK genes from freshly collected specimens and from museum specimens to recognise the evolutionarily independent lineages of the subgenus Alpinobombus as species and to estimate the phylogenetic relationships among these species. We then test whether the colour-pattern polymorphisms are more likely either to be ancestral, or to have arisen more recently and convergently within the daughter species. From our analyses we conclude: (1) that coalescents for the COI gene support nearly twice as many species (nine cf. five) as have been recognised from skeletal morphology and that three of these species, including one taxon not recognised previously, are also supported by PEPCK polymorphisms; and (2) that similarity in colour pattern among parts of these species is likely to have arisen from broadly conserved ancestral polymorphisms, not from convergence. The work reported here builds on a survey of the morphology and colour patterns of the North American Alpinobombus [42] and is part of a revisionary study of Alpinobombus bumblebees world-wide, which will treat type specimens and nomenclatural issues separately.

Materials and Methods

Ethics statement

No humans or vertebrate animals were used in this research. Field permits: (Russia) Bumblebees of the subgenus Alpinobombus are not included in any Red Lists and are not listed under the law ‘On environmental protection’ in Russia. None of the areas where field work was carried out are private. Only Meduza Bay in West Taimyr is a part of the Great Arctic Reserve, where the research was conducted with the permission of the Reserve management. Specimens from the research collection of Novosibirsk passed the State veterinary control before export on loan.

(Norway) No permits were needed for collecting of the Norwegian data.

(Sweden) None of the bees were collected within National Parks or Nature Reserves. They are not subject to any collection restrictions and therefore permits are not required. At present there is no legislation or restrictions that regulates the export of genetic samples from Sweden. There is no other restriction on the transfer of scientific material of this sort between EU member states.

(Denmark / Greenland) The Greenland Ecosystem Monitoring Coordination Group at the National Environmental Research Institute, Aarhus University, approved our research proposals for access and research activities in both 2010 and 2011.

(Canada) No permits were required for collection sites visited.

(USA) Collections from the Rocky Mountain Biological Laboratory in Colorado were made according to an agreement between the RMBL and the United States Forest Service.

We bring together samples of taxa of the subgenus Alpinobombus from throughout their global distributions. Taxonomic and geographical breadth is essential in a taxonomic revision for achieving a representative sample in order to assess the full range of patterns of variation and relationship [47, 48]. This is especially important for reducing the bias that might otherwise result if sampling were restricted to just one or two isolated regions or were to exclude major lineages. To cover all variation, we made a preliminary survey of specimens from museum collections (Fig 1; see the supporting information).

thumbnail
Fig 1. Global Sampling of the Subgenus Alpinobombus.

Grey spots show sample sites for the 4345 specimens examined (Alpinobombus species are unknown from the tropics or from the southern hemisphere; collections listed in Table A in S1 File); black spots show the sites for the 173 samples with COI barcodes; and question marks show the sites for doubtful records (specimens supposedly from sites that combine low latitude with low elevation). Polar projection, North Pole (starred) at the centre of the map. Terrestrial national boundaries and the Arctic Circle are shown as narrow grey lines. Image created in ArcGIS using World_Shaded_Relief basemap which is Copyright 2014 ESRI.

https://doi.org/10.1371/journal.pone.0144544.g001

To recognize species of the subgenus Alpinobombus, we accept the broadly held concept of species as evolutionarily independent lineages (EILs) [1]. To recognize these species in practice, we use statistical tests designed to identify gene coalescents, which provide direct evidence for EILs [49]. Coalescence methods are suggested to be especially useful compared to less direct indicators (including criteria related to mate recognition and reproductive isolation) because they can detect lineage separation even in its early stages [1, 3]. We use one of the fastest evolving, most highly variable mitochondrial genes, coding for COI (cytochrome oxidase subunit I). This is particularly appropriate for detecting species because of its rapid coalescence within lineages [3, 50]. From among recently collected specimens showing the most obvious morphological and colour variation, we sequenced the ‘barcode’ region of the COI gene (Fig 1).

Most of the extraction, amplification, and sequencing work was done at the Biodiversity Institute of Ontario, Guelph, using standard protocols [51] and primers [52]. A minority of specimens were processed at the Institute of Apicultural Research, Beijing. Because the geographic range occupied by Alpinobombus is large and many areas are difficult to access, mitochondrial DNA with many copies per cell has the advantage of being more easily extracted, so that even pinned museum specimens can be sampled. Although 70% of Alpinobombus specimens surveyed are known to be less than a century old, only 8.6% are less than 10 years old. A problem with older museum specimens is that the sequences obtained are often not of full length [53]. Sequence alignments were made using the ClustalW function in MEGA (version 6.06, accessed 2014: megasoftware.net [54]) and amino acid translations were checked with EMBOSS Transeq (accessed 2014: ebi.ac.uk/Tools/st/emboss_transeq/).

To check for possible mismatches between trees for different genes, and hence for possible mismatches between gene trees and species trees [49, 55], we also sequenced the nuclear gene PEPCK (phosphoenolpyruvate carboxykinase; accession numbers for sequences listed in Table A in S2 File, obtained using primers from [32]). This gene has been popular for studies of bumblebees [11, 30, 32].

To recognise species, single-threshold general mixed Yule/coalescent (GMYC) models seek a transition expected in the tree-branching rates between interspecific branches and intraspecific branches [5661]. GMYC analysis requires an ultrametric estimate of the phylogenetic tree that includes only unique haplotypes [59, 62] (T. Barraclough, pers. comm.) in order to avoid spurious inflation of estimates of terminal branching rates. After sorting the sequences by decreasing length, the longest examples of each unique haplotype are recognised using Collapse (version 1.2, accessed 2011: softpedia.com/get/Science-CAD/Collapse.shtml), which ignores sites with missing data. As outgroups for rooting the tree with this fast gene, we include species from closely related bumblebee subgenera [32]: B. (Cullumanobombus) rufocinctus Cresson, B. (Pyrobombus) vagans Smith, B. (Bombus) ignitus Smith, B. (Bombus) terrestris (Linnaeus), and B. (Bombus) cryptarum (Fabricius). The best nucleotide-substitution model for COI according to the Bayesian information criterion (BIC) obtained from MEGA is the general time-reversible model with a gamma-frequency distribution of changes among sites (GTR+Γ). BEAST (version 1.8.0, accessed 2013: beast.bio.ed.ac.uk) was used for Bayesian analysis of multiple model trees [63]. The clock model was set to the uncorrelated lognormal (relaxed clock), the tree-speciation prior was set to a constant-size coalescent process (consistent with the null hypothesis that there is a single species in the data), and the chain length for the Markov-chain Monte Carlo (MCMC) algorithm was set to two billion generations, with sampling of the trees every 200,000 generations. The sample of resulting trees from the MCMC algorithm was examined using Tracer (version 1.6.0, accessed 2013 [63]), which showed that stationarity had been achieved within the first 1% of the total MCMC generations. The large number of MCMC generations was needed with our data to increase effective sample sizes. A maximum clade-credibility tree was obtained from the post burn-in tree sample after rejecting the first 1% of sampled trees using TreeAnnotator (version 1.8.0, as for BEAST). After removal of the outgroups, the GMYC model was applied to the tree using the SPLITS library (version 1.0–11, accessed 2011: r-forge.r-project.org/projects/splits/) running on the R platform (version 2.12.1, accessed 2011: www.r-project.org). The oldest available name [64] is applied to each coalescent group as a prospective species. Multiple runs of all analyses were made and the stability of the results checked.

To check the robustness of the GMYC results by taking into account the uncertainty in estimates of phylogeny, we applied a Bayesian GMYC analysis [59] using the bGMYC library (version 1.0.2, accessed 2015: sites.google.com/site/noahmreid/home/software) on the R platform. bGMYC was applied to the last 100 trees from the sample of 10,000 trees from the BEAST run of two billion MCMC generations (a sample representing the last 20 million trees) that had been checked previously with Tracer for stationarity. bGMYC was run with the default settings (MCMC 50,000 generations, burn-in 40,000, thinning 100).

As another check on the number of species recognised, we apply another related approach for discovering species designed for use with single-locus data, looking for changes in the numbers of nucleotide changes along branches from between- to within-species relationships, based on Poisson-tree processes (PTP: [65]). PTP analysis is less demanding of tree information than the GMYC approach and there is some evidence that it can perform better [65]. We apply the Bayesian implementation on the bPTP server (accessed 2015: species.h-its.org) to a metric maximum-clade-credibility tree from COI barcodes obtained with MrBayes (version 3.1.2, accessed 2011: [66]), using the same data, outgroups, and nucleotide-substitution model, four chains (temperature 0.2), and 100 million generations of the MCMC algorithm with a 1% burn-in selected using Tracer. We used the default bPTP options after removing outgroups from the rooted tree.

The variation in colour patterns for the subgenus Alpinobombus includes many subtle differences [42], which are similar in the two sexes. In all cases colour patterns refer to the colour of the hair on various body regions, not to the colour of the body surface, which is brown or nearly black throughout. Generally, we observe that the presence of non-black hair on the different body segments is not independent, but appears to be conditionally dependent among segments [8] (e.g. there are no Alpinobombus bees with the hair on metasomal tergum 2 black and on tergum 3 yellow, or on tergum 4 orange and on tergum 5 black, although in both cases the reverse conditions are common). Therefore it should be reasonable to begin by coding bumblebee colour patterns for the presence or absence of areas of pale hair and then later adding characters as modifiers for the extent of the pale hair [26, 67]. As a first approximation, we assume that Alpinobombus colour patterns can be reduced to a representation as two principal characters, each with two states. Character (1): whether there are transverse yellow bands in the hair on the body, usually anteriorly and often posteriorly on the thoracic dorsum and/or on metasomal terga 1–2 (state B, banded: minimally with an obvious yellow band on metasomal tergum 2, varying from pale straw yellow to darker orange brown), or whether the hair in these areas is black or nearly black (state U, unbanded). Character (2): whether there is obvious pale hair (most often red or orange, but sometimes yellow or white) on the ‘tail’, consisting of metasomal terga 4–6 (state P, pale: minimally covering tergum 5 or 6), or whether the hair of the tail is black or nearly black (state D, dark). The most widespread colour pattern both among all bumblebees [26] and among the sister-group to the subgenus Alpinobombus, the subgenus Bombus sensu stricto [19], has both the yellow bands and a pale tail. Replacement of the pale bands and the pale tail with black can also be polymorphic within some Bombus s. str. species, although this is relatively rare (occurring in 3/17 and 1/17 species for the two characters respectively). Colour characters were scored for a sample of 1646 Alpinobombus specimens of both sexes from across the distribution range. To reconstruct the history of colour-pattern characters, we estimate the species tree using BEAST linked trees from COI and PEPCK exon and intron sequences. The best (unlinked) nucleotide-substitution models according to the MEGA BIC are HKY for both the PEPCK exon and introns, the tree-speciation prior was set to a birth-death process appropriate for a multi-species tree, and from examination of the results with Tracer, the MCMC algorithm was set to run for 1.5 billion generations. Outgroups and other settings were as for the first analysis. No fossils of Alpinobombus species are known for dating parts of the tree, so the phylogeny was calibrated with a date from a molecular study: Hines [68] estimated the age of divergence between the subgenus Alpinobombus and the subgenus Bombus s. str. to be 13 Ma. Mesquite (version 3.02, accessed 2015: mesquiteproject.org [69]) was used to reconstruct ancestral characters states on this tree. We used parsimony analysis because we wanted to retain information on alternative explanations rather than just the most likely. We coded polymorphic populations as a third intermediate ordered state, allowing state changes in either direction ([70]; W. Maddison, pers. comm.).

Results

We found 46 unique haplotypes among 173 COI-barcode samples of the subgenus Alpinobombus (accession numbers listed in Table A in S2 File). These COI data are information rich, with 92/658 nucleotide sites phylogenetically informative. The GMYC models show a significant change in branching rate through time in the COI barcode tree (likelihood ratio 14.03 between the GMYC multiple-species model and the null model that there is a single species in the group, p = 0.0028). This threshold (at -0.0041 substitutions per nucleotide) leads us to recognise nine candidates for prospective species (with a 95% confidence interval of 8–11 species) within Alpinobombus as in Fig 2, which also shows that each of these species is strongly supported as a monophyletic group.

thumbnail
Fig 2. Recognising Species of the Subgenus Alpinobombus with GMYC.

The 173 sequences reduced to 46 unique COI-barcode haplotypes on a BEAST ultrametric gene tree (outgroups not shown). The single threshold (T = -0.0041) from the GMYC model is shown by the vertical grey bar and the intersecting lineages are interpreted as subtending nine prospective species (black spots show the coalescent node for each species). The tree is the Bayesian ultrametric maximum-clade-credibility tree from a sample of trees after 1% burn-in from 2 billion generations of the MCMC algorithm in BEAST. Values next to the nodes are Bayesian posterior probabilities showing branch support. Each haplotype is represented by one of the longest sample sequences, labelled with a species name, a code that consists of an identifier from the project database and (after the hyphen) from the BOLD database, followed by its geographic origin.

https://doi.org/10.1371/journal.pone.0144544.g002

The bGMYC analysis showed stationarity and a high modal coalescent/Yule ratio consistent with success. The posterior probability distribution is shown against the first sample tree in Fig 3 to provide a ‘heat’ map of the probabilities that the haplotypes are conspecific. Adopting a threshold probability of 0.5 from the posterior mean from the analysis (corresponding to a moderate position, midway between a ‘splitter’ and a ‘lumper’) recognises nine species within the subgenus Alpinobombus on the diagonal (unsurprisingly, the species with more haplotypes have more structure within them). These are the same nine species as are identified in the GMYC analysis in Fig 2.

thumbnail
Fig 3. Recognising Species of the Subgenus Alpinobombus with Bayesian GMYC.

The 173 sequences reduced to 46 unique COI-barcode haplotypes on one BEAST ultrametric gene tree (outgroups not shown). The posterior probability distribution (right, colour scale far right) is plotted against a sample tree from BEAST (left) to provide a ‘heat’ map of the probabilities that haplotypes are conspecific by bGMYC. Black spots show the coalescent node for each species from Fig 2.

https://doi.org/10.1371/journal.pone.0144544.g003

The Bayesian PTP solution with the most support (Fig 4) recognises the same nine prospective species (with a range of 7–30 species), although the individual Bayesian support values for the species are not high (0.48–0.79).

thumbnail
Fig 4. Recognising Species of the Subgenus Alpinobombus with Bayesian PTP.

The 46 unique COI-barcode haplotypes as in Fig 2 on a MrBayes metric gene tree with the Bayesian PTP solution with highest support, showing nine prospective species (black spots show the coalescent node for each species; outgroups not shown). The tree is the maximum-clade-credibility tree after 1% burn-in from 100 million generations of the MCMC algorithm in MrBayes. Values above the nodes are PTP Bayesian support values that all daughter haplotypes belong to a single species population; values below the nodes are Bayesian posterior probabilities showing branch support. Haplotype selection and labels as in Fig 2. The scale bar for branch lengths shows 0.02 substitutions per base position.

https://doi.org/10.1371/journal.pone.0144544.g004

As expected, we found much less variation in the PEPCK gene among species of the subgenus Alpinobombus. Only 4/903 nucleotide sites in our sequences are uniquely diagnostic for species, with unique PEPCK nucleotide changes supporting the species B. alpinus (base position 351), B. neoboreus (bp 755), and two for the unnamed species (bp 360, 533) (Table 2). Additional PEPCK nucleotide differences give further support for the distinction between the species pair B. neoboreus / unnamed (bp 769). The estimate of phylogeny for Alpinobombus species from the combined gene data is shown in Fig 5. This shows strong support for all groups except for the position of B. alpinus, which remains uncertain.

thumbnail
Fig 5. Estimate of phylogeny for species of the subgenus Alpinobombus.

Estimated using a linked-tree BEAST analysis of COI-barcode and PEPCK exon and intron sequences for each species. Values below the nodes are Bayesian posterior probabilities showing branch support. Numbers above the nodes show the estimated dates of divergence events in Ma (millions of years before the present) calibrated with a molecular estimate for the date of divergence between the subgenus Alpinobombus and the subgenus Bombus s. str. from Hines (2008) and grey bars show the 95% confidence limits on the estimated dates of divergence.

https://doi.org/10.1371/journal.pone.0144544.g005

The principal colour patterns recorded for the subgenus Alpinobombus show three of the four possible combinations of the two states for the two characters we survey (Table 3, Fig 6): UP, unbanded pale tails; BP, banded pale tails, and BD, banded dark tails. The fourth state, entirely black (UD), is not recorded here, although one specimen of B. balteatus from Karaginskiy Island (Kamchatka, Russia) comes close, with only a narrow and incomplete yellow band on metasomal tergum 2. Fig 7 shows the distributions of the three principal colour patterns. The UP colour pattern is the most restricted geographically and is concentrated primarily in Europe; UP also occurs in isolated pockets in Russia and in North America but is unrecorded from Greenland. In contrast, the BP and BD patterns co-occur broadly, both taxonomically (Table 3) and geographically (Fig 7). The frequency of the BD colour pattern relative to BP colour pattern may increase towards the north (Fig 7). Nonetheless, the three principal colour-pattern groups are each shared by three of the nine prospective species (Table 3; B. pyrrhopygus, B. polaris, and B. balteatus), as well as being recorded from both the Old World and New World regions, as widespread polymorphisms. While there are often small differences between the most frequent patterns for each species, nonetheless in all cases it is possible to find specimens that resemble other species so closely that diagnosis by colour pattern is not always possible.

thumbnail
Fig 6. Colour Patterns of the Subgenus Alpinobombus.

Simplified diagrams [26] representing the principal variation in colour patterns of the hair on the female dorsum coded as two two-state characters: UP, unbanded pale tail; BP, banded pale tail; BD, banded dark tail; with hair colours: yellow, yellow or yellow-brown hair; orange, orange-red or white hair; dark grey, black hair. Only female patterns are shown because females make up the majority of samples for social bumblebees, although male patterns are similar. Photos show: (left) B. alpinus and (right) B. hyperboreus. Photos by (left) A. Staverløkk and (right) G. Holmström.

https://doi.org/10.1371/journal.pone.0144544.g006

thumbnail
Fig 7. Distribution of Colour Patterns of the Subgenus Alpinobombus.

Pie diagrams (UP/BP/BD) show relative frequencies of colour patterns (black, UP; white, BP; grey, BD) from a sample of 1646 colour-coded specimens with coordinates by regions (clockwise from the top): Greenland (n = 27), Canada and Alaska (n = 257), southern Rockies (n = 25), Russia excluding Murmansk Province (Kola Peninsula) (n = 654), Altai and Sayan (n = 6), Scandinavia and Murmansk (n = 643), and Alps (n = 34). Spot diagrams show the occurrence of each colour pattern: grey circles show all sites represented in the sample, while black spots show those sites in which each colour pattern is recorded: unbanded pale tail (UP n = 329); yellow-banded pale tail (BP n = 1029); yellow-banded dark tail (BD n = 288). Map projection and other symbols as for Fig 1. Image created in ArcGIS using World_Shaded_Relief basemap which is Copyright 2014 Esri.

https://doi.org/10.1371/journal.pone.0144544.g007

thumbnail
Table 3. Principal Colour Patterns for each Species of the Subgenus Alpinobombus.a

https://doi.org/10.1371/journal.pone.0144544.t003

Reconstructions of the most parsimonious ancestral states for the two colour-pattern characters scored for the subgenus Alpinobombus are shown in Fig 8. For the banding pattern, Fig 8 supports either a root polymorphism with two subsequent reversals to a monomorphic banded pattern (to kirbiellus as well as to the ancestor of the group neoboreus-hyperboreus), or alternatively two forward changes to a polymorphic banded/unbanded pattern (to balteatus as well as to the ancestor of the group alpinus-polaris). For tail colour, Fig 8 supports a root ancestor with a polymorphism for tail colour. There is a reversal to a monomorphic pale tail for alpinus and a second forward change to a monomorphic black tail in the ancestor of the group natvigi-hyperboreus. Therefore we infer that both colour-pattern characters are likely to have had polymorphisms within the common ancestors of multiple species, which have then persisted in several of the descendent species as trans-species polymorphisms.

thumbnail
Fig 8. Evolution of Colour Patterns of the Subgenus Alpinobombus.

Reconstruction of ancestral states by parsimony in Mesquite for two colour-pattern characters each with two states (Fig 6) among species based on the estimate of phylogeny in Fig 5, several of which are polymorphic and show both states for both characters. Above: banding colour pattern, with yellow spots showing yellow-banded (B) populations, black spots showing unbanded (U) populations, mixed yellow/black spots showing polymorphic populations, mixed yellow/grey spots showing uncertain polymorphic/monomorphic populations (*for B. alpinus, males from the Alps often have a yellow-banded pattern although this is rare among females). Below: tail colour pattern, with orange spots showing populations with pale (P: orange or white) hair on the tail, black spots showing populations with dark (D: black) hair on the tail, mixed orange/black spots showing polymorphic populations.

https://doi.org/10.1371/journal.pone.0144544.g008

Discussion

Bumblebees have long been used in discussions of the nature of species [8, 21, 23]. Many of the different ideas of what constitutes species [71] have been applied to bumblebees, including the biological species concept [72, 73] and specific-mate recognition systems [67, 74]. Among these ideas, pattern-based and process-based views have sometimes appeared to be irreconcilable [10, 75]. But more recently a consensus has emerged in favour of thinking of species generally as EILs, which provides a framework that is sufficiently broad to accommodate the many different criteria used to recognise species in practice, including both pattern-based and process-based approaches [1].

The ability of the concept of species as EILs to accommodate multiple criteria has led to the idea of integrative taxonomy [2]. The integrative approach is based on the premise that the best way to discover species is to undertake multiple studies of the same set of individuals using different character sets and criteria, to compare the results, and then the best answer is found where different studies agree, i.e. through corroboration. Arguments that integration provides a ‘total evidence’ approach to bring out a shared underlying signal are less well founded, because the different character sets may have phylogenetic histories that are not shared [4]. The integrative approach has been applied previously for recognising bumblebee species [73, 76, 77]. However, our analysis of the subgenus Alpinobombus shows unusually weak agreement for bumblebees among groups from morphological, colour pattern, and genetic data. Morphology has been used in earlier studies to recognise five groups. Genetic data are used here to recognise nine groups, which although they are nested within the morphological groups, are mostly incongruent with them (Fig 9). Colour patterns can be used to recognise three principal groups, although most of the morphological or genetic groups overlap with more than one colour-pattern group (Fig 9). Within the framework of integrative taxonomy, this poor agreement among the groups supported by all three different character sets (the best that is achieved is nested sets: Fig 9) would require invoking special evolutionary explanations for each [2]. In this situation, the explanations might involve time lags between the evolution of some colour, genetic, and morphological characters [1] and this needs further investigation. This lack of congruency among groups based on different character sets may help to explain why previous taxonomic revisions have failed to agree on the number and identity of Alpinobombus species (Table 1).

thumbnail
Fig 9. Disagreement among Groups of the Subgenus Alpinobombus Based on Three Character Sets.

Above, set diagrams for the different but nested groups within Alpinobombus based on exoskeletal morphology (in grey) and based on genes (in black). Below, set diagrams for the different and mostly overlapping groups based on the colour patterns of the hair of the dorsum (in grey) and based on genes (in black). Areas of intersection of sets occupied by bumblebee samples are shown in grey. Abbreviations for colour patterns: BP, yellow-banded pale-tailed; UP, unbanded pale-tailed; and BD, yellow-banded dark-tailed.

https://doi.org/10.1371/journal.pone.0144544.g009

Can we find other character sets that might support groups showing better agreement? Another approach for bumblebees considers characterising chemical extracts of male labial gland secretions. These extracts include sex-attractant or arrestant pheromones [7780], which are expected to be important in specific-mate recognition systems [74]. These extracts can work well as indicators of species [7982]. However, there is also evidence that within some bumblebee species there can be both divergence in these male secretions across different parts of a single species’ distribution, as well as divergence in the female preferences among these secretions [83]. In addition, and unfortunately for the case of Alpinobombus, when closely related bumblebee species have geographical distributions that do not overlap, at least sometimes the mate-recognition pheromones do not differ between these species [73]. We would anticipate that pheromones could be undifferentiated between the three recent east-west hemisphere pairs of sister species recognised in Figs 24, because these species are already prevented from interbreeding and given evolutionary independence through lack of gene flow, but caused in this case by ocean barriers (see below). Consequently there would be no selection for further reinforcement of barriers to interbreeding by causing divergence in sex pheromones.

Interbreeding and hybridisation are often seen as criteria central to species concepts [7173]. Interbreeding is unlikely between the populations that we recognise here as sister species from their gene coalescents (e.g. between B. polaris and B. pyrrhopygus, between B. balteatus and B. kirbiellus, and between B. hyperboreus and B. natvigi) because interbreeding would have disrupted these coalescents. In each case, the distributions of these sister species are now separated so that they are prevented from dispersing and interbreeding by oceans, by 82 km at the Bering Strait and by >1000 km at the Greenland Sea. It is possible that bumblebees might be capable of dispersing across 82 km of sea, but it is unlikely that enough of them would arrive in a suitable physiological condition either to establish colonies or to mate successfully. But more persuasively, we see no evidence of dispersal, interbreeding, or introgression in this group from conflicts between the genetic and geographical data. In these cases, allopatric speciation is likely to have followed the repeated loss of the Bering land connection (which took place over as much as the last 5–7 Ma [84, 85]) according to dates estimated from our phylogenetic tree (Fig 5).

Another response to conflicting evidence has been to seek to accommodate the conflict within a hierarchy, by interpreting the broader groups as species and the narrower groups as subspecies [77, 86] (i.e. by choosing to reject the narrower groups as species). However, not only were gene-coalescent-based methods like GMYC specifically designed to discover species [49, 57, 60], not subspecies, but the justification for recognising subspecies is now far from universally accepted. Subspecies were introduced into bumblebee taxonomy originally as a typological concept, long before the biological species concept, and as a way of giving names to groups of individuals with differing colour patterns within a species [10]. Unfortunately, in practice the concept of subspecies has been applied with even less consistent meaning than the concept of species, to the point where subspecies have been considered to have hindered progress in taxonomy, evolutionary studies, and conservation [8789]. For example, a recent review of the well-known bumblebee B. terrestris defines several subspecies from across Europe that include both discrete colour-pattern groups on islands, as well as dividing segments from more continuous colour-pattern clines across the continent [90]. It is unavoidable that dividing continuous clines is essentially arbitrary. This distinction between islands and the continent in colour-pattern variation is paralleled in evidence from normally highly variable DNA microsatellites, which also support many island groups as distinct, but which show no significant differentiation within the continental population [91]. We wish to avoid sometimes weakly-differentiated colour-pattern groups being over-interpreted in terms of (e.g.) ecology when this has not been justified by a thorough geographical analysis. At the same time, COI-barcode variation shows a single coalescent for the entire species [19] (see below regarding Corsica). To help avoid ambiguous differences in the meaning of subspecies in different cases, we agree that subspecies and trinomials should be replaced, either by direct informal descriptions of the particular colour patterns on which they are usually actually based and which they are used to indicate [67], or preferably by direct descriptions of the underlying genetic patterns from phylogeographic studies [30, 78].

Alternatively, when called upon to make taxonomic decisions in the face of conflicting evidence, we are prompted to consider which criterion within our framework is most appropriate and reliable for recognising species. If we accept the recent move from viewing species as morphotypes towards viewing them as EILs, then gene coalescents, which are very direct evidence of evolutionary independence, may be considered a special case as evidence for recognising species [1, 3, 49]. We still need to consider carefully whether sufficient conditions are likely to have been met to expect that fixed coalescents will have formed in each of these cases, for example whether populations are too large, or whether there is residual gene flow. We also accept that making the choice to view species differently (i.e. not as morphotypes) might alter the shape of taxonomy to some degree [2]. It has been argued that coalescence methods could reduce investigator-driven biases in species delimitation [49] and this should increase the precision in the groups recognised as species.

Both the GMYC and PTP methods have potential pitfalls and have to be applied with care [3, 5961, 65]. For example, both are likely to split any samples isolated on trees by long terminal branches. These branches can become exaggerated as artefacts arising from a variety of causes, including sequencing errors [65], short sequences [59, 92], as well as unrepresentative sampling. For example, when samples from only Corsica and adjacent Europe were analysed with GMYC models, B. xanthopus Kriechbaumer was interpreted as an endemic Corsican bumblebee species separate from the mainland B. terrestris [86]. But when samples from throughout the known global distribution of B. terrestris were analysed [19], including samples from Madeira, the Canary Islands, North Africa, Europe, Russia, Iran, Central Asia, and from as far east as the indigenous eastern limit of the species in Mongolia, then the Corsican ‘xanthopus’ samples were found to be part of the species B. terrestris, closely related to particular subgroups within that species (consistent with the results of an earlier analysis [91]). Sampling can never be ‘complete’, but best practice in morphological taxonomy [48] accepts that it is essential to sample from across the breadth of global ranges of all included taxa when revising any taxonomic group. The same principle holds when sampling genes, as we have attempted to do here (Fig 1). GMYC may also not work well when there has been extremely rapid speciation, either very recently or if populations were separated over a very short period long ago followed by a long slow divergence [59]. There is no reason to believe that these conditions apply in this case (Fig 5).

Some of our prospective species are currently known from few haplotypes, even though they represent larger samples. For the Old World B. hyperboreus, we have only two unique haplotypes shown in our trees (Figs 2 and 4), although these come from a sample of 13 Old World sequences. Six morphologically similar individuals from Greenland share the haplotype of B. natvigi from Canada, the oldest available name [39] for the New World sister species to the Old World species B. hyperboreus. In contrast, more material is needed to clarify the status of the ‘unnamed’ taxon, represented here by sequences from five individuals collected from high elevations in southern Alaska (Alaska range) and the Yukon (Saint Elias Mountains), sequences which differ from one another at two nucleotide positions. These sequences are unlikely to be paralogous nuclear pseudogenes of mtDNA, or ‘numts’ [93], because they have no indels, no stop codons, and have a higher variability with a strong bias towards high A/T frequency at codon position 3 [94]. Heteroplasmy has been found in bees [95], but we examined the trace files for these sequences and found no evidence of double or irregular peaks. Another possible explanation for distinguishing ‘unnamed’ might be incomplete lineage sorting [3], resulting in what might appear in our trees as paraphyly between B. neoboreus and the taxon ‘unnamed’. The counter argument is that the six B. neoboreus sequences (three of them sequenced twice, spanning a range from very pale to very dark specimens, with three unique haplotypes) differ from the five ‘unnamed’ sequences at 19 nucleotide positions (two of them causing amino acid changes at translation). This is sufficient for the GMYC analysis to recognise them as separate species. Independent corroborative support for the unnamed taxon as a separate species is shown by unique diagnostic base changes in the PEPCK gene (Table 2: 454 specimens have been examined for the two taxa combined, although most are too old for us to sequence). The biotic history of the Arctic and of the Beringian refuge may be complex [84, 96]: it is possible that B. neoboreus and the unnamed species could represent relict populations left from successive waves of range expansion and contraction following the cycles of climate change associated with the Pleistocene glaciations. For comparison, high elevations of the Yukon’s Kluane mountains are unusual for having outlying disjunct and genetically divergent populations of an otherwise arctic moth, Gynaephora groenlandica (Wocke) (Erebidae), and of the endemic arctic plants, Oxytropis arctica R. Br. (Fabaceae) and Puccinellia vahliana (Liebm.) Scrib. and Merr. (Poaceae), which have been interpreted as isolated relict populations [97].

Our results stand out compared with other revisionary studies of entire monophyletic bumblebee subgenera using genes [19, 76] in finding evidence that supports nearly twice as many species as had been recognised in morphology-based revisions. Most of the taxa of the subgenus Alpinobombus have been known for more than a century and only the one ‘unnamed’ species may not have been examined for previous revisions (Table 1). The other species recognised here for the first time from gene coalescents differ from the species splits in previous revisions (Table 1), which were based on differences in colour pattern. A large increase in species number is all the more surprising because the Arctic is expected to have relatively low biodiversity [98, 99]. Nonetheless, both the GMYC and PTP results (Figs 24) agree in separating Old World from New World sister populations for B. polaris / B. pyrrhopygus, B. balteatus / B. kirbiellus, and B. hyperboreus / B. natvigi. Discontinuity of populations as geographical disjunction is also a feature of morphologically cryptic species in other genetic studies of insects [100]. These patterns are at odds with ideas in some earlier revisions of the subgenus Alpinobombus, which had concluded that both sister species in some of these species pairs might occur either within the Old World (both of B. polaris and B. pyrrhopygus) [37] or within the New World (both of B. balteatus and B. kirbiellus) [36]. But, for B. kirbiellus, even the rare individuals with black hair on the side of the thorax (resembling many Old World B. balteatus) from as far north in the New World as Ellesmere Island do belong to B. kirbiellus according to our genetic results. Similarly, we infer that both banded and unbanded individuals belong to the Eurasian B. pyrrhopygus. The same nine species were also recognised by the ‘refined single-linkage’ (RESL) clustering procedure (single-linkage clustering followed by Markov clustering) that generates BOLD’s (boldsystems.org) Barcode Index Number (BIN) system from COI barcodes [101]. The RESL procedure is related to the empirical ‘barcoding gap’ idea. Consequently, while it is computationally fast, it lacks the theoretical justification of the coalescence-based approach, which is likely to cause results to diverge in some cases.

Unexpected and morphologically cryptic diversity has also been detected from genetic evidence in other groups of insects that had previously been studied intensively and had been considered taxonomically mature [100, 102]. Among the species recognised here but not in the most recent revisions, B. kirbiellus has been regarded as a separate species by just one reviser [36], although not with quite the same concept (the population from northern Canada was excluded). In contrast, B. pyrrhopygus, B. natvigi, and the unnamed species have no closely corresponding concepts of species in earlier revisions (Table 1). Intensity of sampling may be a factor in the discovery of the unnamed species, but for the others, the geographical breadth of this study is likely to be a key factor (along with the use of genetic data) allowing critical assessment of the variation.

Our results for the subgenus Alpinobombus are the first we know of for bumblebees to demonstrate support for close resemblance among parts of different species populations arising from shared ancestral polymorphisms (Fig 8). Many other cases of resemblance among more distantly related bumblebee species have been ascribed to evolutionary convergence [25, 26, 30]. Whereas spatial segregation (by latitude, longitude, or elevation) is characteristic of many of these convergent bumblebee-colour-pattern groups elsewhere in the world [26], we find broad spatial overlap between the individuals of Alpinobombus with the banded pale-tailed (BP) and the banded dark-tailed (BD) patterns. There may be some variation in the relative frequency of the two patterns with latitude (Fig 7). The unbanded pale-tailed (UP) pattern is more highly concentrated (longitudinally) in Europe, although it is also present in a few locations in both Asia and North America. The majority of species show more than one of these principal colour patterns, and one third of the species show all three (Table 3). This analysis uses our combined COI and PEPCK tree as an estimate of the species tree, which matches a tree for five of the species from five genes (including PEPCK but not COI) [32]. Therefore we feel it is reasonable at present to treat the combined COI and PEPCK tree as an estimate of the species tree. In any case, no estimate of phylogeny from more genes could alter the inference that ancestral polymorphisms are most likely to explain the pattern of extant polymorphisms within Alpinobombus, as long as a similar set of species is supported.

The principal weakness of our analysis of colour-pattern characters is that we do not yet know how colour patterns are inherited for bumblebees [103]. A diversity of possibilities for genes and molecular-development mechanisms is known to exist even among closely related species of flies within the genus Drosophila and among some mimetic butterflies, although within the butterfly genus Heliconius there is also evidence of a broadly conserved genetic basis for mimicry [104]. The simplest explanation for the pattern of states for the two characters of the subgenus Alpinobombus in Fig 8 is that in both cases divergence of the principal colour-pattern states preceded the divergence of many of the species. However, even though this result is more consistent with ancestral polymorphisms for both colour-pattern characters and this result is necessary in order to support the idea, it is not sufficient to prove the case. Ultimately proof will require identification of the genes governing these colour-pattern elements and the demonstration of homology of these genes among the species.

Why would an ancestral polymorphism persist, without one colour pattern becoming fixed in all populations by drift or selection? Local samples of the subgenus Alpinobombus we have examined do appear to show that local populations of a single species are often polymorphic (e.g. within B. polaris and Russian B. pyrrhopygus). Are these different colour patterns selectively neutral, or could there be heterogeneity in selection that actually favours the polymorphism, perhaps depending on temporal or spatial variation in thermal constraints or in the abundance of predators or ‘mimics’ in different microhabitats? Curiously, entirely black individuals of Alpinobombus are not present in our sample. Within Scandinavia, B. balteatus shows a tendency towards an increased frequency of darker colour patterns in the more southern mountains ([37]: her figure 69), which Løken associated with areas that have high humidity near the ground. Pekkarinen [41] recorded ‘almost completely black’ Scandinavian males of B. balteatus and described melanism in this species as ‘typical high-altitude or high-latitude’ melanism. Contrary to this assertion, the highest frequency of bumblebee species with many entirely black individuals world-wide occurs in the tropical lowlands ([26]; see also [105]). Our observations here show that for Alpinobombus, individuals from the highest latitudes (northern Greenland and Novaya Zemlya) are characterised by extensive yellow bands, although often with black tails, and there does appear to be some increase in the frequency of black-tailed individuals towards the north (Fig 7). The thermal properties of bumblebee colour patterns need further study [26].

Conclusions

To examine gene coalescents as at least part of the evidence for discovering bumblebee species, studies should focus on obtaining: (a) homologous gene sequences; (b) comparable-length sequences; (c) sequences from multiple independent genes (which need to be rapidly evolving); and should (d) apply coalescence-based methods such as GMYC analysis. These aims supplement the more general sampling aims long recognised by morphological taxonomists [48] of the essential need to: (e) study a taxonomic group throughout its entire global range so as to include all constituent lineages or taxa; and (f) study many samples from throughout each constituent taxon’s range, both of which are just as important in genetic studies. Much more work is needed to elucidate the genetic control of bumblebee colour patterns and to assess the relative roles of ancestral polymorphisms and convergence in governing polymorphisms within species.

Supporting Information

S1 File. Table A. Depositories.

Collections from which pinned material of the subgenus Alpinobombus has been examined.

https://doi.org/10.1371/journal.pone.0144544.s001

(DOC)

S2 File. Table A. Sequences.

Sequenced samples of bumblebees (genus Bombus) of the subgenus Alpinobombus and outgroups with their accession numbers for GenBank (G) and the public project folder BBAL of the BOLD online database (B).

https://doi.org/10.1371/journal.pone.0144544.s002

(DOC)

Acknowledgments

Thanks to all who generously donated, loaned, or photographed specimens, especially to J. Ascher, Y. Astafurova, C. Buddle, S. Cannings, S. Cardinal, S. Colla, D. Currie, S. Droege, L. Gall, A. Guidotti, B. Harris, K. Martins, L. Packer, M. Proshchalykin, S. Schmidt, V. Scott, D. Sikes, J. Strange, J. Thomas, J. Thomson, J. Weintraub, T. Wheeler, D. Yanega; to the Biodiversity Institute of Ontario for barcoding of North American material; to the Norwegian Barcode of Life Network for barcoding of Norwegian material; to L. Bailey, T. Barraclough, S. Cameron, S. Cannings, M. Kuhlmann, W. Maddison, and A. Vogler for discussion; and to P. Foster for help in setting up the bGMYC software.

Author Contributions

Conceived and designed the experiments: PHW. Performed the experiments: PHW STW. Analyzed the data: PHW. Contributed reagents/materials/analysis tools: PHW AMB BC MVB FØ CR LLR JH CSS. Wrote the paper: PHW.

References

  1. 1. de Queiroz K. Species concepts and species delimitation. Systematic Biology. 2007;56(6):879–86. pmid:18027281
  2. 2. Schlick-Steiner BC, Steiner FM, Seifert B, Stauffer C, Christian E, Crozier RH. Integrative taxonomy: a multisource approach to exploring biodiversity. Annual Review of Entomology. 2010;55:421–38. pmid:19737081
  3. 3. Leliaert F, Verbruggen H, Vanormelingen P, Steen F, Lopez-Bautista JM, Zuccarello GC, et al. DNA-based species delimitation in algae. European Journal of Phycology. 2014;49(2):179–96.
  4. 4. Cardoso A, Serrano A, Vogler AP. Morphological and molecular variation in tiger beetles of the Cicindela hybrida complex: is an 'integrative taxonomy' possible? Molecular Ecology. 2009;18:648–64. pmid:19175505
  5. 5. Michener CD. The bees of the world. 2 ed. Baltimore: John Hopkins University Press; 2007. xiv+913 p.
  6. 6. Linnaeus C. Systema naturae per regna tria naturae, secundum classes, ordines, genera, species, cum characteribus, differentiis, synonymis, locis. Holmiae1758. 823 p.
  7. 7. Radoszkowski O. Révision des armures copulatrices des mâles du genre Bombus. Byulletin’ Moskovskogo obshchestva ispytatelei prirody. 1884;59(1/1):51–92.
  8. 8. Vogt O. Studien über das Artproblem. 1. Mitteilung. Über das Variieren der Hummeln. 1. Teil. Sitzungsberichte der Gesellschaft naturforschender Freunde zu Berlin. 1909;1909:28–84.
  9. 9. Friese H, von Wagner F. Zoologische Studien an Hummeln. I. Die Hummeln der deutschen Fauna. Zoologische Jahrbücher, Abteilung für Systematik, Ökologie und Geographie der Tiere. 1910;29:1–104.
  10. 10. Williams PH. An annotated checklist of bumble bees with an analysis of patterns of description (Hymenoptera: Apidae, Bombini). Bulletin of The Natural History Museum (Entomology). 1998;67:79–152 [updated at www.nhm.ac.uk/bombus/ accessed 2015].
  11. 11. Duennes MA, Lozier JD, Hines HM, Cameron SA. Geographical patterns of genetic divergence in the widespread Mesoamerican bumble bee Bombus ephippiatus (Hymenoptera: Apidae). Molecular Phylogenetics and Evolution. 2012;64:219–31. pmid:22521295
  12. 12. Huang J-X, Wu J, An J-D, Williams PH. Newly discovered colour-pattern polymorphism of Bombus koreanus females (Hymenoptera: Apidae) demonstrated by DNA barcoding. Apidologie. 2014;46(2):250–61.
  13. 13. Gjershaug JO, Staverløkk A, Kleven O, Ødegaard F. Species status of Bombus monticola Smith (Hymenoptera: Apidae) supported by DNA barcoding. Zootaxa. 2013;3716(3):431–40.
  14. 14. Williams PH, Byvaltsev AM, Sheffield CS, Rasmont P. Bombus cullumanus–an extinct European bumblebee species? Apidologie. 2012;44(2):121–32.
  15. 15. Carolan JC, Murray TE, Fitzpatrick U, Crossley J, Schmidt H, Cederberg B, et al. Colour patterns do not diagnose species: quantitative evaluation of a DNA barcoded cryptic bumblebee complex. PLoSONE. 2012;7(1):662–7.
  16. 16. Rasmont P. Les bourdons du genre Bombus Latreille sensu stricto en Europe occidentale et centrale (Hymenoptera, Apidae). Spixiana. 1984;7:135–60.
  17. 17. Krüger E. Phaenoanalytische Studien an einigen Arten der Untergattung Terrestribombus O. Vogt (Hymen. Bomb.). I. Teil. Tijdschrift voor Entomologie. 1951;93:141–97.
  18. 18. Bertsch A, Schweer H, Titze A, Tanaka H. Male labial gland secretions and mitochondrial DNA markers support species status of Bombus cryptarum and B. magnus (Hymenoptera, Apidae). Insectes Sociaux. 2005;52:45–54.
  19. 19. Williams PH, Brown MJF, Carolan JC, An J-D, Goulson D, Aytekin AM, et al. Unveiling cryptic species of the bumblebee subgenus Bombus s. str. world-wide with COI barcodes (Hymenoptera: Apidae). Systematics and Biodiversity. 2012;10(1):21–56.
  20. 20. von Dalla Torre KW. Unsere Hummel- (Bombus) Arten. Der Naturhistoriker. 1880;2:40–1.
  21. 21. Vogt O. Studien über das Artproblem. 2. Mitteilung. Über das Variieren der Hummeln. 2. Teil. (Schluss). Sitzungsberichte der Gesellschaft naturforschender Freunde zu Berlin. 1911;1911:31–74.
  22. 22. Reinig WF. On the variation of Bombus lapidarius L. and its cuckoo, Psithyrus rupestris Fabr., with notes on mimetic similarity. Journal of Genetics. 1935;30:321–56.
  23. 23. Reinig WF. Die Evolutionsmechanismen, erläutert an den Hummeln. Verhandlungen der Deutschen zoologischen Gesellschaft (supplement). 1939;12:170–206.
  24. 24. Tkalcu B. Revision der vier sympatrischen, homochrome geographische Rassen bildenden Hummelarten SO-Asiens (Hymenoptera, Apoidea, Bombinae). Annotationes Zoologicae et Botanicae. 1968;52:1–31.
  25. 25. Plowright RC, Owen RE. The evolutionary significance of bumble bee color patterns: a mimetic interpretation. Evolution. 1980;34:622–37.
  26. 26. Williams PH. The distribution of bumblebee colour patterns world-wide: possible significance for thermoregulation, crypsis, and warning mimicry. Biological Journal of the Linnean Society. 2007;92:97–118.
  27. 27. Tkalcu B. Neue Taxa asiatischer Hummeln (Hymenoptera, Apoidea). Acta entomologica bohemoslovaca. 1989;86:39–60.
  28. 28. Williams PH. Do the parasitic Psithyrus resemble their host bumblebees in colour pattern? Apidologie. 2008;39:637–49.
  29. 29. Richards OW. Parallel colour-variations in humble-bees from the Himalayas. Proceedings of the Royal Entomological Society, London. 1929;3(3):75–6.
  30. 30. Hines HM, Williams PH. Mimetic colour pattern evolution in the highly polymorphic Bombus trifasciatus (Hymenoptera: Apidae) species complex and its comimics. Zoological Journal of the Linnean Society. 2012;166:805–26.
  31. 31. Baum D, Smith S. Tree thinking: an introduction to phylogenetic biology: Roberts and Company; 2012. 496 p.
  32. 32. Cameron SA, Hines HM, Williams PH. A comprehensive phylogeny of the bumble bees (Bombus). Biological Journal of the Linnean Society. 2007;91:161–88.
  33. 33. Rasmont P. Catalogue commenté des bourdons de la région ouest-paléarctique (Hymenoptera, Apoidea, Apidae). Notes Fauniques de Gembloux. 1983;7:1–71.
  34. 34. Reinig WF. Synopsis der in Europa nachgewiesenen Hummel- und Schmarotzerhummelarten (Hymenoptera, Bombidae). Spixiana. 1981;4:159–64.
  35. 35. Hurd Jr PD. Superfamily Apoidea. In: Krombein KV, editor. Catalog of Hymenoptera in America north of Mexico. 2. Washington D.C.1979. p. 1741–2209.
  36. 36. Milliron HE. A monograph of the western hemisphere bumblebees (Hymenoptera: Apidae; Bombinae). II. The genus Megabombus subgenus Megabombus. Memoirs of the Entomological Society of Canada. 1973;89:81–237.
  37. 37. Løken A. Studies on Scandinavian bumble bees (Hymenoptera, Apidae). Norsk entomologisk Tidsskrift. 1973;20:1–218.
  38. 38. Pittioni B. Neue und wenig bekannte Hummeln der Paläarktis (Hymenopt., Apidae). Konowia. 1939;17:244–63.
  39. 39. Richards OW. Some notes on the humble-bees allied to Bombus alpinus, L. Tromsø Museums Årshefter. 1931;50:1–32.
  40. 40. Richards KW. Biology of Bombus polaris Curtis and B. hyperboreus Schönherr at Lake Hazen, Northwest Territories (Hymenoptera: Bombini). Quaestiones entomologicae. 1973;9:115–57.
  41. 41. Pekkarinen A. Morphometric, colour and enzyme variation in bumblebees (Hymenoptera, Apidae, Bombus) in Fennoscandia and Denmark. Acta zoologica fennica. 1979;158:pp. 60.
  42. 42. Williams PH, Thorp RW, Richardson LL, Colla SR. Bumble bees of North America. An identification guide. Princeton New Jersey: Princeton University Press; 2014. 208 p.
  43. 43. Skorikov AS. Die grönländischen Hummeln im Aspekte der Zirkumpolarfauna. Entomologiske Meddelelser. 1937;20(1):37–64.
  44. 44. Skorikov AS. Die Hummelfauna Turkestans und ihre Beziehungen zur zentralasiatischen Fauna (Hymenoptera, Bombidae). In: Lindholm VA, editor. Abhandlungen der Pamir-Expedition 1928. 8. Leningrad: Academy of Sciences of the USSR; 1931. p. 175–247.
  45. 45. Svensson BG, Lundberg H. Distribution of bumble bee nests in a subalpine/alpine area in relation to altitude and habitat (Hymenoptera, Apidae). Zoon. 1977;5:63–72.
  46. 46. Panfilov DV, Shamurin VF, Jurtsev BA. [On the conjugate distribution of bumblebees and legumes in the Arctic]. Bulluten Moskovskogo Obshchestva Ispytateley Prirody, Otdel Biologicheskiy (Moscow). 1960;65(3):53–62.
  47. 47. Maddison WP, Knowles LL. Inferring phylogeny despite incomplete lineage sorting. Systematic Biology. 2006;55(1):21–30. pmid:16507521
  48. 48. Bolton B. How to conduct large-scale taxonomic revisions in Formicidae. In: Snelling RR, Fisher BL, Ward PS, editors. Advances in ant systematics (Hymenoptera: Formicidae): homage to E O Wilson—50 years of contributions. 80: Memoirs of the American Entomological Institute; 2007. p. 52–71.
  49. 49. Fujita MK, Leaché AD, Burbrink FT, McGuire JA, Moritz C. Coalescent-based species delimitation in an integrative taxonomy. Trends in Ecology and Evolution. 2012;27(9):480–8. pmid:22633974
  50. 50. Zink RM, Barrowclough GF. Mitochondrial DNA under siege in avian phylogeography. Molecular Ecology. 2008;17:2107–21. pmid:18397219
  51. 51. Hebert PDN, Ratnasingham S, deWaard JR. Barcoding animal life: cytochrome c oxidase subunit 1 divergences among closely related species. Proceedings of the Royal Society of London (B). 2003;270:S96–S9.
  52. 52. Hebert PDN, Penton EH, Burns JM, Janzen DH, Hallwachs W. Ten species in one: DNA barcoding reveals cryptic species in the neotropical skipper butterfly Astraptes fulgerator. Proceedings of the National Academy of Sciences. 2004;101:14812–7.
  53. 53. Strange JP, Knoblett J, Griswold T. DNA amplification from pin-mounted bumble bees (Bombus) in a museum collection: effects of fragment size and specimen age on successful PCR. Apidologie. 2009;40:134–9.
  54. 54. Tamura K, Stecher G, Paterson D, Filipski A, Kumar S. MEGA6: molecular evolutionary genetics analysis version 6.0. Molecular Biology and Evolution. 2013;30:2725–9. pmid:24132122
  55. 55. Edwards SV. Is a new and general theory of molecular systematics emerging? Evolution. 2009;63(1):1–19. pmid:19146594
  56. 56. Monaghan MT, Balke M, Gregory TR, Vogler AP. DNA-based species delineation in tropical beetles using mitochondrial and nuclear markers. Philosophical Transactions of the Royal Society (B). 2005;360:1925–33.
  57. 57. Pons J, Barraclough TG, Gomez-Zurita J, Cardoso A, Duran DP, Hazell S, et al. Sequence-based species delimitation for the DNA taxonomy of undescribed insects. Systematic Biology. 2006;55(4):595–609. pmid:16967577
  58. 58. Papadopoulou A, Bergsten J, Fujisawa T, Monaghan MT, Barraclough TG, Vogler AP. Speciation and DNA barcodes: testing the effects of dispersal on the formation of discrete sequence clusters. Philosophical Transactions of the Royal Society (B). 2008;363:2987–96.
  59. 59. Reid NM, Carstens BC. Phylogenetic estimation error can decrease the accuracy of species delimitation: a Bayesian implementation of the general mixed Yule-coalescent model. BMC Evolutionary Biology. 2012;12:196. pmid:23031350
  60. 60. Fujisawa T, Barraclough TG. Delimiting species using single-locus data and the Generalized Mixed Yule Coalescent approach: a revised method and evaluation on simulated data sets. Systematic Biology. 2013;62(5):707–24. pmid:23681854
  61. 61. Talavera G, Dinca V, Vila R. Factors affecting species delimitations with the GMYC model: insights from a butterfly survey. Methods in Ecology and Evolution. 2013;4:1101–10.
  62. 62. Monaghan MT, Wild R, Elliot M, Fujisawa T, Balke M, Inward DJG, et al. Accelerated species inventory on Madagascar using coalescent-based models of species delineation. Systematic Biology. 2009;58(3):298–311. pmid:20525585
  63. 63. Drummond AJ, Rambaut A. BEAST: Bayesian evolutionary analysis by sampling trees. BioMed Central Evolutionary Biology. 2007;7:214. pmid:17996036
  64. 64. ICZN. International code of zoological nomenclature. 4 ed. London: International Commission on Zoological Nomenclature; 1999. xxix+306 p.
  65. 65. Zhang J-J, Kapli P, Pavlidis P, Stamatakis A. A general species delimitation method with applications to phylogenetic placements. Bioinformatics. 2013;29(22):2869–76. pmid:23990417
  66. 66. Huelsenbeck JP, Ronquist F. MRBAYES: Bayesian inference of phylogeny. Bioinformatics. 2001;17:754–5. pmid:11524383
  67. 67. Williams PH. The bumble bees of the Kashmir Himalaya (Hymenoptera: Apidae, Bombini). Bulletin of the British Museum (Natural History) (Entomology). 1991;60:1–204.
  68. 68. Hines HM. Historical biogeography, divergence times, and diversification patterns of bumble bees (Hymenoptera: Apidae: Bombus). Systematic Biology. 2008;57(1):58–75. pmid:18275002
  69. 69. Maddison WP, Maddison DR. Mesquite: a modular system for evolutionary analysis. 3.02 ed2015.
  70. 70. Maddison DR, Maddison WP. MacClade version 4: analysis of phylogeny and character evolution. Sunderland, Massachusetts: Sinauer Associates; 2000.
  71. 71. Mallet J. What are species? Trends in Ecology and Evolution. 1997;12:453–4.
  72. 72. Mayr E. Animal species and evolution. Massachusetts1963. xiv+797 p.
  73. 73. De Meulemeester T. Approche intégrative dans la systématique de taxons complexes: bourdons et abeilles fossiles. Mons: Université de Mons; 2012.
  74. 74. Paterson HEH. The recognition concept of species. In: Vrba ES, editor. Species and speciation. Pretoria1985. p. 21–9.
  75. 75. Ackery PR, Vane-Wright RI. Milkweed butterflies, their cladistics and biology, being an account of the natural history of the Danainae, a subfamily of the Lepidoptera, Nymphalidae. London1984. ix+425 p.
  76. 76. Williams PH, An J-D, Huang J-X. The bumblebees of the subgenus Subterraneobombus: integrating evidence from morphology and DNA barcodes (Hymenoptera, Apidae, Bombus). Zoological Journal of the Linnean Society. 2011;163:813–62.
  77. 77. Lecocq T, Dellicour S, Michez D, Dehon M, Dewulf A, De Meulemeester T, et al. Methods for species delimitation in bumblebees (Hymenoptera, Apidae, Bombus): towards an integrative approach. Zoologica Scripta. 2015;44:281–97.
  78. 78. Lecocq T, Dellicour S, Michez D, Lhomme P, Vanderplanck M, Valterova I, et al. Scent of a break-up: phylogeography and reproductive trait divergences in the red-tailed bumblebee (Bombus lapidarius). BMC Evolutionary Biology. 2013;13:263. pmid:24295171
  79. 79. Bertsch A, Schweer H, Titze A. Discrimination of the bumblebee species Bombus lucorum, B. cryptarum and B. magnus by morphological characters and male labial gland secretions. Beiträge zur Entomologie. 2004;54(2):365–86.
  80. 80. Pamilo P, Tengö J, Rasmont P, Pirhonen K, Pekkarinen A, Kaarnama E. Pheromonal and enzyme genetic characteristics of the Bombus lucorum species complex in northern Europe. Entomologica fennica. 1997;7:187–94.
  81. 81. Svensson BG. Pyrobombus lapponicus auct., in Europe recognized as two species: P. lapponicus (Fabricius, 1793) and P. monticola (Smith, 1849) (Hymenoptera, Apoidea, Bombinae). Entomologica scandinavica. 1979;10:275–96.
  82. 82. Terzo M, Urbanova K, Valterova I, Rasmont P. Intra and interspecific variability of the cephalic labial glands' secretions in male bumblebees: the case of Bombus (Thoracobombus) ruderarius and B. (Thoracobombus) sylvarum [Hymenoptera, Apidae]. Apidologie. 2005;36:85–96.
  83. 83. Lecocq T, Coppée A, Mathy T, Lhomme P, Cammaerts-Tricot M- C, Urbanova K, et al. Subspecific differentiation in male reproductive traits and virgin queen preferences, in Bombus terrestris. Apidologie. 2015.
  84. 84. Pringle H. Welcome to Beringia. Science. 2014;343(28 February):961–3. pmid:24578560
  85. 85. Brigham-Grette J. New perspectives on Beringian Quaternary paleogeography, stratigraphy, and glacial history. Quaternary Science Reviews. 2001;20:15–24.
  86. 86. Lecocq T, Brasero N, DeMeulemeester T, Michez D, Dellicour S, Lhomme P, et al. An integrative taxonomic approach to assess the status of Corsican bumblebees: implications for conservation. Animal Conservation. 2014.
  87. 87. Barrowclough GF. Geographic variation, predictiveness, and subspecies. Auk. 1982;99:601–3.
  88. 88. Wilson EO, Brown WL. The subspecies concept and its taxonomic application. Systematic Zoology. 1953;2:97–111.
  89. 89. Zink RM. The role of subspecies in obscuring avian biological diversity and misleading conservation policy. Proceedings of the Royal Society of London (B). 2004;271:561–4.
  90. 90. Rasmont P, Coppée A, Michez D, Meulemeester TD. An overview of the Bombus terrestris (L. 1758) subspecies (Hymenoptera: Apidae). Annales de la Société entomologique de France. 2008;44(1):243–50.
  91. 91. Estoup A, Solignac M, Cornuet J-M, Goudet J, Scholl A. Genetic differentiation of continental and island populations of Bombus terrestris (Hymenoptera: Apidae) in Europe. Molecular Ecology. 1996;5:19–31. pmid:9147693
  92. 92. Williams ST, Ozawa T. Molecular phylogeny suggests polyphyly of both the turban shells (family Turbinidae) and the superfamily Trochoidea (Mollusca: Vestigastropoda). Molecular Phylogenetics and Evolution. 2006;39:33–51. pmid:16483804
  93. 93. Song H, Buhay JE, Whiting MF, Crandall KA. Many species in one: DNA barcoding overestimates the number of species when nuclear mitochondrial pseudogenes are coamplified. Proceedings of the National Academy of Sciences of the United States of America. 2008;105:13486–91. pmid:18757756
  94. 94. Cristiano MP, Fernandes-Salomão TM, Yotoko KSC. Nuclear mitochondrial DNA: an Achilles' heel of molecular systematics, phylogenetics, and phylogeographic studies of stingless bees. Apidologie. 2012;43:527–38.
  95. 95. Magnacca KN, Brown MJF. Mitochondrial heteroplasmy and DNA barcoding in Hawaiian Hylaeus (Nesoprosopis) bees (Hymenoptera: Colletidae). BMC Evolutionary Biology. 2010;10(174):1–16.
  96. 96. Elias SA, Brigham-Grette J. Late Pleistocene glacial events in Beringia. Encyclopedia of Quaternary Science. London: Elsevier; 2013. p. 191–201.
  97. 97. Barro IC, Schmidt BC, Cannings S, Hik DS. First records of the Arctic moth Gynaephora groenlandica (Wocke) south of the Arctic Circle: a new alpine subspecies. Arctic. 2013;66(4):429–34.
  98. 98. Williams PH, Gaston KJ, Humphries CJ. Mapping biodiversity value worldwide: combining higher-taxon richness from different groups. Proceedings of the Royal Society of London (B). 1997;264:141–8.
  99. 99. Gaston KJ. Global patterns in biodiversity. Nature. 2000;405(11 May 2000):220–7. pmid:10821282
  100. 100. Rougerie R, Kitching IJ, Haxaire J, Miller SE, Hausmann A, Hebert PDN. Australian Sphingidae—DNA barcodes challenge current species boundaries and distributions. PLoS ONE. 2014;9(7):12.
  101. 101. Ratnasingham S, Hebert PDN. A DNA-based registry for all animal species: the barcode index number (BIN) system. PLoSONE. 2013;8(8):1–16.
  102. 102. Neumeyer R, Baur H, Guex G-D, Praz C. A new species of the paper wasp genus Polistes (Hymenoptera, Vespidae, Polistinae) in Europe revealed by morphometrics and molecular analyses. Zookeys. 2014;400:67–118. pmid:24843256
  103. 103. Rapti Z, Duennes MA, Cameron SA. Defining the colour pattern phenotype in bumble bees (Bombus): a new model for evo devo. Biological Journal of the Linnean Society. 2014;113:384–404.
  104. 104. Kronforst MR, Barsh GS, Kopp A, Mallet J, Monteiro A, Mullen SP, et al. Unraveling the thread of nature’s tapestry: the genetics of diversity and convergence in animal pigmentation. Pigment Cell Melanoma Research. 2012;25:411–33. pmid:22578174
  105. 105. Williams PH, Ito M, Matsumura T, Kudo I. The bumblebees of the Nepal Himalaya (Hymenoptera: Apidae). Insecta Matsumurana. 2010;66:115–51.