Abstract

Evaluation and conservation of local genetic resources of the domestic honey bee populations is important, especially in regions with high diversity levels as well as high honey bee colony density. Greece is rich in honey bee biodiversity, hosting several subspecies, with the status of them being, though, doubtful. The purpose of the present study was to investigate the genetic relationships of both stationary and movable honey bee populations, originating from many location throughout Greece. Two molecular markers were utilized, namely, the conserved mitochondrial COI gene and the highly variable complementary sex determination (CSD) gene. Samples were collected from nine distant populations: eight populations from colonies that followed the traditional stationary beekeeping model and one following the modern migratory beekeeping model type, where the hives are transferred from place to place according to the season. Regardless the beekeeping model, all populations were characterized by sufficient genetic diversity indicating no signs of inbreeding or any bottleneck effects. Nevertheless, genetic differentiation and phylogenetic analysis in comparison with haplotypes obtained from GenBank revealed a genetic admixture pattern suggesting that movement causes genetic homogeneity, occasionally in the stable reared populations as well. Interestingly, two populations, namely, Kastoria and Protokklisi, belonging to A. m. macedonica population, were significantly differentiated, supporting the maintenance of their genetic integrity. Unfortunately, on the other hand, genetic structure of the populations from Crete (Sasalos population) provided evidence that the indigenous breed from the island, A. m. adami, has probably gone extinct. Future management strategies should focus on the conservation of the local genetic resources in which distinct genetic identity has been sustained.

1. Introduction

Honey bees (Apis mellifera; Hymenoptera: Apidae) constitute one of the most widely studied invertebrate groups, considering their critical role in pollination of crops and wild plants, as well as their contribution in highly valuable products [1]. Also, their ethology is of highly interest since they compose complex social structures [2]. Honey bees are characterized by broad distribution globally, while their distribution has been extensively enhanced by human-mediated translocations on account of apiculture purposes [3]. Particularly, several migration and import events have been reported (e.g., [4, 5]). The nonnative A. mellifera introduced in many countries worldwide makes its conservation a priority [6].

Honey bee populations are generally delineated by high intraspecific and intrapopulation diversity, with distinct phenotypic and genetic features documented among different geographic groups [1, 7]. Based on the variations among geographical groups, numerous different lineages and “ecotypes” have been defined within A. mellifera [8]. This intraspecific diversity presents a fascinating parameter for studying the evolutionary biology, ecological dynamics, and genetic bases of complex traits. At least 30 different subspecies of A. mellifera have been described [8, 9], classified into seven evolutionary different lineages, in which East European (C), West European (M), African (A, L and U), and Western Asia (O and Y) honey bees are included [10, 11]. These subspecies are also known as “geographic races” owing to distributions that correspond to different geographic areas [12].

Apiculture development provides significant benefits for the countries located throughout the Mediterranean Basin, favoured by the climatic conditions as well as the wild plant, herb, and crop biodiversity [13, 14]. In Greece, beekeeping constitutes a sector of high importance since ancient times [15]. Greek apiculture is dominated by nonprofessional apiarists, with a significant increase in the colony density recorded year by year, reaching an average of at least 11.45 colonies per km2 [16]. Most nonprofessional apiarists usually possess no more than 50 hives. Furthermore, in the beekeeping sector in Greece, more elderly farmers are mostly employed, whereas for the one-third of the beekeepers, the beekeeping is not their main occupation [17], a fact affecting the economic, ecological, and social organization importance of the apiculture sector in Greece.

A significant number of Greek apiarists follow the modern migratory apiculture practice and move their hives across the country from spring to autumn [18], whereas there is still a considerable number of beekeepers that rely on the traditional model keeping stationary apiaries. The introgressive hybridization which can occur due to migratory beekeeping and to other beekeeping practices such as queens’ importation results to a disturbance of the domestic subspecies’ biodiversity [19, 20]. Moreover, the imported queens may lead to various unfavourable effects such as bee aggression, strong tendency to swarming, and loss of productivity after the first generation [4, 21]. The introduction of nonindigenous genetic resources may also result to the replacement of native populations by nonnative bee subspecies and could lead to local extinctions and eventually to the reduction of bee biodiversity [22]. Furthermore, the hybridization of A. mellifera local populations with allochthonous populations will affect their survival [23]. Comparisons between bee colonies from different European countries and their environment revealed that the local populations are better adapted to specific environmental conditions (predators, diseases, honey flow, and weather) and occasionally produce higher honey quantities compared to the imported ones [23, 24]. Hybridization could also result to the loss of specific characteristic local populations have, as, for example, the resistance A. m. macedonica bees show to Acarapis woodi mite [25].

Hence, conservation of genetic resources of domestic honeybee populations is of high importance, and knowledge concerning the possibility of domestic honey bee colonies to be genetically improved without the reliance on imported genetic material is crucial [5]. Characterization of the various races can lead to appropriate management of bee colonies for better production performance and their adaptation to local environmental features, to eliminate the dependence on allochthonous genetic resources [4]. Keeping these in mind, it is widely accepted that genetic differences and the origin of honeybee populations in the world are of particular interest [5, 9].

Based on morphometric characteristics, Ruttner [9] concluded that in Greece, there are four subspecies of A. mellifera, i.e., A. m. adami in Crete and the Aegean islands, A. m. macedonica in Macedonia and Thrace, A. m. cecropia in Central and Southern Greece, and A. m. carnica in Ionian Islands [26]. Among these, A. m. macedonica, the Macedonian bee, is the dominant honey bee taxon in Greece, mainly due to the high hybridization rates existing from migratory beekeeping and commercial breeding, a phenomenon which is also observed in other European countries as well [19, 2730]. More recent data shows that high hybridization exists in Greece as A. m. ligustica from Italy and A. m. carnica from the former Yugoslavia can be found in many areas [28, 29].

Genetics of A. mellifera populations have been evaluated based on various genetic markers. Mitochondrial DNA (mtDNA) is one of them and has been used widely because of its simple structure, fast evolutionary rate, and solid genetic information mainly for determining the ancestral origin mtDNA [31, 32]. In a honey bee colony, because all bees (workers and drones) share the same maternal DNA, mtDNA analysis is a good marker for determining the genetic origin [3337].

Honey bees’ sex is under the influence of complementary sex determination (CSD) mechanism [38, 39]. Haploids are characterized by a single CSD allele which gives normal males (drones), while diploids, heterozygous for CSD, turn into females. Additionally, when CSD genes are homozygous, it will lead to a diploid male which will be cannibalized by worker bees at the larval stage [40]. As can be drawn from the above, low allelic diversity for this gene will lead to nonviable diploid drones and as a result will reduce colony strength [41]. Thus, this gene is particularly useful as a genetic marker for population genetic studies of honey bees.

The purpose of the present study was to evaluate the genetic relationships of honey bees from different areas throughout Greece, based on the mitochondrial COI and the nuclear CSD genes. More specifically, the aim of the present study was to investigate the genetic diversity and phylogenetic relationships among these populations as they were collected from long-distant localities and maybe belonged to different subspecies or races. Furthermore, the present study will help to provide an answer to the question whether these populations can be discriminated from each other and assess the level of genetic differentiation. Individuals analyzed in this study originated from stationary apiaries, and we expected to find new haplotypes and to demonstrate that there are still genetically distinct populations in Greece which have to be protected and preserved.

2. Materials and Methods

2.1. Sample Collection and DNA Extraction

Worker bees were collected from nine different regions all over Greece. The 8 regions represented stationary apiaries without any introductions of imported queens, whereas the ninth region, Larissa (number 5 in Figure 1), represented by an apiary with imported genetic material belongs to A. m. carnica (Figure 1 and Table 1). Hence, all the apiaries except one, i.e., Larissa, keep the hives stationary. In all cases, in order to eliminate the possibility of genetically related individuals that may bias the results, all honey bees were collected from different hives in each sampling site. After collection, samples were immersed immediately in 100% ethanol, transferred to the laboratory, and stored at 4°C until DNA extraction.

Before the DNA extraction process, ethanol was totally removed from the samples. Total DNA was extracted from one leg using the NucleoSpin tissue kit (Machinery Nagel, Germany) following the recommendations of the manufacturer. The concentration and the purity of the DNA were evaluated in a Q5000 microvolume spectrophotometer (Quawell, China).

2.2. Molecular Analyses

Extracted DNA from each individual was amplified in separate PCRs using the primer pairs F2 CSD8110/R2 CSD8695 and F1 CSD10984/R1 CSD8657 [42] that target the CSD gene and the primer pair LCO1490/HCO2198 [43] that targets the mitochondrial COI gene. In each 20 μl PCR reaction 10 μl FastGene Taq 2X Ready Mix (NIPPON Genetics, Europe), 0.6 μl of each primer (10 μM), 1 μl of DNA (50 ng/μl), and ultrapure water up to the final volume were included. The conditions of each reaction were 95°C for 3 min, 94°C for 30 sec, annealing temperature for 40 sec, 72°C for 40 sec, and a final amplification step at 72°C for 5 min. The annealing temperatures were 47°C for the COI gene and 51°C for the CSD gene for both primer pairs. After purification of the successfully amplified PCR products using the commercial NucleoSpin Gel and PCR Clean up kit (Machinery-Nagel, Germany) following the manufacturer instructions, the purified products were bidirectionally sequenced applying the Sanger methodology in an ABI 3730xl automatic sequencer. Sequences were aligned in the software MEGAX [44].

2.3. Genetic Diversity Analysis

The aligned DNA segments for both CSD and COI genes were initially analyzed for genetic diversity and genetic differentiation-related parameters in DnaSP [45]. Measurements performed were number of segregating sites (S) haplotype diversity (Hd), average number of nucleotide differences between haplotypes (K), nucleotide diversity (Pi), number of haplotypes excluding and including indel sites h1 and h2, respectively, number of nonsynonymous substitutions per nonsynonymous site (Ka), and number of synonymous substitutions per synonymous site (Ks). The pairwise genetic differentiation between populations was calculated using Fst values.

2.4. Phylogenetic Analysis

Phylogenetic relationships among the novel haplotypes as well as with previously deposited in the GenBank database were evaluated utilizing the MEGAX software. Sequences were aligned using the ClustalW algorithm, and two maximum likelihood (ML) phylogenetic trees were created including sequences retrieved from GenBank with accession numbers for COI gene (MF100917, OP890237, OP763642, OM320635, MT188686, MN714160, AP018404, MF100910, MG298938, KY464957, KU601505, OP890239, OM203307, and OP890236) and for CSD gene (HM588939, HQ915847, KM199798, HQ915823, HQ915857, HQ915780, KY502205, HQ915829, HQ915810, HQ915855, JX915853, JX915858, DQ325102, KY502242, HM588919, HQ622112, KY502245, DQ325082, KY502232, and HQ915790). Only sequences of the same lengths or bigger that were accordingly trimmed were included in the analysis. For ML tree construction with MEGAX, the best-fit substitution models for COI (T93) and CSD (HKY+G+I) phylogenetic trees were determined applying the Akaike information criterion (AIC). Trust in nodes was estimated by 1000 bootstrap repeats. The genetic relationships were visualized in dendrograms created using MEGAX.

Finally, a haplotype network was created to evaluate evolutionary relationships and genetic relatedness of exclusively the Greek haplotypes in CSD, taking into consideration the number of appearances of each haplotype in each population. The median joining network for depicting evolutionary relationships and relatedness was generated in the software PopART [46], using parameters to fit goodness, setting geographic populations.

3. Results

3.1. Polymorphism Analysis

Population genetic diversity indices were calculated for nucleotide data for CSD and COI sequences separately for the nine different populations studied. Concerning specifically CSD sequences, the genetic diversity indices were calculated after the exclusion of indel sites. Genetic variability parameters based on COI and CSD sequences are shown in Tables 2 and 3, respectively.

Seventeen different haplotypes were defined based on the COI sequences obtained, indicating generally low levels of polymorphism in all populations studied. These haplotypes were deposited in the GenBank database, given the accession numbers OR947421-OR947437. The number of segregating sites ranged between 0 and 2, while the average number of nucleotide differences between haplotypes (K) were zero in most of the cases with three exceptions being 0.4 in Ikaria population, 0.6 in Kymi, and 0.93333 in Ioannina. The haplotype diversity (Hd) and nucleotide diversity (Pi) followed the same pattern as K with the values for the same populations being 0.4, 0.6, and 0.733 and 0.00059, 0.00088, and 0.00138, respectively (Table 2). All the above indices confirmed the low genetic diversity in COI gene among the Greek bee populations.

Regarding CSD gene a higher number of haplotypes was defined, reaching 57 in total that were deposited in the GenBank database and given the accession numbers OR972730-OR972776. Furthermore, the number of segregating sites ranged from 3 to 16 while the average number of nucleotide differences between haplotypes ranged among 0.60000 (Ikaria) to 8.80000 (Palamas). The haplotype diversity (Hd) was found high with most populations having one as a value and nucleotide diversity (Pi) followed the same pattern as K with the values ranging from 0.00165 for Ikaria to 0.02424 for Palamas population. All the above indices confirmed the high genetic diversity in CSD gene among the Greek bee populations.

3.2. Genetic Diversity

To evaluate the level of genetic differentiation among populations, studied Fst values were calculated among all pairwise populations. The Fst values were estimated using both genetic markers separately, COI and CSD, and are presented in Tables 4 and 5, respectively.

3.3. Phylogenetic Analysis

Following the results of genetic diversity, we constructed two maximum likelihood phylogenetic trees, one for sequences obtained for each gene (CSD and COI) (Figures 2 and 3, respectively) including some sequences retrieved from GenBank, in order to examine their phylogenetic position as well. The number of sequences retrieved from GenBank were 20 for CSD and 14 for COI. When CSD sequences were analyzed, most of the samples collected from the same location grouped together (Figure 1). Furthermore, genetic distance among most of the international sequences and populations from Larissa is smaller when compared with the genetic distance of the other populations and the cluster with the international sequences.

However, when COI sequences were analyzed, the phylogenetic tree failed to group them according to their origin (Figure 3). A few exceptions noticed, such as population from Protokklisi and Kastoria which seems to have a differentiated genetic identity.

The evolutionary relationships and genetic relatedness of the CSD haplotypes are visualized in the network of Figure 4. It is characterized as a partially star-like haplotype network, with a central widely distributed haplotype in the majority of the populations, surrounded by haplotypes shared by almost all populations. Nevertheless, Protokklisi and Kastoria populations are divided by several mutational steps, grouped outside the star-like pattern. The general pattern of the network, overall, probably indicates, at least partially, the maintenance of genetic integrity of the majority of Greek breeds.

4. Discussion

Introduced species represent an important threat for local ecosystems as they can disrupt the balance and reduce the biodiversity [47]. The increase of hives rearing honey bees outside their native range constitutes a severe factor that negatively influences the biodiversity protection [48].

Greece has been characterized by high honey bee biodiversity hosting four out of the approximately 20 races. In the present study, nine populations, originating from different geographic ranges of honey bee subspecies from Greece, were investigated using a classic conserved mitochondrial DNA marker, as well a highly polymorphic nuclear gene, namely, the COI and CSD genes, respectively. CSD gene, specifically, proved to be a very promising marker in revealing genetic architecture in honey bee populations.

4.1. Genetic Diversity Indices

Animal populations are often exposed to numerous stressful abiotic conditions or to pathogenic microorganisms with harmful impacts, and genetic variation can be beneficial towards their survival [49]. Analysis of CSD sequences revealed that genetic variability values by the means of average number of nucleotide differences between pairwise sequences as well as of total nucleotide diversity per site indicate considerable levels of genetic diversity in all examined populations, within the range observed in previous studies based on the same marker [42, 5054]. These inferences reveal a general pattern of sufficient welfare levels in reared honey bee populations in Greece.

High levels of genetic variability are not uncommon in reared honey bees elsewhere [42, 50, 51, 53, 55, 56]. For instance, the vast majority of honey bee populations from the USA have been characterized by high levels of haplotype diversity, ranging from 0.236 to 0.763 [10]. Similarly, honey bee populations from Iraq indicated considerable genetic diversity, ranging from 0.39 to 0.47 based on COI and 16S rDNA genes, with an average value of 0.44 [57]. Native Western African honey bee populations from Nigeria showed an overall even higher level of genetic diversity at the level of 0.702 [58]. These values are not only in agreement with our results but also demonstrate a well-preserved genetic variation despite the fact that concern reared populations, in which phenomena such as inbreeding and bottleneck effects may result to detrimental consequences for the fitness of populations.

Yu et al. [59] suggested that plant nature reserve creates a stable habitat that protects A. mellifera populations from the bottleneck effect, resulting in high genetic diversity values as well as the emergence of novel haplotypes. Interestingly, in our case, Kastoria population exhibited the highest genetic variation indices. Kastoria is located in close proximity with the Prespa National Park, wherein the greatest plant biodiversity in Greece has been identified. Worth noting that while Greek flora hosts more than 6,000 species and subspecies, consisting of the one-fourth of the Mediterranean plant diversity [60], the wide area of Lake Prespa hosts the one-third of Greek flora [61]. Keeping this in mind, the population of Kastoria, in which the highest genetic diversity was observed, may be affected by the high diversity of plants that in line with the suggestion of Yu et al. [59], favours the maintenance of genetic variation.

On the other hand, previous studies from Greece demonstrated reduced genetic diversity [28, 30] in comparison to the current findings. These contradictory results may be probably attributed to the nature of the genetic markers examined. The COI gene serves as a very useful tool for species identification and phylogeny reconstruction; however, according to our results in accordance with other recent studies (e.g., [42]), we can propose CSD as a by far more appropriate marker for genetic diversity estimations, on account of its high polymorphism level. It should be also noted that the maternal inheritance of mtDNA, a principle for the majority of animal species including honey bees, denotes that all the workers and drones in a colony share the mtDNA haplotype of the queen. Variation in the mtDNA of honey bees has been used to provide insights into their biogeography and, like morphometric data, reveals several main evolutionary lineages of Apis mellifera.

Furthermore, all collected samples except from Larissa originated from beekeepers that follow the traditional beekeeping method, without moving their hives. The latter is reflected in the genetic diversity of the population from Larissa, which was the lowest observed. Thus, it cannot be excluded that in commercial beekeeping, operations for the management of the colonies have resulted in genetic diversity loss. Even based on mitochondrial diversity, native-range population genetic variation is often higher, on the contrary with introduced or populations of moving hives, which occasionally carry the same COI haplotype [62]. As a consequence, the establishment of proper management practices is of major importance for sustaining this genetic diversity, both on a biodiversity point of view and for designing genetic improvement programs.

4.2. Population Structure and Phylogeny

The majority of the investigated populations were classified according to their geographic origin based on the CSD gene analyses. The only exception was the population from Sasalos, Crete, which was admixed with A. m. cecropia instead of A. m. adami that was expected based on the morphological assignment proposed previously. This finding, due to the high conservation priority and importance, is discussed in a separate section, namely, in Local Extinctions and the Case of A. m. adami.

According to Ruttner [9], based on classical morphometric markers, the subspecies inhabiting Greece are A. m. adami in Crete, whereas in Peloponnese, Palamas and Kymi bees belong to A.m. cecropia and in Kastoria and Protokklisi to A. m. macedonica, while in Ioannina and Larissa, both A. m cecropia and A. m. macedonica could be found. Interestingly, among Greek haplotypes geographically assigned in accordance with their origin, a number of haplotypes from various distant origins were interfered, including China, Russia, and the USA. This admixture pattern can be probably attributed either to the introduction of allochthonous queen bees in Greece or less likely to the introduction of Greek indigenous queens to other countries. Genetic relationships of Larissa populations provide evidence that probably the former is the case. Indeed, migratory beekeeping had been observed since ancient times, specifically at certain times of the year, when beekeepers moved their apiaries to areas with new blooms in order to collect larger quantities of honey or to suitable areas for overwintering [63]. More specifically, migratory beekeeping from ancient times until today is practiced all over Greece, carried out over long distances [64].

These findings indicate that especially for the populations of Larissa in which admixture is probably evinced, they probably have not maintained their genetic identity and hence they were genetically closer to various sequences obtained from GenBank. Thus, it is proposed that modern apiculture practices are responsible in a great extent for genetic introgression and loss of genetic structuring of local honey bee populations. The alternative scenario, i.e., the reverse genetic introgression of Greek indigenous honey bees to the Americas or towards Anatolia, could be easily precluded, since although occasionally indigenous honey bees have been introduced worldwide as well, these introductions are rather of smaller scale.

On the other hand, the inferences based on the CSD gene are not completely reflected in the analysis of the COI sequences, where in some populations, phylogenetic analysis failed to reveal the structure according to all different Greek populations. Only in distinct populations, subdivision was observed, particularly in cases where this subdivision was more robust. Previous genetic studies based on 16S, COI, and ND5 genes in Greece revealed that A. m. macedonica subspecies could be discriminated from the other populations studied while it was characterized by the higher genetic distance from other populations [30]. This is also in line with our results where populations from Kastoria and Protokklisi are supposed to have A. m. macedonica subspecies according to morphometric classification by Ruttner [9].

Overall, based on the different results of each of the two markers applied, CSD gene may provide a very good mean for population genetics in hone bees, enlightening the relationships among local populations. This is in agreement with the general assumption, suggesting that particularly in invertebrates, mtDNA-conserved regions such as COI gene serve as appropriate markers for species identification, but instead other genomic regions, usually nuclear, may be more informative for resolution of genetic composition. Although mitochondrial DNA is a very useful tool in phylogeographic and genetic analysis, it is characterized by the absence or low recombination and relatively high evolutionary rate, conserved structure, and maternal inheritance [6567]. This is probably the reason why clustering by mitochondrial DNA did not give a clear clustering by geographic groups.

4.3. Local Extinctions and the Case of A. m. adami

Population genetic studies constitute the gold standard means for the detection of local extinctions either from the wild, or in reared populations as well [68]. The current genetic study did not demonstrate the expected genetic distinction of Cretan honey bees, based on classic morphometric analysis [9]. Keeping in mind that the analyzed population was collected from a stable apiary geographically isolated, in the area of Sasalos, reared for several decades with the traditional technique, our results support and confirm that the subspecies adami has been probably gone extinct.

A. m. adami was considered to be the endemic breed of Crete, named so in honour of Brother Adam, who extensively studied highlighting its specific features. It was classified as a distinct subspecies by Ruttner in 1975, just a few years before it was subjected to great demographic pressure, to the point of near extinction, during the 1980s from the Varroa mite, which in those years invaded the Greek area, and to which it proved to be quite vulnerable. Added to this blow was the uncontrolled importation of swarms of other breeds from mainland Greece (mainly A. m. macedonica and A. m. ligustica) and the inevitable admixture with them, to the extent that many claim that this breed has probably now disappeared or concluding that there seemed to be no pure populations of A. m. adami left on the island [69]. It should be noted, however, that later research examining the phylogenetic relationships among honey bee populations from different geographic regions of Greece proposed the existence of A. m. adami in Kasos, Kythira, and Ikaria [30]. Particularly, in the islands of the Eastern Aegean (Karpathos, Kasos, Rhodes, Kos, Chios, and Lesvos), the closely related populations of honey bees indicated a great similarity with A. m. adami, without being officially classified in it, though. This has been questioned ever since, as the populations of the aforementioned islands are clearly distinguishable from A. m. anatoliaca of the neighbouring Asia Minor coasts. Nevertheless, this was not the case in the population from Ikaria according to our results, which was clustered together with A. m. macedonica, supporting the scenario of genetic admixture in this island as well.

The populations of indigenous honey bees from Crete have been afterwards proposed to be probably composed of introgressed populations by A. m. ligustica and A. m. macedonica [70, 71]. The population investigated in the present study supports this polyhybridization pattern, with honey bee haplotypes from Sasalos mostly genetically introgressed with A. m. ligustica, which is in disagreement with the study by Momeni et al. [71]. Thus, interestingly, the investigated populations from the two islands (Ikaria and Crete) have not been influenced by the same admixture patterns and probably followed separate routes. However, we need to keep in mind that Crete especially is a very large island and the populations studied in different studies originated from different areas of the island.

4.4. Future Conservation Practices Proposed

Modern apiculture is facing the pressure of the increasing consumer demand for honey bee products, leading to the requirement of increased productivity. In order to meet these requirements, beekeepers, based on their reliance that local breeds are of lower productivity, introduce allochthonous queens in an effort to improve the quantity of honey bee products [15]. These introductions may have detrimental effects on the local honey bee populations, ranging from loss of their distinct genetic identity that aids to their adaptation to the local habitats up to the occasional complete extinction of some breeds.

Based on the current investigation, although genetic diversity indices of the examined populations did not indicate any sign of fitness loss or reduction of such bottleneck effects, the genetic admixture patterns observed may set in risk the genetic integrity of honey bee population in Greece. Indeed, the general trend of replacing local populations with improved ones, in combination with the movement of hives, may lead to the loss of the genetically adapted populations. Thus, keeping in mind the urgent need for maintenance of local genetic resources in Greek apiculture, appropriate management strategies should be followed.

More specifically, based on our results, at least two investigated populations are genetically differentiated, namely, Protokklisi and Kastoria, maintaining their genetic distinction. In an attempt to prevent beekeepers from admixing these populations and at the same time support the economic and productive efficiency of local honey bee breeds, alternative options may be promising. Particularly, European Union has established various labels, such as protected geographical indication (PGI) and protected designation of origin (PDO), in order to protect quality and highlight and assure authenticity. Our data demonstrate that at least the two aforementioned populations may fulfil the requirements for these classifications. Nevertheless, more studies are needed, examining various genetic markers and other biological aspects to confirm such an assumption.

5. Conclusions

The samples collected in the present study are covering a vast area of Greece, ranging from Alexandroupoli, Evros in the North-East to Kastoria, and Western Macedonia in the North-West to Crete in the South. The main findings presented here were that all honey bee populations investigated, covering a vast area of Greece, ranging from Evros in the North-East to Kastoria, in the North-West, and up to Crete in the South, are characterized by sufficient genetic diversity indicating no signs of inbreeding or any bottleneck effect. However, a genetic admixture pattern was observed probably due to human interference. Although Kastoria and Protokklisi populations, belonging to A. m. macedonica subspecies, were significantly differentiated, Sasalos population from Crete provided evidence that the indigenous breed from the island, A. m. adami can be considered extant. Sustainable conservation practices as well as development of breeding protocols will contribute to the preservation of genetic heritage of Greek A. mellifera. Apart from the contribution of the present study to future conservation practices, trade control, certification of subspecies can also help towards the authentication of the valuable bee products.

Data Availability

Genomic sequences derived in this study were submitted to the GenBank database and given the accession numbers OR947421-OR947437 and OR972730-OR972776.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

The authors are grateful to the beekeepers for their valuable help and contribution in sample collections. This research is supported by the Administrative Region of Western Macedonia, Greece (Special Account for Research and Funds Project Number: 80753).