Genetic diversity and differentiation among insular honey bee populations in the southwest Indian Ocean likely reflect old geographical isolation and modern introductions

With globalization the Western honey bee has become a nearly cosmopolitan species, but it was originally restricted to the Old World. This renowned model of biodiversity has diverged into five evolutionary lineages and several geographic “subspecies.” If Apis mellifera unicolor is indubitably an African subspecies endemic to Madagascar, its relationship with honey bees from three archipelagos in the southwest Indian Ocean (SWIO) hotspot of biodiversity is misunderstood. We compared recent mtDNA diversity data to an original characterization of the nuclear diversity from honey bees in the Mascarenes and Comoros archipelagos, using 14 microsatellites, but also additional mtDNA tRNALeu-cox2 analysis. Our sampling offers the most comprehensive dataset for the SWIO populations with a total of 3,270 colonies from 10 islands compared with 855 samples from Madagascar, 113 from Africa, and 138 from Europe. Comprehensive mitochondrial screening confirmed that honey bees from La Réunion, Mauritius, and Comoros archipelagos are mainly of African origin (88.1% out of 2,746 colonies) and that coexistence with European lineages occurs only in the Mascarenes. PCA, Bayesian, and genetic differentiation analysis showed that African colonies are not significantly distinct on each island, but have diversified among islands and archipelagos. FST levels progressively decreased in significance from European and African continental populations, to SWIO insular and continental populations, and finally among islands from the same archipelago. Among African populations, Madagascar shared a nuclear background with and was most closely related to SWIO island populations (except Rodrigues). Only Mauritius Island presented clear cytoplasmic disequilibrium and genetic structure characteristic of an admixed population undergoing hybridization, in this case, between A. m. unicolor and A. m. ligustica, A. m. carnica and A. m. mellifera-like individuals. Finally, global genetic clustering analysis helped to better depict the colonization and introduction pattern of honey bee populations in these archipelagos.

Introduction varied from an island to another [37] reaching an exclusive level in the eastern Rodrigues population [39]. These findings confirmed that multiples introductions occurred in the Mascarenes but their effect on La Réunion and Mauritius populations have never been characterized using nuclear markers. Therefore, it is unknown whether these populations are undergoing hybridization. Mitochondrial DNA sequencing indicated that Comoros Archipelago might act as a contact area between the Africa coast and Madagascar, and require to be carefully examined using nuclear markers.
The present study characterizes for the first time, the nuclear genetic diversity of honey bee populations within La Réunion and Mauritius (Mascarenes Archipelago) and populations found in Grande Comore, Mohéli, Anjouan, and Mayotte Islands (Comoros Archipelago). We compared the mtDNA polymorphism with nuclear diversity and structuration help in the detection of ongoing hybridization between subspecies on La Réunion and Mauritius. Since multiples lineages coexist in these islands, processes shaping the genetic diversity are hard to disentangle without a large and diverse dataset, capable of discerning between African and European populations. For that reason, we implemented additional sampling from continental native populations and previous datasets from Madagascar [40], Seychelles [38], and Rodrigues [39], to assess the relationship among insular and continental populations. Using the most comprehensive genetic dataset of A. mellifera in the SWIO, we propose an interpretation of the intraspecific phylogeographic patterns in the three archipelagos.

Honey bee population sampling
Genetic diversity in A. mellifera populations from the SWIO islands was assessed by de novo genotyping of worker honey bees from 2,860 colonies from both insular and continental areas. In order to have a comprehensive understanding of SWIO honey bee phylogeography, a dataset containing 1,528 individuals formerly described in Madagascar, Seychelles, Rodrigues, South Africa, and Italy was also employed [37] ( Table 1). The sample (N = 4,388) covered honey bee populations in the SWIO islands and different habitats throughout Africa and Europe (Fig 1, S1 and S2 Figs).
Sampling efforts focused on previously undescribed populations using nuclear markers from La Réunion, Mauritius, Grande Comore, Mohéli, Anjouan, and Mayotte Islands. To increase the probability of obtaining samples representative of these populations, to the extent possible, colonies were collected in different habitats across each island. For La Réunion and Mauritius, collection sites encompassed urban areas as well as virgin tropical forest in National Parks (Parc National de La Réunion and Black River Gorges National Park, respectively). 15.8% of the known managed honey bee populations in La Réunion (13,000 managed colonies; GDS, 2014 personal communication) and 13.6% of those in Mauritius (2,700 colonies; Jhumun [41]) were sampled. Beekeeping is poorly developed in the Comoros Archipelago, so that honey bee foragers were collected every 5 km, whenever it was possible. For each managed or wild colony, workers were collected at the entrance or inside each colony. Immersion in 95% ethanol immediately killed workers, and they were stored at -20˚C until laboratory processing.
Sampling of honey bee colonies in continental Africa and Europe provided reference native populations for comparison. A total of 113 colonies were sampled at 28 A. m. cecropia). Finally, 12 colonies to a strictly known introduced insular population from Tahiti were collected. Table 1. Location details for populations sampled in islands of the southwest Indian Ocean, Africa, and Europe. N de novo the number of honey bee colonies newly sampled, dataset μsat & tRNA Leu -cox2: number of colonies for which microsatellites and mtDNA tRNA Leu -cox2 dataset were previously described, with associated references in brackets, dataset tRNA Leu -cox2: number of colonies for which only mtDNA (tRNA Leu -cox2 dataset) was previously described, but microsatellites were analyzed only in this study. All maps depicting sampling lovations were generated and using QGIS software [42] and derived from open source layers OSM and the world border country polygon shapefile freely available from https://github.com/petewarden/openheatmap/tree/master/mapfileprocess/test_ data/TM_WORLD_BORDERS-0.3.

Population genetic analysis
Genetic diversity was analyzed at both fine and coarse scales by considering different population levels with at least five individuals per apiary. Intra-island/country genetic variation was estimated using each site as a population unit (398 sites) while inter-insular/continental analysis considered all sites from the same island/country (11 islands and 21 countries). Null allele frequency per locus (A null ) was estimated for each population unit with FREENA [49]. The mean number of alleles (N all ), observed heterozygosity (H obs ), unbiased expected heterozygosity (H nb ), and F IS per population unit were estimated using GENETIX 4.05 [50]. Distributions of alleles within and among insular and continental populations was calculated using ADZE [51]. Allelic richness was computed and tested with R and the PopGenReport package [52] using the rarefaction method (for each island/country with a minimum of 6 diploid individuals). Population unit pairwise F ST s were estimated using FSTAT 2.9.3.2 software [53]. GENEPOP 3.4 was used to test deviation from Hardy-Weinberg (HW) equilibrium and genetic differentiation between population pairs [54]. Regarding comparisons among pairs of insular and continental populations, only French sites were kept as separated populations due i) to divergence of mtDNA and ii) the nuclear background shown by significant F ST values. Principal Component Analysis (PCA) was used to further explore genetic differentiation between populations, using R 3.0.2 software [55] and adegenet 1.3-9.2 package [56]. Genetic structure among populations was additionally investigated using STRUCTURE 2.3.3 [57]. A total of 10 6 simulations using 10 5 burn-in steps and MCMC (Markov Chain Monte Carlo algorithm) steps were run for all samples (N = 4,388) simultaneously, considering a K interval  with ten iterations each. The optimal number of clusters was estimated using the ΔK method of Evanno [58]. In parallel, INSTRUCT software [59], which takes into account inbreeding, was run with the same parameters to confirm results from STRUCTURE. Discriminant Principal Components Analysis (DAPC) was also used to support population structure results [60]. Subsequent runs were performed to evaluate intra-island structure in i) La Réunion sites (N = 2,050), ii) Mauritius (N = 367), iii) islands of the Comoros Archipelago (N = 120) and iv) continental populations (N = 263). STRUCTURE HARVESTER [61], CLUMPP [62] and DISTRUCT 1.1 [63] were used to develop the graphical output.

Results
Full sample details including sampling geo-coordinates, multilocus genotypes using the mtDNA tRNA Leu -cox2 intergenic region, and 14 loci microsatellites are available in S1 Table. Distribution of mitochondrial evolutionary lineages Successful analysis of the tRNA Leu -cox2 intergenic region identified 19 restriction profiles in 4,252 colonies from SWIO, Africa, Europe, and Tahiti honey bee populations. The African lineage was characterized by A1, A4, A6, A11, A14, A16, A48, A49, A50, A51, Z2 and Z7 profiles while the European M lineage was distinguished by M3, M4, M6, M7, M7' and M65. On gels, C1 and C2 profiles were difficult to discriminate (1bp difference), so those were coded as C1/C2, but both indicated the European C lineage.
In the SWIO, the distributions and proportions of honey bee mtDNA evolutionary lineages were similar to those reported previously based upon smaller sampling [37]. In Madagascar, Mahé, Praslin, La Digue, Grande Comore, Mohéli, Anjouan and Mayotte, all colonies had tRNA Leu -cox2 DraI profiles characteristic of the African lineage. The Mascarenes Archipelago was the only region in the SWIO to exhibit three different lineages: i) in La Réunion, 95.5% lineage A, 4.3% C, and 0.3% M; ii) in Mauritius, 54.2% A, 44.7% C and 1.1% M; iii) in Rodrigues, 100% C. Proportions of evolutionary lineages formerly reported for the Mascarenes using tRNA Leu -cox2 sequencing (N = 130 for La Réunion and N = 239 for Mauritius [37]) was largely confirmed, despite a massive difference in sample size (N = 2,050 for La Réunion and N = 367 for Mauritius). As for continental populations, distribution of tRNA Leu -cox2 restriction profiles shifted from A maternal lineages in African colonies to M and C maternal lineages in European colonies ( Table 2).

Nuclear genetic diversity in insular and continental honey bee populations
All samples (N = 4,388) were considered successfully genotyped (with less than 40% missing data). Preliminary analysis showed negligible low null allele frequencies for all insular and continental populations listed in Table 2 (site details are given in S2 Table). The asymptotic trend observed in allele accumulation curves may indicate that the majority of alleles at the 14 studied microsatellites loci were captured in the La Réunion and Mauritius populations (S3 Fig). The mean number of alleles per site (n ! 5) showed that African and Mascarenes honey bees have the highest genetic diversity of the populations tested (Fig 2A). However, when looking more closely at La Réunion and Mauritius populations, the mean number of alleles seems to be related to site sampling size (S4 Fig). Allelic richness using rarefaction was more appropriate to confirm this pattern due to the sample size difference. The SWIO insular population allelic richness was not significantly different from continental European populations, but was significantly lower than continental African ones. This lower allelic richness was likewise observed for another insular population in São Tomé (Ar = 2.94 ± 0.71). In the Mascarenes Archipelago, only La Réunion (Ar = 3.18 ± 0.71) and Rodrigues (Ar = 3.09 ± 0.67) displayed allelic richness levels comparable to those of the Cameroon and Malawi populations (Fig 3).

One genetic cluster exists in La Ré union, but Mauritius has two
In La Réunion, low pairwise F ST values between sites (from 0.004 to 0.034) were detected, but only 48 pairs out of 3,655 (considering n ! 5) were significantly different (P < 0.00037). Among these significant pairs, 25 involved two sites from remote environments located in a geological cirque (REU126 and 124). Bayesian clustering (STRUCTURE and INSTRUCT) and multivariate methods detected no genetic structuration among honey bee samples from La When comparing nuclear diversity to mtDNA and location data, one might see that structuration was most likely linked to the presence of African or European lineage in each sampling site from Mauritius ( Fig 4B). The pattern also indicated presence of "hybrid" as some individuals presenting African lineage mtDNA haplotypes were assigned (through a gradient of probabilities) to the same nuclear cluster as European C lineage individuals. The reciprocal situation was also true. Considering a probability threshold of 50%, 70.3% of the colonies from Mauritius were assigned to cluster 1 and the rest of cluster 2 ( Fig 4B). Among the 109 colonies assigned to cluster 2, i) 105 showed mtDNA tRNA Leu -cox2 characteristic of the European C lineage and ii) 83 of these colonies belonged to sites MUS20 and 21 (same beekeeper, Fig 4B).

Honey bee populations from the Comoros Archipelago are structured by island
Pairwise differentiation test using F ST values within Grande Comore, Mohéli, Anjouan, and Mayotte islands could not be performed due to limited sample sizes. Nonetheless, differentiation indices among islands of the archipelago were low, but significantly differentiated (Table 3). Mohéli was not significantly differentiated from all other populations, which might reflect the smaller sample size in the SWIO (N = 11) associated with the highest null allele frequency (A null = 8.3%). Genetic differentiation among the Comoros islands was even visible with the distinction of four groups on PCA (Fig 5) supported by best model estimated by Bayesian clustering method (ΔK4 = 532.1, Fig 4C). At K = 2, individuals from Grande Comore  Genetic diversity and structure of honey bees from the southwest Indian Ocean islands Global dataset structure analysis showed genetic distinction between SWIO and outgroup populations The complexity of honey bee phylogeography in the SWIO was represented in the Fig 6 by STRUCTURE bar plots and distribution of the detected genetic clusters. For the Bayesian method, the uppermost value of K (K = 2) was not coherent with results of other analyses such as PCA (Fig 5B), DAPC, or with STRUCTURE runs computed on each island. To avoid interpretation problems, the choice of K followed recommendation for such case [66] by considering the Ln(K) standard deviation (S9 Fig). DAPC method gave similar results that are available in supplemental S10 and S11 Figs. Here, we concluded that K = 5 is the most appropriate model based upon the present sampling.
Considering a probability threshold of 70%, at K = 5, colonies sharing same African mtDNA haplotypes in Madagascar, La Réunion, Mauritius, and Seychelles were assigned to different genetic clusters (Fig 6B). At K = 5, Mascarenes Archipelago colonies were grouped by island i) 84.0% of La Réunion in red cluster 3, ii) 54.8% of Mauritius in green cluster 5 and iii) 99.6% of Rodrigues colonies assigned to blue cluster 1. As SWIO reference, 99.4% of Madagascar colonies clustered to yellow cluster 2. Mauritius shared a nuclear genetic diversity (alleles) with Seychelles colonies mostly assigned to the same cluster (96.4% to green cluster 5 for Mahé, Praslin, La Digue confounded). Compared to all other populations, Mauritius presented the most mixed population with four clusters (Fig 6). Global STRUCTURE analysis showed Genetic diversity and structure of honey bees from the southwest Indian Ocean islands that one Mauritius cluster is highly similar to individuals from Europe and Rodrigues. Finally, the Comoros Archipelago presented an interesting mixture of cluster 4 (pink) derived mainly from continental Africa and cluster 2 mainly found in Madagascar. A geographic paradox appeared as 100% of the colonies from the most distant island to Madagascar, Grande Comore, was assigned to the same cluster. This SWIO island structure was also supported by significant differentiation estimated by pairwise significant F ST values, but still lower than those of the continents (Table 3). In addition, PCA (Fig 5) showed that SWIO island point scatters were all distinct apart for the case of Mauritius overlapping with several populations.
The global population structure view and the progressive increment of K allow to show multiple conclusions: i) Madagascar honey bees were genetically distinct from continental African colonies, ii) La Réunion bees differentiated early from Madagascar genetic group, iii) Mauritius and Seychelles bees share a genetic background with La Réunion, but were still differentiated, iv) the second intra-cluster found on Mauritius was associated to one European cluster, v) Comoros Archipelago bees share a genetic background with continental African populations and Madagascar, but still formed a different population, and vi) Rodrigues shares a genetic background with European populations, yet formed a new population, unlike Tahiti.

Discussion
This study shows that honey bee (A. mellifera) populations established in the Mascarenes and Comoros archipelagos present high nuclear microsatellite polymorphism for insular populations. Despite having generally similar African mtDNA backgrounds, different levels of nuclear genetic differentiation were apparent within and between archipelagos. No genetic structure was found at La Réunion, while in Mauritius, genetic data indicated the coexistence of two clusters. Both populations have a close relationship to Madagascar honey bees rather than to native African continental populations. Ongoing hybridization between African and European lineages was evident in Mauritius.
In the SWIO, several insular populations occur, and each one possesses a singular genetic pattern in terms of evolutionary lineages and nuclear genetic diversity and structure. Such a complex evolutionary pattern will need to be broken down for each island, in the future.

The Comoros Archipelago as a contact area between Africa and Madagascar
Comoros Archipelago honey bee colonies analyzed are exclusively descended from African lineages, and more particularly the mitochondrial A I sub-lineage observed both in Madagascar and the main part of Africa [15,37,40]. The previous sequencing of two mtDNA non-coding and coding regions showed that Comoros honey bees mainly share a common haplotype with A. m. unicolor in Madagascar, but also exhibit a private haplotype [37]. These colonies also display shared nuclear alleles with Madagascar and African continental populations, putting them in an intermediate position. The Comoros Archipelago is equidistant between Africa and Madagascar (300 km) and could represent an exchange pathway between landmasses.
Each island possesses specific genetic clusters clearly differentiated from the closest neighboring island. Such a pattern is not surprising for insular honey bee populations, as it has been previously reported for the Seychelles [38] and other archipelagos, such as the Canary Islands [67], Madeira Islands, and the Azores [68,69]. The four Comoros islands are separated by an oceanic barrier of approximately 40 km (Grande Comore-Mohéli) to 190 km (Grande Comore-Mayotte), which may have been sufficient to restrict gene flow among the islands. Genetic structure within the archipelago also suggests no recent human-mediated introductions of honey bees from neighboring islands. Unlike the Seychelles Archipelago, beekeeping is not developed in the Comoros Archipelago, so human-mediated exchanges are unlikely to have occurred among the four islands. Paradoxically, Grande Comore, the closest island to the East African coast, was the most genetically similar to Madagascar. The other three islands indicate progressive colonization and gene flow from the African coast and Madagascar, likely via island hopping, as observed for chameleons [70]. On Grande Comore, the unexpected highest proportion of individuals assigned to the Madagascar cluster, probably indicates a recent unreported human introduction.
The African honey bee genetic background predominates in La Réunion. La Réunion was dominated by A. m. unicolor colonies (95.2%), but mtDNA sequences characteristic of A. m. mellifera, A. m. carnica, and A. m. ligustica were also found. However, if multivariate analysis indicated that this population was similar to, but distinct from Madagascar, it also showed a certain proximity to European colonies. The presence of European honey bee colonies confirmed previous reports stating that introduction of exotic subspecies occurred before import prohibitions established in 1982 [71]. Unlike Rodrigues or Mauritius, from the same archipelago, no honey bees with European mtDNA differed from African colonies in their nuclear identities. Conversely no individual with African mtDNA clustered with European colonies on the basis of nuclear genetic diversity. Such patterns could result from low importation levels compared to the large African pool preexistent in La Réunion. In addition, asymmetric introgression in favor of African over European lineages is a known phenomenon in honey bees defined as Africanized [72][73][74][75][76]. Differences in sperm-related genes were proposed to explain this reproductive advantage of African vs. European honey bees [12].
Despite the topology of La Réunion, reaching an elevation of 3069 m [77], no genetic structure was observed across the island. This indicated that gene flow is sufficient to maintain genetic homogeneity. Exchanges between distant locations may have been facilitated by intraisland beekeeping movements following resource cycles. Beekeepers from La Réunion move their hives two or three times per year from the lowlands to the highlands, and from the humid east coast to the drier west coast, following the availability of flowering plants. Homogenization of the genetic pool by such practices has already been described in continental honey bee populations [78].
Interestingly, La Réunion has the highest proportions of African colonies and undegraded native forest habitats, estimated at 25% of the original cover [77]. This difference in environmental conditions compared to Mauritius (2% original forest cover) and Rodrigues (0%) has undoubtedly influenced the established populations. Natural selective pressures exerted on La Réunion colonies may be advantageous to tropically adapted subspecies, such as A. m. unicolor, more than to temperate, introduced subspecies. However, it would be interesting to investigate this in wider genomic and coding regions to determine whether this adaptation is due only to the African genome. As an example, it has recently been shown that a positive selection signal appeared in Chromosome 11 with an excess of European over African ancestry which may confer an adaptive advantage to Africanized honey bee populations found in Brazil [79].
Admixed population and recent hybridization undergoing in Mauritius. For Mauritius, the genetic composition of honey bee populations was the most complex case among all the SWIO islands, highlighting the impact of human honey bee importation. In our sampling, similar proportions of local and exotic lineages were detected. However, one apiary highly contributed to the European frequency as a beekeeper imported a huge batch of A. m. ligusticalike queens [41]. As in La Réunion, European lineages present on Mauritius corroborates records stating that regular introductions were carried out [41,80]. Nuclear analysis uncovered on-going hybridization on this island between African and European genetic groups. Mixing of these two divergent lineages was revealed by cytoplasmic disequilibrium as several individuals with A. m. unicolor mtDNA haplotypes were poorly differentiated and were assigned to the same genetic cluster as native European colonies, and the reverse occurred. In both mtDNA African and European colonies found in Mauritius, a continuum in nuclear assignment toward the opposite cluster suggested ongoing introgression in both directions. Such phenomena have already been reported on other European islands and between the African and European C lineage in A. m. iberiensis [69] and A. m. siciliana [81].
As at Rodrigues, C and M colonies were genetically similar to continental populations of A. m. ligustica and carnica. In Mauritius, two apiaries (MUS20 and 21) showed "pure" European colonies due to relatively recent introductions confirmed by field investigation and government reports [41]. Whether these colonies originated directly from native continental areas in Europe or were introduced from another exotic population could only be determined with wider sampling. Heavy deforestation, leaving only relictual native vegetation in Mauritius, and replacement with high proportions of exotic plants may have also influenced survival of European colonies [82,83].
Genetic diversity and structure in the SWIO reflects island effects and influences of beekeeping practices. Geographic and climatic barriers have played an important role in the evolution of the Western honey bee (A. mellifera) into five lineages with up to 31 subspecies [13,84,85]. The divergence signal can be weak in the case of landscape continuum as for A. m. iberiensis [86,87] or sharper, as in the genetic differentiation of A. m. unicolor from Madagascar [40]. However, a comprehensive picture of indigenous A. mellifera population histories may be difficult to develop due to the long relationship with humans and modern global transportation of honey bee colonies [20,21]. The occurrence of a novel mitochondrial SWIO African sub-group and private haplotypes in the Mascarenes (except Rodrigues), Comoros and Seychelles suggested ancient colonization events [37]. The larger screening of the tRNA Leu -cox2 intergenic region done here in the SWIO populations confirmed that these insular populations are mainly of African origin. Since mtDNA indicates only maternal lineage, the nuclear genetic diversity observed confirmed that the SWIO islands (excepted Rodrigues) are more closely related to Madagascar than to any African populations. Nevertheless, the differentiation index (F ST ), Bayesian, and multivariate analyses showed that SWIO populations present nuclear genetic diversity distinct from Madagascar. Yet this observed divergence between Madagascar and others SWIO islands was still lower than that observed between continental African A and European M and C populations.
Nuclear genetic structuration among SWIO populations could be explained by several nonexclusive hypotheses. First, to obtain such genetic differentiation levels among islands, geographic isolation by ocean barriers should have been maintained long enough to restrict the gene flow, resulting in homogenization. By comparison, the local honey bee from the Balearic Archipelago experienced introduction probably in the XVIII and XIX centuries, and the F ST values among populations (0.04-0.27 [88]) were of the same order of magnitude as in the Mascarenes (F ST = 0.04 to 0.19) or the Comoros (F ST = 0.14 to 0.29). Evidence of gene flow between La Réunion, Mauritius, and the Seychelles Archipelago is highlighted by global Bayesian clustering analysis. This may have resulted from past connections via the exchange of colonies, possibly during the colonial period through the route to India in the XVIII century.
Secondly, genetic differentiation among islands could be due only to divergent lineages assembly. Nevertheless, in native continental areas, F ST values among populations were !2x higher than among SWIO islands with different evolutionary lineages coexistence. In comparison, Rodrigues, which is likely from a European genetic background (mtDNA and microsatellites), showed differentiation values higher than the admixed population of Mauritius, La Réunion, or any other islands from the Seychelles or Comoros archipelagos. This could mean that the genetic diversity pool created by European colonies brought to La Réunion and Mauritius might not fully explain their differentiation from other local populations. If not European, the genetic differentiation may have resulted from colonization by or introduction of African lineages. However, Madagascar was the closest native African population the SWIO islands and none of the continental African populations showed similar proximity (congruence of all analyses: PCA, STRUCTURE, and genetic differentiation test).
At the time when A. mellifera started to diverge throughout its native range, all islands of the SWIO were completely formed, and the Mascarenes were already colonized by angiosperms. For example, the Dombeyoideae family (~25 to 35 Ma [89]) and the Acacia heterophylla species, visited by honey bees now, are believed to have reached the island around 1.4 Ma ago [90] from Hawaii. Consequently, these islands already possessed habitats suitable for generalist pollinators that require pollen and nectar for survival. Madagascar has been identified several times as a base for colonization and radiation into neighboring archipelagos, whether for flora [91][92][93] or fauna [94][95][96][97]. The relatively lower genetic diversity observed in the 11 islands of the SWIO compared to African continental populations suggested progressive colonization [9]. A loss of genetic diversity could be associated with founder events where only a sub-sample of African diversity reached Madagascar. After diversification on that island, a similar evolutionary process likely occurred in the Mascarene, Seychelles, and Comoros archipelagos. Despite being the potential source population for SWIO islands, Madagascar had lower levels of heterozygosity and alleles per locus than any other island. One hypothesis is that Madagascar colonies experienced bottleneck events, possibly due to loss and fragmentation of original habitats caused by deforestation [98]. A similar observation was made on African honey bee colonies, deforestation being identified as the major threat to wild African colonies [99]. A second possible explanation is that part of the high genetic diversity observed in the SWIO archipelago in regard to Madagascar is the result of admixture, which reduces the negative effects of a bottleneck [100] and has been shown to increase diversity levels in honey bee populations [23]. In all the SWIO archipelagos, hybridization occurred or is ongoing with dissimilar assemblies African A I -Malagasy in Comoros, African Z-Malagasy in Seychelles [38], and African-European in the Mascarenes. All these elements, combined with significant population differentiation among SWIO insular honey bee populations, seem consistent with "natural" colonization. Yet, this requires further investigation, as several factors are unknown for each of these populations, complicating interpretation.

Conclusions
Genetic diversity and structure of honey bee populations of SWIO islands suggests ancient colonization events of A. m. unicolor from Madagascar to the Mascarenes and Seychelles archipelagos, old enough to detect population differentiation within the sub-lineage. The use of nuclear and mitochondrial markers uncovered the presence of exotic subspecies and different levels of hybridization with indigenous populations in the archipelagos. The numerous interactions recorded between A. mellifera and endemic species [33,34,36,[101][102][103][104][105] with some remarkably benefits [106], stress the importance of preserving this species.
Apart from its ecological role, these populations with singular genetic diversity deserve particular attention, especially against the global loss of honey bee colonies [107][108][109]. Now that whole-genome sequencing has become more affordable, it would be interesting to investigate the effects of hybridization between African and European lineages in the SWIO islands using a genomic approach. Such data could offer better resolution for estimating times of divergence and would allow us to better retrace the demographic history of these insular populations.  Table. Mitochondrial COI-COII intergenic region diversity (based on DraI restriction profiles) and nuclear diversity indices for each SWIO, African, and European sampling site. N: number of colonies per site; NA COI-COII : number of individuals with missing COI--COII data, N all : mean number of alleles; H nb and H obs : unbiased expected and observed heterozygosity, respectively; F IS ( Ã significant at P < 0.05) and A null : mean allele null frequency. (DOCX) S3 Table. Pairwise F ST values among sites from Zimbabwe (ZWE), France (FRA), and Italy (ITA) with N ! 5, based on 14 microsatellites. After Bonferroni corrections, permutations tests were only significant among French sites (in bold P < 0.000549). Colors as in Table 3. (DOCX) S4 Table. Pairwise F ST values among sites at Mauritius Island with N ! 5 based on 14 microsatellites. Statistical significance for the permutation tests after Bonferroni corrections is indicated in bold (P < 0.000476). Colors as in Table 3. (DOCX) study. We warmly thank both entomological teams of Mauritius and Rodrigues and the Seychelles Agricultural Agency for helping in sampling in the Seychelles archipelago. We wish to thank Steve D. Aird, OIST technical editor, for his time and detailed comments on an earlier version of the manuscript. We thank the editor and two reviewers for their thorough review and constructive comments, which helped us to improve the manuscript.