ABSTRACT
Genetic diversity of four populations of honey bee, Apis mellifera from two vegetation zones in Nigeria namely; southwest rainforest and northcentral derived savanna was analysed using fifteen morphometric characters and five microsatellite loci. Discriminant function analysis (DFA) of the morphometric data revealed a considerable variation of morphological characters between the sampled localities while Principal Component Analysis (PCA) and Canonical Variate Analysis (CVA) produced overlapping clusters of the populations sampled indicating lack of separation between the various populations. The genetic diversity (FST) revealed low differentiation among populations suggesting that geographic distance was not an impediment to gene flow among populations. The overall FIT value indicated that the four populations have a deficiency of heterozygotes suggesting the presence of inbreeding within populations. The analysis of molecular variance (AMOVA) showed that 91% of the total molecular variance existed within the populations while 9% existed among populations, indicating low inter population genetic variation. It is suggested that there is an apparent loss of genetic diversity in the populations of A. mellifera studied in the two vegetation zones of Nigeria. This could have implication for the health and stability of these bee populations.
Key words: Biodiversity, genetic diversity, honey bees, morphometry, microsatellite, population genetics.
The true honey bee, Apis mellifera L. is known to be one of the most economically valuable insects because of its honey production and pollinating activities (Lawal and Banjo, 2010). The services of bees and other pollinators to agriculture is estimated to be between $235- $577 billion per year worldwide (FAO, 2018). The A. mellifera originated from Africa (Whitfield et al., 2006), and is naturally distributed to Europe and Asia (Howpage, 1991; Nedic et al., 2011). It was introduced into America and Australia by humans (Tunca and Kence, 2011). The species is found on every continent except Antarctica, that is, all the habitats on the planet that contain insect-pollinated flowering plants. Approximately forty-three (43) subspecies based on geographic variations are recognized (Engel, 1999). The subspecies are divided into four major branches, based on work by Ruttner (1988) and confirmed by mitochondrial DNA analysis (Smith, 1991; Garnery et al, 1992; Palmer et al., 2000) and Single Nucleotide Polymorphism (Whitfield et al., 2006).
Recently, honey bee populations all over the world have been reported to be on the decline due to Colony Collapse Disorder (CCD). Akinwande et al. (2013) reported a decline in the number of honey bee colonies from selected apiaries in Southwestern Nigeria. There were 58.34, 44.84 and 40.61 average percentage declines in colony establishment in Lagos, Ogun and Osun States respectively. The presence of pests and diseases, pesticide poisoning, lack of queen rearing, poor hive and seasonal management were suggested to be the major factors responsible for the annual decline in honey bee colony establishment. Oyerinde and Ande (2009) reported that 15.01% of the 2000 installed bee hives in Kwara State, Nigeria had established bee colonies. Some biotic factors were suggested to have been responsible for the rather low bee colonizing record of the state. There is little information about the genetic diversity of honey bee populations in Nigeria. This study provides information on the genetic diversity within and between populations of A. mellifera which will be useful for taxonomic re-evaluation of the species for subsequent conservation efforts.
Sampling
Samples of workers of A. mellifera were collected randomly from 28 colonies in two vegetation zones in Nigeria namely: the Tropical Rainforest (Edunabon in Osun State and Ijebu-Ode in Ogun State) and the Derived Savanna (Kabba in Kogi State and Malete in Kwara State). Figure 1 is the map of the study area showing the sampling sites with the map of Nigeria inset. Sampling was carried out in apiaries which do not practice migratory beekeeping, and the hives sampled were stationary during the sampling period. Seventy specimens from each location were sacrificed in ether vapor, and then preserved in 90% ethanol for morphometric studies. Seven specimens from each location were preserved in sample bottles containing 90% ethanol for genomic DNA extraction. Identification was done with a dissecting binocular microscope, using bee identification keys of Michener (2007).
Morphometric analysis
A total of 28 colonies were subjected to morphometric analysis. Ten worker bees per colony were dissected and measured for 15 morphometric characters according to Ruttner et al. (1978) using a dissecting binocular microscope and vernier calliper. The details of the characters measured are shown in Table 1.
Microsatellites analysis
Genomic DNA was extracted from the thorax of seven (7) worker bees per population using the CTAB (Cetyl Trimethyl Ammonium Bromide) method. The isolated DNA was analysed in a thermocycler using five microsatellite primers (A024, A028, A043, A088 and A113) selected according to Genebank which had been previously used by Franck et al. (2001). Polymerase chain reactions were carried out in standardized 10 μL reaction mixture containing 0.1 mM of each deoxyribonucleoside triphosphate (dATP, dCTP, dGTP, dTTP), gel loading buffer, stabilizers, 0.3 µl each of the forward and reverse primer sets, 1.5 mMMgCl2, (NH4)2SO4, 0.5 µl Taq polymerase in 1X buffer, 5.9 µl of PCR grade water, and 20 ng total genomic DNA. The mixture was incubated in a GeneAmp PCR thermocycler programmed as follows: 30 cycles each of denaturation, annealing and extension temperature at 95°C for 20 s, 58-60°C for 25 s and 72°C for 45 s respectively and a primer extension temperature of 72°C for 60 s followed by final extension temperature at 72°C for 10 min to complete the amplifications. The amplicons generated were then subjected to electrophoresis on 1.4% Agarose gel and visualized by staining with ethidium bromide.
Statistical analysis
Principal Component Analysis (PCA), Canonical Variate Analysis (CVA), Discriminant Function Analysis (DFA) and Cluster Analysis on morphometric data of the honeybee populations were performed using the software, PAST (Hammer et al., 2006). Data generated from microsatellite studies were analysed using GenAlex 6.502 Software (Peakall and Smouse, 2006, 2012). The total number of alleles, allele frequencies, average number of alleles per locus, observed (Ho) and expected heterozygosity (He) for each population across the loci, were estimated. Analysis of Molecular Variance (AMOVA), fixation indices (FST, FIT and FIS), degrees of heterozygosity and polymorphism, mean gene flow and Hardy-Weinberg Equilibrium (Nei, 1978) were also estimated. Phylogenetic tree was constructed using PHYLIP-3.695 (Felsentein, 2014).
Morphometric studies
The average values, range and standard deviation of all the morphometric characters measured are shown in Table 2. Wide range of sizes (Standard deviation, SD) especially on proboscis length (PL) and right forewing length (RFL) were found in samples from Edunabon and Ijebu-Ode.
Principal Component Analysis (PCA) of the 15 morphometric measurements of A. mellfera from the four study areas (Figure 2) showed overlapping of all the four clusters produced. The CVA plot (Figure 3) showed overlapping of clusters of specimens from the different populations studied. Figure 4 shows the respective morphometric characters and their loadings on PC1, which indicated that proboscis length is the main characteristic responsible for variation among the four populations of A. mellifera studied (loading, 0.5294), while right forewing length has the second (0.4913) heaviest loading. Proboscises were shown to be longer in the Rainforest than in the Derived Savanna. Similarly, size of tibia and metatarsus were longer in the Rainforest than in the Derived Savanna. Inter-locality variations in wing characteristics were observed which have also been reported by Sharma (1983), and Tahmasebi et al. (2002).
Discriminant Function Analysis (DFA) showed no significant difference between specimens from Edunabon and Ijebu Ode with individuals from both locations overlapping along the discriminant function plot. Also, only 84.3% of the specimens concurred with their a priori classification showing that the discriminant function did not recognize any significant difference among the specimens based on their locations. Significant differences were observed between specimens collected from Edunabon and Kabba, Edunabon and Malete, Ijebu-Ode and Kabba, Ijebu-Ode and Malete, and Kabba and Malete with individuals from each paired location clearly separated along the discriminant function plots. Moreover, 98.57, 99.29, 97.14, 93.57 and 98.57% of the specimens concurred with their a priori classification. Similarly, DFA showed significant difference between specimens collected from the Rainforest and Derived Savanna zones of Nigeria with individuals from both zones clearly separated along the discriminant function plot. Moreover, most of the specimens (95.36%) concurred with their a priori classification. The Unweighted Pair Group Method with Arithmetic mean (UPGMA) dendrogram revealing the clustering pattern of the specimens of A. mellifera across the two vegetation zones is shown in Figure 5.





Microsatellite DNA studies
The mean heterozygosity for the samples across all loci was 50% (Table 3) while the average number of alleles observed (Na) for the total population was 3.450. Unbiased expected heterozygosity (uHe) ranged from 0.830 to 0.997 with a mean value of 0.902±0.118. With a 95% threshold, the percentage of polymorphic loci at the population level was 100%. Observed heterozygosity (Ho) per locus had values of 0.500 each while the expected heterozygosity (He) per locus ranged from 0.770 in locus A043 to 0.836 in locus A024 (Table 3). Shannon’s Information Index (I) ranged from 0.624 in locus A024 to 0.862 in locus A028 with an average value of 0.731±0.072.
Genetic diversity parameters based on allelic frequencies are also shown in Table 3. Observed heterozygosity (Ho) per population had values of 0.500 each while the expected heterozygosity (He) per population ranged from 0.766±0.080 for population 4 to 0.838±0.140 for population 3. The mean Ho and mean He were 0.500±0.000 and 0.800±0.146 respectively. The heterozygosity level within a subpopulation (FIS), the heterozygosity level in total populations (FIT) and the degree of genetic differentiation of subpopulations (FST) are presented in Table 4. All 5 loci illustrated deficiency of heterozygotes in the four populations. The mean FIT amounted to 0.407 ± 0.103 (from 0.347 to 0.494) and the mean FIS across loci was 0.375 ± 0.103 (from 0.351 to 0.402). The fixation coefficients of subpopulations for the loci studied within the total populations, measured as an FST value, varied from 0.030 (A088) to 0.062 (A028), with a mean value of 0.049 ± 0.241. This signified that 4.9% of the total diversity existed among populations while the remaining 95.1% existed within populations. The mean gene flow (Nm) among populations which gives information about genetic divergence or genetic similarity of subpopulations due to gene flow was 5.340 ± 0.432. In other words, gene exchange between populations was low. Analysis of Molecular Variance (AMOVA) showed that 91 % and 9% of the total molecular variance was within and among populations respectively (Table 5). This implies that the populations were not significantly different from each other.


A summary of the test for departure from Hardy-Weinberg (H-W) equilibrium across loci and populations showed that at P <0.05, three (A024, A028 and A113) out of the five loci (60%) studied in the Rainforest and four (A024, A028, A043 and A113) for the Derived Savanna zones (80%) were in H-W equilibrium. The Chi-square (χ2) test (P<0.05) indicated significant departures from H-W equilibrium in many cases (55%). All the deviations were primarily attributed to heterozygotes deficit.
The UPGMA cluster analysis based on Nei’s unbiased genetic distances (GD) is shown in Figure 6. The dendrogram separates the four populations into two (2) major clusters with three sub clusters. The first cluster consists of only population 1 while the other cluster consists of populations 2, 3 and 4. Within the second cluster, populations 2 and 4 are clustered together.
The mean values of the set of morphometric characters measured agrees with those reported for the subspecies of A. mellifera in Nigeria (Dukku, 2016) and sub-Saharan Africa (Ruttner, 1988; Yu et al., 2012), though not in agreement with those reported by Oyerinde et al. (2012) and Ajao et al. (2014), whose values fell outside the range reported for all subspecies of A. mellifera. This agreement validates the correctness of the measurements taken in this study. The differences in measured wings characteristics could result from the use of the wings for flight during foraging and thermal regulation of comb. The length of the proboscis was considered a very important character because it showed the geographical variability more than all the other studied characters (Marghitas et al., 2008). Inter-locality variations of the proboscis, tibia and metatarsus are in line with Allen’s rule: appendages of the body relatively shorter in the North than in the South (Ruttner, 1988). Although gradual variation was established along the Rainforest-Derived Savanna continuum, no morphometric differentiation has yet been found, in spite of the geographic distance and prominent differences in humidity and altitude.
The Inbreeding coefficient (FIS) value indicated that overall, the four populations had heterozygotes deficit suggesting the presence of inbreeding within populations which could lead to subsequent loss of unexploited genetic potential. The mean gene flow (Nm) among populations, which gives information about genetic divergence or genetic similarity of subpopulations due to gene flow, indicated that there was small genetic differentiation among the populations. In other words, gene exchange between populations was low. AMOVA revealed that most of the variability (91%) was observed in individuals within populations. Measurements of genetic distance (GD) revealed that the Rainforest populations were genetically more diverse (0.600) than the Derived Savanna populations (0.021).
Analysis of the genetic diversity of the A. mellifera populations suggests a possible loss of variability. This loss could be attributed to inbreeding depression and/or any of population restructuring, loss of habitat through deforestation, hunting for honey involving killing of wild colonies, natural selection, genetic drift and introduction of exotic honeybees and the parasitic mite, Varroa destructor (Akinwande et al., 2013).
The result of this study reveals that A. mellifera populations studied are morphometrically similar. There is a need to maintain a healthy level of genetic variability in A. mellifera populations, therefore efforts should be made to protect bees from the threats to their abundance, diversity and health. This may be achieved by monitoring and curtailing the effect of inbreeding depression, population restructuring, deforestation, poaching, natural selection, genetic drift and introduction of exotic honeybees and the parasitic mite, V. destructor. Also, enforcement of legislations aimed at protecting honeybees, in particular, and the ecosystem in general should be put in place.
The authors have not declared any conflict of interests.
REFERENCES
Ajao AM, Oladimeji YU, Idowu AB, Babatunde SK, Obembe A (2014). Morphological characteristics of Apis mellifera L. (Hymenoptera: Apidae) in Kwara State, Nigeria. International Journal of Agricultural Sciences 4(4):171-175.
|
|
Akinwande KL, Badejo MA, Ogbogu SS (2013). Challenges Associated with the Honey Bee (Apis Mellifera Adansonii) Colonies Establishment in South Western Nigeria. African Journal of Food, Agriculture, Nutrition and Development 13(2):7467-7484.
|
|
|
Dukku HU (2016). Evaluation of morphometric characters of honeybee (Apis mellifera L.) populations in Nigeria. PeerJ PrePrints.
Crossref
|
|
|
Engel MS (1999). The taxonomy of recent and fossil honey bees (Hymenoptera: Apidae: Apis). Journal of Hymenoptera Research 8:165-196.
|
|
|
Food and Agriculture Organization (FAO) (2018). Food and Agriculture Organization of the United Nations. Why bees matter: the importance of bees and other pollinators for food and agriculture. Available at:
View
|
|
|
Felsentein J (2014). PHYLIP (Phylogeny Inference Package) Version 3.695. Department of Genome Sciences, University of Washington, Seattle.
|
|
|
Franck P, Garnery L, Loiseau A, Oldroyd BP, Hepburn HR, Solignac M, Cornuet JM (2001). Genetic diversity of the honeybee in Africa: microsatellite and mitochondrial Data. Heredity 86:420-430.
Crossref
|
|
|
Garnery L, Cornuet JM, Solignac M (1992). Evolutionary history of the honeybee A. mellifera inferred from mitochondrial DNA analysis. Molecular Ecology 1:145-154.
Crossref
|
|
|
Hammer O, Harper D, Ryan P (2006). PAST- Paleontological Statistics Software Package for Education and Data Analysis version 1.58. Paleontologica Electronica 4(1):9.
|
|
|
Howpage D (1991). The apiculture development project of Sri Lanka. Journal of Beekeeping Development 19:10-11.
|
|
|
Lawal OA, Banjo AD (2010). Appraising the bee keeping knowledge and perception of pests' problem in beekeeping business at different ecological zones in South-Western Nigeria. World Journal of Zoology 5(2):137-142.
|
|
|
Marghitas LA, Paniti-Teleky O, Dezmirean D, Mărgăoan R, Bojan C, Coroian C, Laslo L, Moise A (2008). Morphometric differences between honey bees (Apis mellifera carpatica) populations from transylvanian area. Zootehnie si Biotehnologii 41(2):309-315.
|
|
|
Michener CD (2007). The Bees of the World. 2nd edition The Johns Hopkins University Press, Baltimore.
|
|
|
Nedic N, Jevtić G, Jež G, Anđelković B, Milosavljević S, Kostić M (2011). Forewing differentiation of the honey bees from Serbia. Biotechnology in Animal Husbandry 27(3):1387-1394.
Crossref
|
|
|
Nei M (1978). Estimation of average heterozygosity and genetic distance from a smaller number of individuals. Genetics 89:583-590.
|
|
|
Oyerinde AA, Ande AT (2009). Distribution and impact of honey bee pests on colony development in Kwara State, Nigeria. Journal of Agriculture and Social Sciences 5:85-88.
Crossref
|
|
|
Oyerinde AA, Dike MC, Banwo OO, Bamaiyi LJ, Adamu RS (2012). Morphometric and landmark based variations of Apis mellifera L. wings in the savannah agro-ecological zone of Nigeria. Global Journal of Science Frontier Research 22:33-41.
|
|
|
Palmer MR, Smith D, Kaftanoglu O (2000). Turkish honeybees: genetic variation and evidence for a fourth lineage of A. mellifera mtDNA. Journal of Heredity 91:42-46.
Crossref
|
|
|
Peakall R, Smouse PE (2006). GENALEX 6: Genetic analysis in Excel. Population genetic software for teaching and research. Molecular Ecology Notes 6:288-295.
Crossref
|
|
|
Peakall R, Smouse PE (2012). GENALEX 6.5: Genetic analysis in Excel. Population genetic software for teaching and research- an update. Bioinformatics 28:2537-2539.
Crossref
|
|
|
Ruttner F, Tassencouyt L, Louveaux J (1978). Biometrical statistical analysis of the geographic variability of Apis mellifera L. Apidologie 9:363-381.
Crossref
|
|
|
Ruttner F (1988). Biogeography and Taxonomy of Honeybees. Berlin, Heidelberg, New York: Springer-Verlag.
Crossref
|
|
|
Smith DR (1991). Mitochondrial DNA and honeybee biogeography. Diversity in genus Apis, Westview Press, Oxford. pp. 131-176.
Crossref
|
|
|
Sharma PC (1983). Morphometric studies on Apis florea F. and Apis dorsata F. of Himachal Pradesh and Punjab. Thesis. Himachal Pradesh University, India.
|
|
|
Tahmasebi G, Ebadi R, Tajabadi N, Akhoundi M, Faraji S (2002). The effects of geographical and climatic conditions on the morphological variation and separation of Iranian small honeybee (Apis florea F.) populations. Journal of Science and Technology of Agriculture and Natural Resources 6(2):169-176.
|
|
|
Tunca RI, Kence M (2011). Genetic diversity of honey bee (Apis mellifera L.: Hymenoptera: Apidae) populations in Turkey revealed by RAPD markers. African Journal of Agricultural Research 6(29):6217-6225.
Crossref
|
|
|
Whitfield CW, Behura SK, Berlocher SH, Clark EG, Johnston JS, Sheppard WS, Smith DR, Suarez AV, Weaver D, Tsutsui ND (2006). Thrice out of Africa: ancient and recent expansions of the honey bee, Apis mellifera. Science 314(5799):642-645.
Crossref
|
|
|
Yu L, Xie W, Wu H, Zou Y, Nan Q, Zhu, L, Ji H, Wu Q (2012). Morphological characteristics and microsatellite DNA genetic diversity of Nigeria African honey bee, AnhuiApis mellifera and their hybrid generation II. Acta Ecologica Sinica 32:3555-3564.
Crossref
|
|