Open access journal: http://periodicos.uefs.br/ojs/index.php/sociobiology ISSN: 0361-6525 DOI: 10.13102/sociobiology.v64i1.1179Sociobiology 64(1): 33-41 (March, 2017) Morphometric diversity and phylogenetic relationships among Iranian honey bee (Apis mellifera meda Skorikow, 1829) populations using morphological characters 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). Sociobiology An international journal on social insects A Rahimi1,2, A Mirmoayedi1, D Kahrizi3, L Zarei3, S Jamali1 Article History Edited by Yi Juan-Xu, South China Agricultural University, China Received 07 August 2016 Initial acceptance 23 October 2016 Final acceptance 10 March 2017 Publication date 29 May 2017 Keywords Dendrogram, Honey bee, Morphometric, Morphometric diversity, Phylogenetic relationships. Corresponding author Ataollah Rahimi Department of Plant Protection, College of Agriculture, Campus of Agriculture and Natural Resources, Razi University, Kermanshah, Iran E-Mail: Rahimi.ata.1@gmail.com 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 1- Department of Plant Protection, College of Agriculture, Campus of Agriculture and Natural Resources, Razi University, Kermanshah, Iran 2- Department of Animal Science, College of Agriculture - Kifri, Garmian University, Kalar, As Sulaymaniyah, KRG of Iraq 3-Department of Plant Breeding and Biotechnology, College of Agriculture, Campus of Agriculture and Natural Resources, Razi University, Kermanshah, Iran RESEARCH ARTICLE - BEES A Rahimi et al. – Morphometric diversity and phylogenetic relationships in Iranian honey bee34 parts of Middle Eastern group (O), subspecies from North-East Africa and Eastern part of Ethiopia group (Y) (Ruttner, 1988, 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 and Apis 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 km2. 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. Fig 1. The geographical locations of honey bee populations collected in this study. Sociobiology 64(1): 33-41 (March, 2017) 35 Materials and Methods 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. (1985, 1988). To measure morphometric Table 1. Sampling localities, geographical positions and number of honey bees used for morphometric analysis. Province The number of apiary The number of colonies The number of bees Latitude Longitude Altitude Kermanshah 5 7 140 34° 29′ 15.62 N 047° 05′ 57.65 E 1664 Hamadan 5 8 160 34° 43′ 24.64 N 048° 31′ 01.27 E 2464 Ilam 5 5 100 33° 34′ 31.97 N 046° 27′ 26.44 E 1548 Lorestan 5 6 120 33° 32′ 01.81 N 048° 14′ 12.03 E 1802 Esfahan 5 10 200 32° 36′ 56.09 N 051° 34′ 04.20 E 1618 Chaharmahal and Bakhtiari 5 7 140 32° 28′ 19.40 N 050° 06′ 30.53 E 2582 Fars 5 10 200 28° 59′ 19.33 N 053° 39′ 56.47 E 1563 Sistan and Baluchestan 5 5 100 29° 30′ 53.89 N 060° 50′ 17.17 E 1407 Kerman 5 6 120 29° 18′ 19.25 N 057° 08′ 29.72 E 2838 Kurdistan 5 7 140 34° 42′ 48.06 N 046° 51′ 34.16 E 1775 West Azerbaijan 5 10 200 37° 34′ 04.14 N 044° 44′ 17.19 E 2138 East Azerbaijan 5 10 200 38° 06′ 13.12 N 046° 11′ 10.99 E 1345 Ardabil 5 9 180 38° 13′ 58.28 N 048° 10′ 05.22 E 1455 Zanjan 5 6 120 36° 42′ 37.00 N 048° 22′ 25.99 E 1546 Tehran 5 6 120 35° 43′ 44.11 N 052° 07′ 00.70 E 2233 Razavi Khorasan 5 7 140 36° 25′ 56.06 N 059° 31′ 26.22 E 922 North Khorasan 5 6 120 37° 29′ 16.88 N 057° 09′ 56.89 E 1454 Golestan 5 5 100 36° 50′ 47.46 N 054° 40′ 12.88 E 181 Mazandaran 5 10 200 36° 25′ 22.19 N 052° 14′ 39.71 E 135 Gilan 5 10 200 36° 50′ 23.28 N 049° 25′ 57.83 E 329 characters of the wing, the right side wings were separated and put into a mixture of alcohol and honey, then transferred 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. 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, A Rahimi et al. – Morphometric diversity and phylogenetic relationships in Iranian honey bee36 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) (Table 2). G18D7A4HWWHWLFWWFWLdfS.O.V 188.970 **197.502 **10.185 **0.245 **0.228 **0.043 **1.817 **19Province (A) 30.321 n.s32.846 n.s1.977 n.s0.056 n.s0.061 n.s0.016 n.s0.072 n.s80City (Province) (B) 61.888 **61.936 **2.601 **0.200 **0.124 **0.022 **0.091 **50Error 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). Table 2. Continue TTCHLLTFTLSCSICIPLdfS.O.V 20.436 **1.089 **0.176 **3.887 **0.007 **0.839 **0.728 **19Province (A) 1.989 n.s0.179 n.s0.035 n.s1.266 n.s0.002 *0.256 n.s0.150 n.s80City (Province) (B) 3.616 **0.385 **0.047 **0.880 n.s0.001 n.s0.225 **0.214 **50Error 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 *: significant at p<0.05, **: significant at p<0.01, n.s: non – significant, df: degree of freedom, S.O.V: source of variation PL, TFTL and HLL (mm). Table 2. Analysis of variance of 14 morphometric traits in Apis mellifera meda. 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. Fig 2. Dendrogram following WARD method for clustering populations of honey bees from Iranian provinces. Sociobiology 64(1): 33-41 (March, 2017) 37 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 2nd 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 3td 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 4th factor covered 15.2 % of the data variance and had positive loading. The 5th factor with 10.7 % of the data variation showed negative loading for the width of hind wing and A7 angle. Finally, the 6th factor covered 9.2 % of the variation and had a positive loading for sternite index (Table 5). TTCHLLTFTLSCSICIPIG18D7A4HWWHWLFWWFWL 1FWL 10.225 *FWW 10.228*0.355 **HWL 10.236*0.0040.358 **HWW 10.616**0.626**0.241*0.452 **A4 10.575**0.298*0.532**0.1960.575 **D7 10.636**0.665**0.259*0.579**0.1320.549 **G18 10.261*0.206*-0.1470.1210.1460.1240.335 **PI 10.148-0.114-0.252*-0.144-0.226*-0.1200.144- 0.327 **CI 10.141-0.142-0.106-0.1690.1390.1270.116-0.1840.301**SI 10.115-0.583**-0.389**0.0890.1120.026-0.136-0.1770.238*-0.010SC 1-0.122-0.278*-0.297*0.257*-0.1280.0510.1110.1020.1190.249*0.054TFTL 10.250*0.600**-0.260*0.609**0.467**-0.224*0.544**-0.180-0.138-0.1570.1960.384**HLL 10.651**0.519**0.660**0.245*0.258*-0.427**-0.1840.1240.309**-0.157-0.1560.1810.264 *TCC Table 3. Correlation coefficient between fourteen morphometric traits in honey bee Apis mellifera meda. Fig 3. Two dimensional biplot composed of the 1st and the 2nd Principal Components of morphometric characters of Iranian A.mellifera meda honey bee populations. A Rahimi et al. – Morphometric diversity and phylogenetic relationships in Iranian honey bee38 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 Third Component FWL -0.309 0.158 -0.258 FWW 0.004 0.385 0.437 HWL -0.268 0.249 0.210 HWW -0.270 0.176 -0.360 A4 -0.250 0.354 -0.015 D7 -0.357 0.152 0.132 G18 -0.319 0.239 0.179 PI 0.053 0.375 -0.396 CI 0.324 0.148 0.275 SI -0.260 -0.174 -0.322 SC 0.264 0.301 0.067 TFTL 0.136 0.391 -0.039 HLL 0.337 0.190 -0.248 TTC 0.302 0.237 -0.339 Eigenvalue 5.89 4.16 1.33 Proportion 42.1 29.8 9.6 Cumulative variance (%) 42.1 71.8 81.4 Table 4. Principal component analysis of fourteen morphometric traits in populations of Apis mellifera meda of different parts of Iran. Traits Factor 1 Factor 2 Factor 3 Factor 4 Factor 5 Factor 6 FWL 0.821 -0.142 -0.117 -0.199 -0.268 0.240 FWW 0.249 -0.209 -0.698 0.432 0.213 -0.301 HWL 0.312 0.066 -0.881 -0.078 -0.230 0.107 HWW 0.396 -0.031 -0.205 -0.098 -0.850 0.130 A4 0.474 -0.141 -0.682 0.101 -0.477 0.055 D7 0.869 0.190 -0.359 -0.129 -0.106 0.109 G18 0.857 0.147 -0.410 0.058 -0.185 -0.016 PI 0.174 -0.912 -0.278 0.085 -0.105 0.000 CI -0.439 -0.205 -0.048 0.697 0.331 -0.223 SI 0.223 0.176 0.002 -0.328 -0.107 0.879 SC -0.009 -0.537 0.032 0.640 0.224 -0.360 TFTL -0.006 -0.401 -0.205 0.730 -0.353 -0.238 HLL -0.458 -0.723 0.096 0.432 0.125 -0.063 TTC -0.251 -0.787 0.223 0.363 -0.054 -0.288 Proportion 23.2 18.9 16.2 15.2 10.7 9.2 Cumulative variance (%) 23.2 42.1 58.3 73.5 84.2 93.4 Table 5. Matrix of factors rotation (Varimex rotation) for fourteen traits of Apis mellifera meda in various populations of different locations in Iran. Fig 4. Phylogenetic tree based on UPGMA method to compare subspecies of Apis mellifera meda with different subspecies of A.mellifera. Sociobiology 64(1): 33-41 (March, 2017) 39 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-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, A Rahimi et al. – Morphometric diversity and phylogenetic relationships in Iranian honey bee40 1988, 1992; Garnery et al., 1993; Sheppard et al., 1997; Engel, 1999; Sheppard & Meixner, 2003; Arias & Sheppard, 2005; Meixner et al., 2011). Acknowledgements The authors appreciate the kindness of those beekeepers who allowed us to collect honey bee samples from their apiaries and equally thank the authorities of Plant Protection and Biotechnology departments of Razi university who let us a free access to laboratory facilities which led to fulfill of the present study. References Arias, M. & Sheppard, W. (2005). Phylogenetic relationships of honey bees (Hymenoptera: Apinae: Apini) inferred from nuclear and mitochondrial DNA sequence data. Molecular Phylogenetics and Evolution, 37: 25–35. Dupraw, E.J. (1965). Non-Linnean taxonomy and the systematic of honey bees. Systematic Zoology Journal, 14: 1-2. Engel, M.S. (1999). The taxonomy of recent and fossil honey bees (Hymenoptera: Apidae). Journal of Hymenoptera Research, 8: 165-196. Farshineh Adl, M.B., Vasfi Gencer, H., Firatli, C. & Bahreini, R. (2007). Morphometric chatacterization of Iranian (Apis mellifera), Central Anatolia (Apis mellifera anatolica) and Caucasian (Apis mellifera caucasica) honey bee populations. Journal of Apicultural Research, 46: 225-231. Fetayeh, A., Meixner, M. & Fuchs, S. (1994). Morphometrical investigation in Syrian honey bee. Apidologie, 25: 396 – 401. Garnery, L., Solignac, M., Celebrano, G. & Cornuet, J.M. (1993). A simple test usingrestricted PCR amplified mitochondrial DNA to study the genetic structure of Apismellifera L. Experienta, 49: 1016-1021. Harvey, W.R. (1990). Users guide for lsmlmw and mixmdl, po-2 Version. The Ohio State University, U.S. Kandemir, I., Ozkan, A., & Mohammad, G.M. (2004). A scientific note on allozyme in persian honey bee (Apis mellifera meda) from the Elburz mountains in Iran. Apidologie, 35: 521-522. Kence, M., Jabbari Farhoud, H. & Tunca, R.I. (2009). Morphometric and genetic variability of honey bee (Apis mellifera L.) populations from northern Iran. Journal of Apicultural Research, 48: 247 – 255. Maganus, R. & Tripodi, A. (2011). Mitochondrial DNA diversity of honey bees, Apis mellifera L. (Hymenoptera: Apidae) from queen breeders in the United States. Journal of Apicultural Science, 55: 5-17. Meixner, M.D., Leta, M.A., Koeniger, N. & Fuchs, S. (2011). The honey bees of Ethiopia represent a new subspecies of Apis mellifera - Apis mellifera simensis ssp. Apidologie, 42: 425-437. Meral, K. (2010). Honey bee (Apis mellifera L.) in Turkey: Biodiversity using Geometric Morphometrics Analysis. 4th European Conference of Apidologie, September 2010, Metu- Ankara, Turkey. Page: 134. Rahimi, A. & Asadi, M. (2010). Morphological characteristics of Apis mellifera meda (Hymenoptera: Apidae) in Saghez (west of Iran). Nature Montenegro, 10: 101-107. Rahimi, A. & Mirmoayedi, A. (2013). Evaluation of morphlogical characteristics of honey bee Apis mellifera meda (Hymenoptera: Apidae) in Mazandaran (North of Iran). Technical Journal of Engineering and Applied Sciences, 3: 1280-1284. Rahimi, A., Asadi, M. & Abdolshahi, R. (2014) A. Genetic diversity of honey bee (Apis mellifera meda) populations using microsatellite markers in Jiroft. Journal of Science and Tactic of Honey Bee, 6: 26-33. Rahimi, A., Miromayedi, A., Kahrizi, D., Abdolshahi, R., Kazemi, E. & Yari KH. (2014) B. Microsatellite genetic diversity of Apis mellifera meda skorikov. Molecular Biology Reports, 41: 7755–7761. Rahimi, A., Hasheminasab, H. & Azati, N. (2015). Predicting honey production based on morphological characteristics of honey bee (Apis mellifera L.) using multiple regression model. Ecology-Environment-Conservation, 21: 29-33. Ruttner, F., Tassencourt, L. & Louvaux, J. (1978). Biometrical – statistical analysis of the geographic variability of Apis mellifera L. Apidologie, 9: 363 – 381. Ruttner, F., Pourasghar, D. & Kauhausen, D. (1985). Honey bees of Iran .2. Apis mellifera meda Skorikoow. The Persian bee. Apidologie, 9: 363-381. Ruttner, F. (1988). Biogeography and Taxonomy ofHoneybees. Springer-Verlag. Berlin., Germany, 285 pp. Ruttner, F. (1992). Naturgeschichte der Honigbienen. Ehrenwirth Verlag. Münich. Germany. 357 pp. Sheppard, W.S., Arias, M.C., Greech, A. & Meixner, M.D. (1997). Apis mellifera ruttneri, a new honey bee subspecies from Malta. Apidologie, 28: 287-293. Sheppard, W.S. & Meixner, M.D. (2003). Apis mellifera pomonella, a new honey bee subspecies from Central Asia. Apidologie, 34: 367–375. Shykon, J.A., Poschmid, H. & Stort, A. (1991). Genetic relatedness and eusociality. Parastie mediated selection on the genetic composition of groups. Behavioral Ecology and Sociobiology, 28: 371-376. SPSS. (2004). SPSS for Windows, Release 13.0 Standard Version, SPSS Inc., 1989- 2004. Sociobiology 64(1): 33-41 (March, 2017) 41 Tahmasebi, GH. (1996). Morpological and biochemical survey of honey bee (Apis mellifera L) populations in Iran. Ph.D. Thesis, Tarbit Modares University, Tahran, Iran. Tahmasebi, GH., Abadi, R., Esmaili, M. & Kambozia, J. (1998). Morphological study of honey bee (Aips mellifare L) in Iran. Journal of Agricultural Sciences and Natural Resources, 2: 89-100.