Pliocene Origin , Ice Ages and Postglacial Population Expansion Have Influenced a Panmictic Phylogeography of the European Bee-Eater Merops apiaster

Oscillations of periods with low and high temperatures during the Quaternary in the northern hemisphere have influenced the genetic composition of birds of the Palearctic. During the last glaciation, ending about 12,000 years ago, a wide area of the northern Palearctic was under lasting ice and, consequently, breeding sites for most bird species were not available. At the same time, a high diversity of habitats was accessible in the subtropical and tropical zones providing breeding grounds and refugia for birds. As a result of long-term climatic oscillations, the migration systems of birds developed. When populations of birds concentrated in refugia during ice ages, genetic differentiation and gene flow between populations from distinct areas was favored. In the present study, we explored the current genetic status of populations of the migratory European bee-eater. We included samples from the entire Palearctic-African distribution range and analyzed them via mitochondrial and nuclear DNA markers. DNA data indicated high genetic connectivity and panmixia between populations from Europe, Asia and Africa. Negative outcomes of Fu’s Fs and Tajima’s D tests point to recent expansion events of the European bee-eater. Speciation of Merops apiaster started during the Pliocene around three million years ago (Mya), with the establishment of haplotype lineages dated to the Middle Pleistocene period circa 0.7 Mya. M. apiaster, which breed in Southern Africa are not distinguished from their European counterparts, indicating a recent separation event. The diversification process of the European bee-eater was influenced by climatic variation during the late Tertiary and Quaternary. Bee-eaters must have repeatedly retracted to refugia in the Mediterranean and subtropical Africa and Asia during ice ages and expanded northwards during warm periods. These processes favored genetic differentiation and repeated lineage mixings, leading to a genetic panmixia, which we still observe today.


Introduction
Abiotic and biotic conditions have an important influence on the genetic structure of resident and long-distance migratory birds.When physical and/or ecological barriers limit gene flow between populations, phylogeographic differentiation can occur.Without barriers, gene flow will take place between populations.When differentiated lineages meet in refugia, lineages can mix, leading to genetic panmixia.These scenarios have been observed in phylogeographic studies of several migratory birds [1][2][3].Philopatry, habitat fragmentation, specific migratory flyways, and climatic conditions can influence the speciation of migratory taxa by limiting gene flow between populations and promoting genetic divergence [4].Species with differing migratory routes and destinations are usually genetically distinct, while species with similar breeding and wintering areas show weak genetic structure [2,5].
In the last two million years, a cyclical change of climate was observed in the northern hemisphere, with long warm and long cold periods.The Pleistocene climatic oscillations have led to compressions and expansions of the distribution range of many Palearctic bird taxa.During ice ages, populations retracted to distinct refugia in the Mediterranean, African and Asian subtropics.During warm periods, they left the refugia and distributed themselves over wide areas of the Palearctic again.These frequent cycles promoted genetic differentiation, gene flow and even panmixia [6].These environmental conditions have had important impacts not only on the bird richness and distribution but also on their migration systems [7].Thus, the migration systems of birds developed in response to long-term climatic cycles.They determined the phylogeographic pattern in many species [8].
The bee-eaters (family Meropidae) consist of widely distributed species in Africa, Asia, Europe, Australia, and Indonesia, inhabiting a variety of habitats from arid areas to wet forests [9].Phylogenetic relationships between the bee-eaters have been explored on the basis of DNA data [10] and morphological diagnosis [11].The phylogenetic analyses of all bee-eaters, based on mitochondrial and nuclear data, revealed two main clades, which differ in their migratory behavior: one lineage comprising basically sedentary species and one including mostly migratory species [10].
The long-distance migrant European bee-eater Merops apiaster exhibits a wide distribution in the Palearctic.Its breeding range covers a large area in Europe, extending to Western Asia, as well as the northwest and the extreme south of Africa [12].Although this species usually inhabits mainly dry areas, populations of bee-eaters have expanded towards Central/Southeast Europe, north of the Alps and Great Britain during the last decades, probably due to the occurrence of warmer summers [13].
The phylogenetic position and the genetic variability of the European bee-eater have been studied before [9][10][11]14].A low level of genetic inter-population divergence in M. apiaster had been observed [10].However, the detailed evolutionary history of M. apiaster and the clarification of the phylogenetic relationships between the bee-eaters is still a demand.Therefore, in this study we investigated sequences of mitochondrial marker genes and nuclear microsatellite markers: (i) to explore the evolutionary history of the migratory European bee-eater, (ii) to trace the genetic connectivity between breeding populations in Europe, Asia and Africa, and (iii) to estimate the speciation time of M. apiaster.Special emphasis was placed on populations in Germany (Rhineland-Palatinate, Baden-Württemberg, Saxony and Saxony-Anhalt), which have been established during the last two decades, and on populations which breed in Southern Africa.
The DNA was extracted from tissue, feathers, blood and swab samples by the standard phenol-chloroform protocol [15].The mitochondrial cytochrome b (Cyt b) and cytochrome oxidase I (COI) genes were amplified using the primers published from previous studies [3,16,17].PCR conditions were set as follows: 94 • C denaturation for 5 min, 52-62 • C for annealing temperature for 45 s, and 72 • C for extension for 2 min, running in 38 cycles.A final extension step was conducted at 72 • C for 10 min.The reaction was executed in a 40 µL final volume, each reaction containing 60 µg of DNA, 38.7 µL of H 2 O, 100 µM of dNTP, PCR Solution Buffer Complete (10x), 10 pmol of forward and reverse primers, and 2 units of Taq Polymerase.
The PCR product was purified by precipitation in 4 M NH 4 Ac and absolute ethanol (1:1:10), centrifuged at 13,000 rpm for a period of 30 min, followed by centrifugation with 70% ethanol using the same settings and later dissolved in 25 µL of sterile H 2 O.For sequencing, 1.0 µL of sequencing primer (10 pmol/µL) was combined with 7.0 µL of PCR products.The Sanger sequencing was executed on an ABI 3730 automated capillary sequencer (Applied Biosystems, Carlsbad, CA, USA) with the ABI Prism Big Dye Terminator Cycle Sequencing Ready Reaction Kit 3.1 (carried out by STARSEQ GmbH, Mainz, Germany).The nucleotide sequences were aligned using the program CodonCode aligner version 5 (http://www.codoncode.com/aligner).Additionally, we refined the alignment by eye and checked for stop codons in the sequences to exclude the possibility that nuclear copies had been amplified.

Demographic Analysis
DNAsp 5.0 [18] was used to access the demographic parameters for each population: number of haplotypes (H), number of polymorphic (segregating) sites (s), haplotype diversity (Hd), nucleotide diversity (Pi), Fu and Li's F, Tajima's D and Fu's Fs statistics.To test for signals of sudden demographic expansion [19], the mismatch distribution was estimated using the constant population size model and the evaluation of the data was verified through a goodness-of-fit test between the observed and simulated data (Raggedness statistic r) [20].

Phylogenetic Analyses and Molecular Dating
To employ the best evolutionary model for each gene we used jModelTest 2.1.7[21].HYK + InvGamma model is selected by the Akaike Information Criterion (AIC) to the Cyt b gene.The molecular timing tree estimation was assessed based on the evolutionary rate of the Cyt b gene using its known evolutionary rate (2.1% substitutions/site/lineage/million years) [22], set to strict model and Yule tree, the length of chain set to 10 7 , and sampling every 1000th tree in Beauti 1.8.0 [23].The phylogenetic tree was performed in Beast 1.8.0 [23] and the stationarity and convergence of runs were checked using Tracer 1.5 [24].The tree with the best likelihood was summarized using TreeAnnotator 1.8.0 [23].Other species belonging to the Merops genus were also included in the molecular dating phylogenetic tree: M. albicollis, M. boehmi, M. bullockoides, M. gularis, M. hirundineus, M. leschenaulti, M. malimbicus, M. nubicus, M. nubicoides, M. orientalis, M. persicus, and M. viridis [10].The analysis was conducted on the Cipres platform [25].

Genotyping and Microsatellite Analyses
A total of eleven polymorphic microsatellite loci were selected [26].The corresponding primers were used to genotype 279 specimens of M. apiaster.The thermal cycling was done under the following conditions: 95 • C for 5 min, 29 cycles of 95 • C denaturation for 30 s, 51-61 • C annealing temperature for 1.5 min, 72 • C for 30 s of extension and 4 • C for 15 min.We used the Type-it PCR Kit for Microsatellites (Qiagen, Cat.No. 206243, Hilden, Germany), see Reference [3] for the reaction details.
PCR products were analyzed on an ABI gel electrophoresis sequencer conducted by GATC Biotech AG (Konstanz, Germany).The software PEAK Scanner TM version 2.0 was used to determine the allele sizes.The alleles were arranged by size in classes and plotted in the most probable class of allele size.The presence of null alleles and genotyping mistakes (amplification mistakes, recording of stutter peaks and short allele dominance) were verified using the software Microchecker version 2.2.3 (University of Hull, Hull, United Kingdom) [27].Departures from the Hardy-Weinberg proportion and linkage equilibrium were verified for all loci using Genepop version 4.2 (Laboratiore de Genetique et Environment, Montpellier, France) [28], with significance estimated by Fisher exact, running the Markov chain algorithm with dememorization of 10,000.
The following parameters were calculated for all populations in Genetix version 4.05.2 (Université de Montpellier II, Montpellier, France) [29]: number of alleles (Na), allele richness (AR), observed heterozygosity (Ho), expected heterozygosity (He), gene diversity (HS), total gene diversity (HT).Fstat version 2.9.3 (UNIL, Lausanne, Switzerland) [30] was employed to determine the genetic polymorphism for each microsatellite locus (Na-number of alleles, Ar-average allelic richness, Ho-observed heterozygosity, Hs-gene diversity, HT-overall gene diversity, and Fst-genetic differentiation between populations).The levels of genetic differentiation between populations were determined using the Fst test [31] and the analysis of molecular variance (AMOVA) was conducted in Arlequin version 3.1 (Institute of Ecology and Evolution, University of Bern, Switzerland) [32].The discriminant analysis of principal components was performed using R package adegenet; this approach uses the k-mean clustering algorithm according to the Bayesian Information Criterion to detect genetic clusters in the population [33].
The genetic structure was tested with a Bayesian clustering method and implemented in Structure version 2.3.4 using the admixture model without prior information about the populations [34].MCMC parameters were fixed at 50,000 and 20,000, with run lengths of 10 6 iterations and K varying from 1 to 18 clusters.The calculations were repeated five times to confirm convergence among estimated settings.The most probable number of clusters (K) was based on the harmonic mean estimator [34] and on the second order rate of variation of the likelihood function referring to K (∆K) [35].
The program Bottleneck [36] was used to test the heterozygosity excess in each population to identify population bottlenecks.We used the two-phase model of microsatellite evolution [37].By means of GeneClass version 2.0 [38], we determined the direction of recent gene flow of M. apiaster between Europe and Africa by setting the parameters proposed by Paetkau et al. 2004 [39].

Sequence Data
In total, we sequenced the mtDNA Cyt b and COI of 175 samples.We included samples of M. apiaster available in Genbank as well.Overall, we analyzed Cyt b of 299 samples.We found 48 polymorphic sites and 34 parsimony informative sites in Cyt b.The dataset shows substantial genetic variability.We identified 56 haplotypes, 22 of them with more than one sample.Two dominant haplotypes came from several localities in Europe, Asia, and Southern Africa: H1 with 77 samples from Greece, Spain, Macedonia, Germany, Turkey, Southern Africa and China; and H2 with 57 samples from Greece, Spain, Macedonia, Germany, Slovakia, Tunisia, Southern Africa and China (Figure 1, Table 1).
The overall haplotype network based on Cyt b yielded two main haplotypes (Figure 1) and was congruent with the dataset of COI (Table S2).Various unique haplotypes of low frequency derived from them and formed a star-like structure.The unique haplotypes were distinguished from the main haplotype by only one or two mutation steps.Overall, bee-eaters from all localities shared both main common haplotypes (H1 and H2), with the exception of bee-eaters from Turkey, which showed only haplotype H1, and those from Slovakia, Tunisia and Georgia, which had only haplotype H2 (Figure 1, Table 1).It is remarkable, that bee-eaters, which breed in Southern Africa, did not reveal unique haplotypes, but shared haplotypes with European populations, indicating that they probably derived from them.
Additionally, the haplotype network structure (Figure 1) revealed recent expansion events, which have been confirmed by the negative results of the demographic population parameters, test Fu's Fs and Tajima's D (Table 1).Data on genetic diversity, quantified by haplotype diversity (Hd) and nucleotide diversity (Pi), is summarized in Table 2.The overall values demonstrated high haplotype diversity 0.878 (Hd) and low nucleotide diversity 0.0043 (Pi).To test the hypothesis of population expansion in M. apiaster, the distribution of pairwise differences was calculated using the DnaSP version 5.0.Unimodal mismatch distribution and the goodness-of-fit of the observed data (Raggedness statistic r) revealed a model of sudden population expansion, subsequently to the last glacial bottleneck event (Figure 2).This finding was also supported by the Raggedness statistic with r = 0.0722 between observed and simulated data.population expansion in M. apiaster, the distribution of pairwise differences was calculated using the DnaSP version 5.0.Unimodal mismatch distribution and the goodness-of-fit of the observed data (Raggedness statistic r) revealed a model of sudden population expansion, subsequently to the last glacial bottleneck event (Figure 2).This finding was also supported by the Raggedness statistic with r = 0.0722 between observed and simulated data.

Microsatellite Data
Out of the eleven microsatellite loci, six loci produced well-defined alleles, which were in Hardy-Weinberg and linkage equilibrium (BE129, BE72, BE192, BE216, BE233, and BE231).Therefore, we restricted the following analysis to these six microsatellite loci.Additionally, we did not find evidence for large allele dropout, null alleles and stuttering signals using Microchecker 2.2.3.The diversity estimates differed between the microsatellite loci (Table 3) and between populations analyzed (Table 4).The total number of alleles per locus ranged from 7 to 17 alleles and the average allelic richness (A r ) varied from 1.477 to 1.874.The average gene diversity (H s ) of all loci was 0.717 (range 0.456-0.867),whereas the overall gene diversity (H T ) range was 0.458-0.885(Table 3).Additionally, the lowest observed heterozygosity (H o ) and expected heterozygosity (H e ) were found for the population of Northern Africa, respectively 0.472 and 0.495, whereas the highest scores (H o = 0.750) were observed in Scandinavia and Northern Asia (Table 4).The pairwise differences Fst varied from 0 to 3.6% (Table 5).Interestingly, the Fst values from the population of Southern Africa were the lowest.The AMOVA results demonstrated that most of the molecular variance was distributed within populations (92%) rather than among groups (2.5%) or among populations (1.2%) (Table 6).Furthermore, the PCoA analysis grouped the populations into two clusters, however, samples from the distinct groups appear in all clades, evidencing a lack of genetic structure and indicating panmixia (Figure 3).Furthermore, the PCoA analysis grouped the populations into two clusters, however, samples from the distinct groups appear in all clades, evidencing a lack of genetic structure and indicating panmixia (Figure 3).Table 6.Analysis of molecular variance (AMOVA) for six microsatellite loci to test the variation among the following breeding populations of the bee-eaters: Group 1 (Spain, Germany (RLP, BW, Saxony-Anhalt), Greece); Group 2 (Denmark, Sweden); Group 3 (Slovakia); Group 4 (Turkey, Georgia, Kuwait); Group 5 (Egypt, Tunisia); Group 6 (Kazakhstan, Russia) and Group 7 (Southern Africa, Malawi).The estimated number of clusters obtained through the Bayesian cluster analysis, in a structure using the ∆K method, was six clusters.Overall, the different sampled bee-eater populations share genetic information; however, a slight differentiation was perceived in the populations from Western Asia, East Europe, Scandinavia and Southern Africa with K = 5 (Figure 4).No signal of a potential bottleneck was detected, based on the L shape of the distribution of allele frequencies observed in all the populations (Table S3).The estimated number of clusters obtained through the Bayesian cluster analysis, in a structure using the ΔK method, was six clusters.Overall, the different sampled bee-eater populations share The direction of the population movements based on the gene flow among the populations and estimated throughout the software GeneClass version 2.0 [38], detected the following routes: RLP/BW (western Germany) to Northern Asia, Saxony-Anhalt (eastern Germany) to West Asia and to RLP/BW, Northern Africa to East Europe, Southern Africa to Greece, Western Asia to East Europe and Greece to Western Asia.

Divergence Time Estimation
Molecular clock estimation, based on the evolutionary rate of Cyt b, revealed splitting events within the Meropidae family taking place during the Miocene/Pliocene/Pleistocene periods.The speciation of M. apiaster apparently took place in the mid-Pliocene with a highly supported sister relationship with M. malimbicus.The population establishment of the European bee-eater probably occurred during the late Pleistocene, about 0.6 Million years ago (Mya) (Figure 5).
Molecular clock estimation, based on the evolutionary rate of Cyt b, revealed splitting events within the Meropidae family taking place during the Miocene/Pliocene/Pleistocene periods.The speciation of M. apiaster apparently took place in the mid-Pliocene with a highly supported sister relationship with M. malimbicus.The population establishment of the European bee-eater probably occurred during the late Pleistocene, about 0.6 Million years ago (Mya) (Figure 5).

Historical Demography of the European Bee-eater
Our data reveal a weak phylogeographic structure between the populations of European beeeaters, which suggests high gene flow between populations from different localities, resulting in a genetic panmixia.These results are in agreement with previous DNA studies from other migratory birds of the Western Palearctic [2,3,40], which revealed low phylogeographic structures.The lack of genetic structure observed for the populations of M. apiaster is likely to be explained by shared ancestral polymorphism due to the recent divergence in this group.

Historical Demography of the European Bee-eater
Our data reveal a weak phylogeographic structure between the populations of European bee-eaters, which suggests high gene flow between populations from different localities, resulting in a genetic panmixia.These results are in agreement with previous DNA studies from other migratory birds of the Western Palearctic [2,3,40], which revealed low phylogeographic structures.The lack of genetic structure observed for the populations of M. apiaster is likely to be explained by shared ancestral polymorphism due to the recent divergence in this group.
The low genetic differentiation of Palearctic birds can be explained by the frequent climatic oscillation during the Pleistocene, which repeatedly forced the birds to retreat to more suitable regions in the Mediterranean, subtropical and tropical Africa and Asia during glacial periods [8].During warm periods, populations moved north again and differentiated locally.Migratory birds can present low genetic structure because they share reproductive sites; consequently, a reduction in the amount of exclusive of alleles within the populations occur.Thereby, populations from different localities in Europe likely mix with each other during the winter season and colonize other points after they return to Europe, which has probably resulted in the observed gene flow between birds from different populations [41].It is also possible, that the lineages mixed earlier when they retracted to refugia.Furthermore, some bird species have the ability to change the migratory pathway because of social or local environmental circumstances [42].This plasticity was necessary for the adaptation of birds in different regions and different climatic conditions [43], however it has not yet been recognized in European bee-eaters.
The outcome of the demographic tests of Tajima's D, Fu's Fs, and Fu and Li's F can explain the demographic forces that have defined the historical population dynamics of European bee-eaters.The negative values herein accessed reflect an excess of rare polymorphisms in all populations of M. apiaster, which is linked to population expansion and positive selection.Indeed, historical signs of population expansion northwards toward central Europe and Great Britain have been recorded for European bee-eaters, which were related to global warming and the appearance of new suitable habitats [13].
The slight genetic differentiation observed in the populations from Eastern Europe and Western Asia (Figure 4) might be associated with the distinct migration routes used by them.A similar pattern of the slight distinction between eastern and Western Europe was previously observed from microsatellite data [14].European bee-eaters use two migration routes to Africa: Bee-eaters breeding in Western Europe move across the Straits of Gibraltar to overwintering areas in Western Africa, between Senegal and Nigeria, and Eastern European populations travel down to Southern Africa past the Eastern Mediterranean Sea and the Sahara Desert.These two migratory routes seem to drive the populations to differentiate genetically into two distinct groups [14].

Lack of Differentiation of the African Populations from Merops apiaster
The population, which breeds in the western half of Southern Africa and Southern Namibia and spends the non-breeding season in Malawi, was expected to show unique alleles and to form a distinct cluster in the structure analyses.Nonetheless, they in fact share haplotypes and belong to the same genetic pool as the Eurasian migratory populations.As the African populations share the same haplotypes and microsatellite variation with Eurasian populations, they might have derived from migrants which remained in Southern Africa for breeding.

Molecular Dating of the European Bee-eater
The Neoaves appeared in the Cretaceous but rapidly radiated in the early Tertiary, with intensive speciation processes during the Eocene.The origin of the Coraciiformes, to which the Meropidae belong, is dated to the Late Eocene circa 55 Mya [43].Contributions of Woodruff 2010 [44] on the biogeography of Southeast Asia highlighted that landscape modifications of this region have resulted in multiple sea-level oscillations, successive biotic compressions, and expansions.Intensive diversification of taxa in the genus Merops was dated to the Pliocene and Pleistocene periods.Many bird species of the rainforests have also differentiated since the late Miocene and mid-Pliocene [45].
The evolution of the genus Merops partially agrees with the speciation process linked to dispersal routes for passerine birds proposed by Johansson et al., 2008 [46].The ancestry of this family might have started in Asia from where bee-eaters dispersed in two directions during the Miocene: One route to Africa may have been through Asia and the second route over Madagascar to continental Africa.
In conclusion, molecular time estimations indicate that the radiation of bee-eaters took place mainly during the Pliocene.Populations of M. apiaster evolved about 0.6 Mya, during a period, which is characterized by climatic oscillations between warm and cold periods in the northern hemisphere.Although European bee-eaters breed from Western Europe to the Far East, they do not show phylogeographic differentiation.Haplotypes show a panmictic pattern, which probably resulted from population mixings in the winter quarters and in refugia during the last glaciation periods.European bee-eaters have spread and now breed in Southern Africa, but they remain migrants.

Figure 1 .
Figure 1.Haplotype network of Merops apiaster using statistical parsimony based on the cytochrome b gene.Colors correspond to the specific origins.

Figure 2 .
Figure 2. Mismatch distribution of frequencies of the pairwise number of differences among Merops apiaster populations.The dotted line represents the observed values in a stable population with a constant population size.The solid line represents the observed frequency of the pairwise mismatch distribution of mtDNA.

Figure 1 .
Figure 1.Haplotype network of Merops apiaster using statistical parsimony based on the cytochrome b gene.Colors correspond to the specific origins.

Figure 1 .
Figure 1.Haplotype network of Merops apiaster using statistical parsimony based on the cytochrome b gene.Colors correspond to the specific origins.

Figure 2 .
Figure 2. Mismatch distribution of frequencies of the pairwise number of differences among Merops apiaster populations.The dotted line represents the observed values in a stable population with a constant population size.The solid line represents the observed frequency of the pairwise mismatch distribution of mtDNA.

Figure 2 .
Figure 2. Mismatch distribution of frequencies of the pairwise number of differences among Merops apiaster populations.The dotted line represents the observed values in a stable population with a constant population size.The solid line represents the observed frequency of the pairwise mismatch distribution of mtDNA.

Figure 3 .
Figure 3. Discriminant analysis of principal components based on six microsatellite loci of Merops apiaster.The dots correspond to the individuals while the ellipses represent the clusters.Colors represent each different cluster.RLP/BW-Rhineland-Palatinate, Baden-Württemberg and Saxony-Anhalt (Germany).

Figure 3 .
Figure 3. Discriminant analysis of principal components based on six microsatellite loci of Merops apiaster.The dots correspond to the individuals while the ellipses represent the clusters.Colors represent each different cluster.RLP/BW-Rhineland-Palatinate, Baden-Württemberg and Saxony-Anhalt (Germany).

Figure 3 .
Figure 3. Discriminant analysis of principal components based on six microsatellite loci of Merops apiaster.The dots correspond to the individuals while the ellipses represent the clusters.Colors represent each different cluster.RLP/BW-Rhineland-Palatinate, Baden-Württemberg and Saxony-Anhalt (Germany).

Figure 4 .
Figure 4. Graphic representation of Bayesian cluster analysis conducted in a structure based on six microsatellite loci performed with a K from 2 to 6. Population names are provided under the graphics.Each individual is represented by one vertical line and the colors point toward the probability of belonging to each of the genetic clusters.

Figure 4 .
Figure 4. Graphic representation of Bayesian cluster analysis conducted in a structure based on six microsatellite loci performed with a K from 2 to 6. Population names are provided under the graphics.Each individual is represented by one vertical line and the colors point toward the probability of belonging to each of the genetic clusters.

Figure 5 .
Figure 5. Phylogenetic relationships and molecular dating of Merops apiaster inferred throughout Bayesian analysis using Cyt b sequences (in total 1059 bp).Posterior probabilities are placed at the nodes.

Figure 5 .
Figure 5. Phylogenetic relationships and molecular dating of Merops apiaster inferred throughout Bayesian analysis using Cyt b sequences (in total 1059 bp).Posterior probabilities are placed at the nodes.

Table 2 .
Haplotype table based on cytochrome b of Merops apiaster.Single haplotypes were not included.N = haplotype frequency.

Table 3 .
Genetic diversity at six microsatellite loci in European bee-eaters.
Weir and Cockerham 1984) of alleles, A r -average allelic richness, H o -observed heterozygosity, H s -gene diversity, H T -overall gene diversity, Fst-genetic differentiation between populations (estimated based onWeir and Cockerham 1984).

Table 4 .
Population size (N), the average number of alleles per locus (Na), observed heterozygosity (Ho) and expected heterozygosity (He) at six microsatellite loci in the European bee-eater.
RLP/BW-Rhineland-Palatinate and Baden-Württemberg (Germany), Scandinavia-Denmark and Sweden, Northern Africa-Egypt and Tunisia, Western Asia-Turkey, Georgia and Kuwait, Northern Asia-Russia and Kazakhstan, Sub-Saharan Africa-Southern Africa and Malawi.

Table 5 .
Pairwise Fst comparison matrix based on the microsatellite data of the European bee-eater.Palatinate and Baden-Württemberg (Germany), Scandinavia-Denmark and Sweden, Western Asia-Turkey, Georgia and Kuwait, Northern Africa-Egypt and Tunisia, Northern Asia-Russia and Kazakhstan, Sub-Saharan Africa-Southern Africa and Malawi.In bold p < 0.05.