Morphometric diversity and phylogenetic relationships among Iranian honey bee ( Apis mellifera

Honey bee (Apis mellifera L.) as a natural component of the local biota and beneficial insects make an important contribution to honey, wax, propolis, royal jelly, bee poison production and they are most important pollinators, too by helping plants reproduce and increases in agricultural productivity and environmental conservation (Maganus & Tripodi, 2011; Rahimi et al., 2014 B). Honey bees expanded to Northward to the South of Scandinavia, Southward to the Cape of Good Hope, Westward to the Dakar in Senegal, Eastward to Mashad in Iran and Oman in the gulf of Oman with human-assisted migration. In sixteenth century by the human migrations to Oceanic islands and to North and South American continents, Abstract In this study, the morphometric diversity and phylogenetic relationships of Iranian honey bee populations, were investigated using 14 morphometric characteristics. A total of 2250 young adult worker bees from 20 different populations in 20 different provinces of Iran were collected during June to October 2014. The results of nested analysis of variance showed that there were significant differences (P<0.01) between the provinces for all analyzed morphometric traits indicating the existence of a diversity among them. Correlation coefficient analysis showed a high degree of association among the most of the traits. This correlation coefficient should be a putative mean to improve of certain characters in breeding of honey bee. Principal component analysis revealed three principal components explained 81.5% of the total variation. Cluster analysis using WARD method classified honey bee populations into two main groups. The first group includes the honey bees collected from North, Northwest and West portions of Iran. The second group was represented by the honey bees from Eastern North, Central and Southern regions of Iran. The phylogenetic tree based on UPGMA method divided 29 subspecies of honey bee to 5 distinct clusters. The Iranian subspecies honey bee composed of a shared clade with subspecies of Eastern Mediterranean, Near East and Eastern parts of Middle East (O branch). Sociobiology An international journal on social insects


Introduction
Honey bee (Apis mellifera L.) as a natural component of the local biota and beneficial insects make an important contribution to honey, wax, propolis, royal jelly, bee poison production and they are most important pollinators, too by helping plants reproduce and increases in agricultural productivity and environmental conservation (Maganus & Tripodi, 2011;Rahimi et al., 2014 B).Honey bees expanded to Northward to the South of Scandinavia, Southward to the Cape of Good Hope, Westward to the Dakar in Senegal, Eastward to Mashad in Iran and Oman in the gulf of Oman with human-assisted migration.In sixteenth century by the human migrations to Oceanic islands and to North and South American continents,

Abstract
In this study, the morphometric diversity and phylogenetic relationships of Iranian honey bee populations, were investigated using 14 morphometric characteristics.A total of 2250 young adult worker bees from 20 different populations in 20 different provinces of Iran were collected during June to October 2014.The results of nested analysis of variance showed that there were significant differences (P<0.01) between the provinces for all analyzed morphometric traits indicating the existence of a diversity among them.Correlation coefficient analysis showed a high degree of association among the most of the traits.This correlation coefficient should be a putative mean to improve of certain characters in breeding of honey bee.Principal component analysis revealed three principal components explained 81.5% of the total variation.Cluster analysis using WARD method classified honey bee populations into two main groups.The first group includes the honey bees collected from North, Northwest and West portions of Iran.The second group was represented by the honey bees from Eastern North, Central and Southern regions of Iran.The phylogenetic tree based on UPGMA method divided 29 subspecies of honey bee to 5 distinct clusters.The Iranian subspecies honey bee composed of a shared clade with subspecies of Eastern Mediterranean, Near East and Eastern parts of Middle East (O branch).
the honey bees were transported into the new continents.In such a distribution over vast area, many factors such as geographical isolation and ecological adaptations caused the appearance of different populations and emergence of new subspecies (Ruttner et al., 1978;Ruttner, 1988;Tahmasebi et al., 1998;Rahimi & Asadi, 2010).Twenty nine subspecies and many ecotypes of honey bees recorded from different locations in the world.These subspecies were characterized based on a set of characteristic (morphological, behavioral and molecular) and divided into five evolutionary groups; subspecies from West Mediterranean and North-West European group (M), subspecies from South and Central Africa group (A), subspecies from Central Mediterranean and South Eastern European group (C), subspecies from East Mediterranean, Near East and Eastern parts of Middle Eastern group (O), subspecies from North-East Africa and Eastern part of Ethiopia group (Y) (Ruttner, 1988(Ruttner, , 1992;;Garnery et al., 1993;Sheppard et al., 1997;Engel, 1999;Sheppard & Meixner, 2003;Arias & Sheppard, 2005;Meixner et al., 2011).In the Middle East there are eight subspecies, and between them, Apis mellifera meda exist in Iran (Ruttner et al., 1985;Kandemir et al., 2004;Rahimi et al., 2014 A).The Iranian subspecies honey bee (Apis mellifera meda) are spread over Northern part of Iraq, South Eastern part of Turkey (Ruttner, 1988), North part of Syria (Fetayeh et al., 1994) and in all parts of Iran except Eastern desert regions (Tahmasebi et al., 1996).The Iranian subspecies honey bee has similar morphological characters with Italian subspecies honey bee (Apis mellifera ligustica Spinola, 1806).In other cases, however, biological, molecular and behavioral traits seem to be different between them for example: a high rate of reproduction, a very good tolerance for hibernation and high aggressiveness (Ruttner, 1988).
Morphological characters in zoology is routinely used to differentiate different subspecies from each other, still is the common way to characterize honey bee populations in many countries of the world (Ruttner et al., 1978).Ruttner et al. (1978) collected subspecies of honey bees from different countries and showed that the European subspecies (Apis mellifera carnica Pallmann, 1879 and Apis mellifera caucasica Pallmann, 1889) with a bigger body size located the right part of the axis using Principal Component Analysis (PCA), the African subspecies (Apis mellifera jemenitica Ruttner, 1976 andApis mellifera lamarckii Cockerell, 1906) formed the left side of axis and the Iranian subspecies honey bee occupied the central part of the axis.Kence et al. (2009) grouped seven populations of Iranian and Caucasian subspecies of honey bee existing in Turkey based on morphological characters.The Iranian subspecies honey bee formed a distinct group whereas Caucasian and other populations formed two other groups.In the Meral (2010) study of genetic diversity among Turkish honey bee populations, different Turkish populations formed three distinct groups.Tahmasebi et al. (1998) studied morphological characters of different populations of Iranian honey bees and showed that the most of honey bee populations in Iran belonged to Iranian subspecies honey bee which had significant differences from introduced European subspecies.In another study, Iranian subspecies honey bee have some overlapping characters with Turkish honey bee, but the samples collected from different regions of Iran (Urmia, Tabriz and Tehran regions) formed a distinct group (Farshineh-Adl et al., 2007).Rahimi and Mirmoayedi, (2013) reported that the populations of honey bee from Mazendaran province (North of Iran) formed three groups, and they concluded that the climatic factors should be blamed for the population separation.
Iran is a vast country comprising a land area about 1648195 km 2 .From the geographical point of view, Iran has a very diverse climatic condition with diverse patterns of plant communities forming four climatic zones simultaneously; therefore it has a high potential for beekeeping and honey production.Therefore, genetic structure maintenance and preservation of Iranian honey bees is the first step to explore of this potential.The objective of this study was to determine the morphometric variation and phylogenetic relationships among honey bee samples collected from different locations of from Iran using morphometric characters analysis.

Sampling
A total of 2250 young adult worker bees from 150 colonies of Apis mellifera meda were sampled from 100 different localities distributed in 20 Iranian provinces during June to October 2014 (Fig 1, Table 1).Samples were taken from honey bee colonies in most active apiaries from five cities in each province, one apiary in each city and one to three hives per apiary were selected for sampling.Young adult worker bees directly were collected from the brood areas on the combs, using an open mouth jar containing the cotton soaked to the ether and poured into Pumpel solution (30 parts distilled water, 6 parts formaldehyde, 15 parts ethanol 65% and 2 parts acetic acid).

Selected morphometric characters and morphometric measurements
From every McCarthy glass tube, which contained 20 bees, fifteen bees were randomly selected to measure the morphometric traits.We have chosen fourteen morphometric characters of honey bees to distinguish the local subspecies from widely used forty characters.The measured characters were as follows; FWL -forewing length, FWW -forewing width, HWL -hind wing length, HWW -hind wing width, A4, D7 and G18: angle in the forewing, PI -proboscis length, CI -cubital index, SI -sternite index, SC -scutellum color, TFTL -third & forth tergite length, HLL -hind leg length, TTC -third tergite color.All measurements were done based on Ruttner et al. (1985Ruttner et al. ( , 1988)).To measure morphometric

Statistical analysis
A nested analytical design was used.The analysis of variance and means comparisons were done with SAS statistical package V.14 (Harvey, 1990).Multivariate analytical techniques including correlation analysis, principal component analysis, factor analysis and cluster analysis were performed using SPSS V. 13.0 (2004).To draw a phylogenetic tree between Apis mellifera meda and other subspecies, morphological information of other honey bee subspecies obtained from the World Bank database.

Results
In order to study morphometric diversity and investigate phylogenetic relationships, 2250 worker honey bees from different Iranian honey bee populations were used to analyze fourteen characters.The results were checked with keys provided by Ruttner et al. (1985).The results showed that there were five subspecies of honey bee in our collected samples from different regions of Iran as follows; Apis mellifera meda, Apis mellifera ligustica, Apis mellifera carnica, subspecies hybrid midnite and subspecies hybrid Star-Line.

Analysis of variance
The analysis of variance showed a significant difference among all different Iranian provinces for all characters (P<0.01).There were no differences between cities in the same province for all character, except for the sternite index (P<0.05).Besides, analysis of variance revealed that there were significant differences among hives in each cities for the following characters: forewing length, forewing width, hind wing length, hind wing width, A4 angle, D7 angle, G18 angle, proboscis length, cubital index and third & forth tergite length.No significant differences were observed among hives within cities for sternite index and scutellum color (P>0.05) ( FWL: forewing length, FWW: forewing width, HWL: hind wing length, HWW: hind wing width, A4, D7 and G18: angle in the forewing *: significant at p<0.05, **: significant at p<0.01, n.s: non -significant, df: degree of freedom, S.O.V: source of variation FWL, FWW, HWL and HWW(mm), A4, D7 and G18 (degree).

Correlation coefficient
Pearson correlation coefficients between fourteen characters were calculated (Table 3).Length of fore wing showed a positive correlation with width of fore wing, length of hind wing, angle A4, angle D7, angle G18, length of proboscis, sternite index, hind leg length and third tergite color, but the forewing length had a significant negative correlation with cubital index.

Cluster analysis
Dendrogram constructed using WARD method divided the provinces into two groups.Twelve provinces including Kermanshah, Hamadan, Mazandaran, Ardabil, Golestan, Kurdistan, Ilam, East Azerbaijan, West Azerbaijan, Gilan, Zanjan and Tehran clustered were grouped together, and other provinces including Fars, Sistan and Beluchistan, Esfahan, Chaharmahal and Bakhtiari, Razavi Khorasan, Lorestan, Kerman and North Khorasan provinces were placed together in the second group (Fig 2).The results of cluster analysis confirmed by discriminant analysis.

Principal Component Analysis (PCA)
PCA is a variable reduction technique which maximizes the amount of variance accounted for in the observed variables by a smaller group of variables called components.The first three axes of the principal component analysis explained 81.5 % of the total variation (Table 4).The first two components suggested that diversity in honey bee populations were influenced by forewing length, D7 angle, G18 angle, cubital index, hind leg length, third tergite color, forewing width, A4 angle, proboscis length, scutellar color and third & forth tergite length (Fig 3).

Factor analysis
Principal factors were determined as the variables with the highest projection scores on the principal components (Table 5).Results showed that six main factors accounted for 93.4 % of the total variability.The first factor which explained 23.1 % of total variation emphasized the forewing length, D7 angle and G18 angle, and can be named as the flight factor (Table 5).The results of correlation analysis also confirmed the existence of negative correlation between the length of fore wing, the index of cubital and the length of hind leg which are in accordance with the results obtained by the factor analysis.The 2 nd factor covered 18.9 % of the data variation and had negative loadings for length of proboscis, scutellar color, length of hind wing and the third tergite color.The 3 td factor justified 16.2 % of the data variability and had the maximum of negative loading for the characters the width of fore wing, length of hind wing and A4 angle.The 4 th factor covered 15.2 % of the data variance and had positive loading.The 5 th factor with 10.7 % of the data variation showed negative loading for the width of hind wing and A7 angle.Finally, the 6 th factor covered 9.2 % of the variation and had a positive loading for sternite index (Table 5).

Phylogenetic tree
The phylogenetic trees are generally plotted to find the genetic distances between the collected samples of honey bees.Phylogenetic tree based on morphometric character was drawn using UPGMA method to compare phylogenetic relationships Iranian subspecies honey bee with other subspecies.Phylogenetic tree showed 29 honey bee subspecies could be divided into 5 distinct clusters (Fig 4 ).

Discussion
There is one native subspecies of honey bee in Iran, Apis mellifera meda, which during the past thousands years of Iranian civilization, are adapted to the different climate conditions and plant flora varieties in this area.Also, certain resistance against diseases and pests which attack to honey bees in Iran have been acquired by this subspecies (Ruttner, 1978;Tahmasebi et al., 1998).We have compared morphometric

Traits
The First Component The Second Component The  characters of Iranian honey bee samples with Keys provided for study of subspecies of honey bees by Ruttner et al. (1978), five subspecies of honey bee have been identified in different regions of Iran including Apis mellifera meda, Apis mellifera ligustica, Apis mellifera carnica, subspecies hybrid midnite and subspecies hybrid Star-Line.The only native honey bee subspecies in Iran is Apis mellifera meda (Ruttner et al., 1985;Kandemir et al., 2004;Rahimi et al., 2014 A) and other identified subspecies in this study are not native to Iran.Probably, these subspecies had entered Iran as a result of trafficking queen.
The existence of other subspecies in Iran can be a warning to beekeeping community because it eliminates racial purity Iranian honey bee.Moreover the results of other studies (Tahmasebi et al., 1998;Rahimi & Mirmoayedi, 2013;Rahimi et al., 2014 B) support the particular identity of Iranian subspecies honey bee.Conservation of genetic purity in Iranian subspecies honey bee as genetic resources and most expensive natural heritage should be considered in order to preserve them against admixture with other subspecies.Probably, the admixture of Iranian subspecies honey bee with other regional subspecies is occurring spontaneously in border regions; however we should be worried about the illegal trafficking of queens from other honey bee subspecies into Iran.
The results of analysis of variance showed that there was a significant difference between the provinces for all the honey bee characters (P<0.01) which could support the existence of genetic diversity among bees from Iranian different provinces.Moreover, most of the fourteen morphometric characters had significant positive correlations to each other.Shykon et al. (1991) showed significant positive correlations between length of fore wing and characters such as length of proboscis, length of hind leg and third and forth tergite length.Tahmasebi et al. (1996) observed significant positive correlations between length of fore wings with characters such as length of proboscis, length of hind leg and third and forth tergite length.Rahimi et al. (2015) have found positive correlations between the length of fore wings and the honey production of the Iranian honey bee.
Cluster analysis classified honey bee populations into two groups; a group consisted of populations from North, Northwest and West provinces of Iran while the other group comprised of populations from Eastern North, central and the Southern Iranian provinces.In the study of Kence et al. (2009) on five honey bee populations from North, Northwest and West of Iran, all the populations composed single group in agreement with our results.Tahmasebi et al. (1998) have found different three groups by cluster analysis of Iranian honey bee populations; the first group consisted of the populations from North and Northeast, the second included the populations from Northwest and West and the third consisted of the populations from central parts of Iran.The different results in our findings is likely due to the event of dryness which Iran encountered during last twenty years (1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014), this is because of scarcity of raining and the obligation of beekeepers to displace their hives to new sites.These hive migrations influenced the trade of queens among beekeepers and induced genetic introgression of the bee populations and disrupting the previously groups detected by Tahmasebi et al. (1998).Principal Component Analysis (PCA) revealed three principal components accounted for 81.5 % of the total variance.Most of the variation could be explained by the two first components (71.9 %) supporting the existence of two distinct groups.The first group composed the provinces Tehran, Mazandaran, Gilan, East and West Azerbaijan and Kermanshah had the bees with higher cubital index, the longer hind leg and the lighter color of scutellum.The second group comprised populations belonged to the provinces Kurdistan, Ilam, Zanjan, Kerman and Sistan and Beluchistan with the wide forewing, longer hindwing, wide A4 and G18 angles, longer proboscis, longer third & forth tergite length and the biggest abdomen.
As an exploratory tool to reduce many dependent variables to a few principal ones, we used factor analysis for fourteen morphometric characters.The results revealed six hidden factors that explained 93.4 % of variation between studied populations.Phylogenetic tree drawn based on UPGMA method divided 29 honey bee subspecies into 5 distinct clusters (Fig 4).Dupraw (1965) has grouped the honey bee subspecies based on geographical distribution into five groups.Iranian honey bee subspecies grouped with East Mediterranean, Near East and Eastern Middle East subspecies (O branch).Ruttner (1988), Garnery et al. (1993), Engel (1997), Sheppard et al. (1997), Meixner (2003), Arias and Sheppard, (2005) and Meixner et al. (2011) have used cluster analysis to grouping 29 honey bee subspecies based on morphological, behavioral and molecular traits; they detect five evolutionary groups.The Iranian subspecies honey bee formed a group with East Mediterranean, Near and Middle Eastern subspecies (O branch).
In the current study, the phylogenic tree drawn based on measured morphological characters along with other subspecies data, grouped 29 honey bees subspecies into five groups.Iranian honey bee subspecies with subspecies Apis mellifera ligustica, Apis mellifera anatolica, Apis mellifera syriaca, Apis mellifera cyprica, Apis mellifera adami, Apis mellifera armeniaca, Apis mellifera caucasica and Apis mellifera pomonella were placed in the same group.These subspecies except for Apis mellifera ligustica are the same subspecies of evolutionary group (O branch).The results of this study were in accordance with previous studies.Many authors have speculated about a close relationship for the morphological characters between Iranian honey bee subspecies and Italian subspecies (Apis mellifera ligustica).Whereas these subspecies were clustered in the same group, but the Italian honey bee subspecies has biological, behaviorial and molecular characters different from Iranian honey bee subspecies, so honey bee researchers have grouped Apis mellifera ligustica with central Mediterranean and South East European subspecies into one group (Ruttner,

Fig 1 .
Fig 1.The geographical locations of honey bee populations collected in this study.

Fig 2 .
Fig 2. Dendrogram following WARD method for clustering populations of honey bees from Iranian provinces.

Fig 3 .
Fig 3. Two dimensional biplot composed of the 1 st and the 2 nd Principal Components of morphometric characters of Iranian A.mellifera meda honey bee populations.

Fig 4 .
Fig 4. Phylogenetic tree based on UPGMA method to compare subspecies of Apis mellifera meda with different subspecies of A.mellifera.

Table 1 .
Sampling localities, geographical positions and number of honey bees used for morphometric analysis.
to two layers glass slides and arranged in a numbered regular manner.The wings were attached firmly to glass slides following the evaporation of alcohol.Slide projector was used to measure the wings dimensions in each sample analyzed.Stereomicroscope with calibrated eyepiece was used to measure other body characters.Color characters were measured according Ruttner's international index.

Table 2 .
Analysis of variance of 14 morphometric traits in Apis mellifera meda.

Table 3 .
Correlation coefficient between fourteen morphometric traits in honey bee Apis mellifera meda.

Table 4 .
Principal component analysis of fourteen morphometric traits in populations of Apis mellifera meda of different parts of Iran.

Table 5 .
Matrix of factors rotation (Varimex rotation) for fourteen traits of Apis mellifera meda in various populations of different locations in Iran.