Genetic Relatedness of Low Solitary Nests of Apis dorsata from Marang, Terengganu, Malaysia

Knowledge on the population of genetic structure and ecological behaviour of Apis dorsata from Peninsular Malaysia is needed for effective management and conservation of this species since unsustainable whole solitary low nest cutting for product harvesting is the current common practice here. The analysis of 15 single locus DNA microsatellite markers on samples from 20 solitary nests of A. dorsata showed that while these markers were polymorphic, high intracolonial relatedness existed. Furthermore, in general, slightly negative values of intercolony relatedness (R) among the nests of A. dorsata were found. However, positive values of mean intercolony relatedness were observed between 54 pairs of nests out of 190 possible combinations. The R values among nest pairs 3–4 and 3–5 was higher than 0.50 showing that their queens were half siblings, whereas nest pair 19–20 showed relatedness of 0.95 indicating that the same queen was sampled. The results that we obtained could not conclusively support the hypothesis of this study that the honey hunters in Marang district of Malaysia repeatedly harvest the same nest located at a different site and at a different time during the same honey harvesting season. However, our finding of an appreciable level of intercolonial relatedness between several pairs of nests in this pioneer study indicated that a comprehensive study with a larger sample size of solitary nests found throughout the region would be necessary to provide concrete proof for this novel idea.


Introduction
The Indo-Malaya region is known as the centre of honeybee diversity, whereby seven out of the nine honeybee species of the world are sympatric and endemic to this region [1]. Among the Apis species, the nests of A. dorsata are known as the sources of wild honey in Malaysia [2]. Many colonies of the single-comb, opennesting and sedentary of A. dorsata are found to nest either singly or low to the ground, or high in aggregates on tree limbs of tall bee trees along the coastal, submerged Melaleuca forest in Marang, Terengganu, Malaysia [3]. Of the two nest types, it is the low solitary nests that are usually harvested in Peninsular Malaysia and thus is the subject of this study.
Seasonal migration and aggregation found densely on a bee tree are known as unique characteristics in A. dorsata, which may differentiate this particular species from other Apis species [4]. The seasonal migration between alternative nesting sites is done to find available forage [5] and control the levels of Tropilaelaps clarea (parasitic mite) [6]. While the aggregation structure (up to 50 or 100 nests of A. dorsata on a bee tree high in the air [7], usually from 5 to 40 metres above the ground in the rainforest) is a unique feature of this honeybee species, abundance of low aggregations (2 to 10 nests on a tree) and solitary nests (only one or rarely two nests on a tree) are more common in this species [7], especially in the secondary forests in the district of Marang, Terengganu, Malaysia ( Figure 1). Both low aggregations and low solitary nests are found on trees of less than 5 metres in height, which makes it easy for honey hunters to harvest these nests. The honey hunters in Marang district climb the tree and remove the whole comb of this type of A. dorsata nests easily during the days of harvesting season (Personal observation).
In Malaysia, the interaction between man, the forest and A. dorsata has been established for several hundreds of years [8]. Most honey hunters are found in the states of Kedah, Terengganu and Negeri Sembilan in the Peninsular Malaysia [8]. Collection of wild honey is lucrative and generates an income of about RM6,000 (USD2000) per harvesting season in about three months [8].
Because of the inherent organic properties, high medical and nutritional values of wild honey, its price is much higher than the commercial honey, which may be produced by A. mellifera and A. cerana [9]. Due to the large size of the A. dorsata nest, a considerable amount of wild honey (up to 45 kg) [10] may be stored by a nest. This amount of honey tempts local honey hunters to harvest the nests of A. dorsata during flowering season frequently. Therefore, Malaysian honey hunters harvest the whole nest of A. dorsata by cutting and taking the different parts of a nest for selling. This method of nest harvesting is especially done on solitary nests, which are spread throughout the rainforest of Malaysia. This common harvesting method is deemed unsustainable as the bees may avoid this area as their nesting site in the future causing a decrease of the A. dorsata population here.
Honey hunters in Marang claim that they harvest about 600 or more nests of A. dorsata per harvesting season (June, July and August) in this area (Personal communication). This figure is seemingly huge, and based on the A. dorsata biology and plant flora, it is not possible to have this huge number of A. dorsata nests in this area. Furthermore, A. dorsata needs specific plant species for providing nectar, pollen and establishment of its nest. Saberioon et al. (2010) reported that in Marang, A. dorsata constructs its nest on the branches of available nesting support of Acaia auriculiformis and Melaleuca cajuputi trees which are the sources of nectar and pollen for producing strong flavour and weak density of wild honey in the forested areas of Marang district [11].
We hypothesize that the honey hunters in Marang harvest the same nests, in which the queen of each nest heads away from the original nest site together with her nest members after the first harvest to construct a new nest at another site in the same season and in the available nesting support in the vicinity of the same geographical location. Information on the genetic structure and relatedness among the solitary nests of A. dorsata could provide us with the estimation of honey harvesting dates from the new nests to help rectify the current unsustainable nest harvesting practices.  Paar et al. (2004) studied the genetic structure and relatedness among and within the A. dorsata aggregations in northeast India [5]. Their results revealed significant genetic differentiation among the aggregations. They also found positive relatedness among nests within aggregations which indicated that the colonies within aggregations were more highly related than random colonies. However, they did not study the low solitary nests it was necessary for us to do so since this nest type is mainly harvested in Malaysia.
There is currently no data on the genetic relatedness and population structure of A. dorsata solitary nests from Marang district (Terengganu, Malaysia). The knowledge on the genetic structure of A. dorsata is needed for effective management and conservation strategies of this species in Malaysia, especially to prevent over harvesting of honey. Overall, the objectives of this study were to estimate the level of intercolonial and intracolonial relatedness and to identify the population structures of Malaysian A. dorsata low solitary nests using single locus DNA microsatellites. According to Insuan et al. (2007), microsatellites are efficient markers to reveal the detailed genetic structure of A. dorsata population [12].

Microsatellite Amplification and Genetic Variation
All the 15 loci were successfully amplified with the expected sizes of PCR products ( Table 1).The inferred genotypes at 15 loci of the queen of each solitary nest are presented in Table 2. The number of the effective alleles displayed per locus (allelic diversity) acts as a useful measure of genetic variability within a population [13]. A total of 134 alleles were found at the 15 microsatellite loci in the honeybees of Marang district population. The number of alleles per locus ranged from 6 alleles for locus Ap273 to 11 alleles for locus Ag005a. All loci were found to be polymorphic. H o ranged from 0.30 to 0.60 across the 15 loci, with mean of 0.4160.025 and mean of H e 0.8260.007 (Table 3).
Higher PIC values of a molecular marker indicate higher heterozygote frequency in a population, as well as more genetic information [14]. In this research, all loci showed high polymorphism. It was found that all PIC values were less than their related heterozygosity. The obtained results indicated that the loci with more alleles contained higher rates of heterozygosity and PIC values, but this was not absolute because of the effect of allele frequencies. For example, Ag005a had 11 alleles, while its related heteozygosity and PIC values were lower (H e = 0.829, PIC = 0.646) than locus Ap049 with 10 alleles (H e = 0.854, PIC = 0.718). Thus, it was concluded that at Ap049 locus, the rare alleles contributed little effect on heterozygosity. Furthermore, the heteozygosity at this locus exceeded the expected heterozygosity computed from the number of alleles as a result of mutation drift equilibrium [15]. The mean effective number of alleles (5.73860.223) across the 15 microsatellite loci indicated that the sample size of this study was enough [16].
A polymorphic locus must have at least 0.10 of heterozygosity to reflect the variation of genetic structure [17]. The mean heterozygosity (0.82260.007) across the 15 microsatellite loci was sufficient to show genetic variation in the studied population.

Genetic Structure
The deviations from HWE were assessed using the method of Guo and Thompson (1992) for each locus-nest combination using a Markov chain of 100 000 steps and 1000 dememorization steps [18]. For social insects, the colony structure may show significant deviations from HWE and linkage disequilibrium tests whereas no deviation exists in a population if the workers' genotypes are used [5]. So, in this research, these analyses were carried out using the queen's alleles only. The results showed that no nest had significant deviation from HWE (P.0.05). There was also no significant deviation from genotypic linkage disequilibrium (P.0. 34). The results also demonstrated that the low solitary nests which were sampled for this study were mating randomly.

Population Genetic Differentiation
Among closely related populations, genetic drift is assumed to be the main factor in genetic differentiation [19]. Hence, F ST index is suggested to be used for estimating the differences between populations [20,21]. Theoretically, F ST vales change between 0 and 1 [22]. All F ST values between the solitary nests of A. dorsata were significantly different from zero (P,0.01). A moderate level of genetic differentiation (F ST ) among the sampled solitary nests (mean of F ST = 11.6%) was found. In addition, the AMOVA indicated that this percentage of differentiation was due to the differences among the solitary nests (11%), an appreciable amount of genetic variation among workers of all nests (59%), and within workers of each nest (30%), as shown in Figure 2. The low percentage of differentiation (11%) among the solitary nests could  support this research hypothesis, which states that the honey hunter may be harvesting the same nest at a different site of a same geographical area and at a different time. The genetic variation among A. dorsata workers of all nests and within workers of each nest were resulted as the reproductive system in honeybees and high levels of polyandry [23].

Genetic Relatedness among the Solitary Nests of A. dorsata
We followed Oldroyd et al.'s (1995) analytical approach to determine whether the spatial distributions of the solitary nests of A. dorsata were random or otherwise [16]. The Poisson method and negative binomial method were used. An organism under Poisson distribution is dispersed randomly, whereas an aggregated organism follows negative binomial distribution in the environment (Dale, 2000 cited in [16]). Both types of distribution were tested, because two types of nests (the solitary and aggregations) were found in the study site. This paper is focused on the results of our work on the more frequently harvested solitary nests. The observed and expected distributions of nests under these models showed that the observed distribution differed significantly (P,0.01) from the negative binomial distribution, indicating that the solitary nests were distributed randomly throughout the plots (Figure 3).
The degree of intercolonial genetic relatedness (R) among the solitary nests of A. dorsata revealed positive relatedness between some of the nest pairs based on the queen's heading of them. The mean relatedness coefficient among the solitary nests of A. dorsata was negative (R = 20.05360.016). The 95% confidence interval was determined by bootstrap resamplings and limits of 20.092 and 0.094 were found. Under Kin selection hypothesis, the majority of relatedness coefficients among pairs of queens should be positive. However, positive values of mean intercolonial relatedness were observed only between 54 pairs of nests out of 190 possible combinations ( Table 4). The R values among nest pairs 3-4 and 3-5 was higher than 0.50 showing that their queens were half sibling, whereas nest pair 19-20 showed a relatedness of 0.95 indicating that the same queen was sampled ( Table 2 and  Table 4).
Despite the tendency of A. dorsata to aggregate on a single tree which suggested that the colonies within the aggregations might show higher degree of relatedness than random colonies, Paar et al. (2004) found that the long distance migration between the original and alternative sites might minimize the genetic differentiation among the geographical locations as well as within the aggregations of this species [5]. Hence, the low genetic differentiation among and within the A. dorsata aggregations was due to gene exchange among them. The slightly negative values of intercolony relatedness within the aggregations of A. dorsata were found using microsatellites [5]. In their study, no pairs of colonies within the aggregations showed high relatedness. Only one pair of the queens carried at least one identical allele at every locus. These queens were considered to have formed mother and daughter colonies [5]. Paar et al. (2004) further explained that the limited number of bee workers per colony and the loci which they used affected their results.
Based on the objectives of the present study, the pair wise nests numbered 19 and 20 showed the highest intercolonial relatedness (R = 0.95) among all the pairs of nests (Table 4) while queen pair carried out two same identical alleles at 10 loci and at least one allele for the remaining five loci. On the contrary, the nest number 16 had the lowest degree of relatedness (R = 0.016) with nest numbers 19 and 20. Despite the distance between plot B and plot C, and the different date of sampling at these sites, 95% of genetic relatedness was found between one pair of nests (numbered 19 and 20). This high level of relatedness between this pair of nests raised the possibility that the honey hunters might have harvested the same nest after the nest shifted to a new site in the same season and at the same geographical area. Hence, it was concluded that the plundering or destruction of a particular A. dorsata nest for honey harvesting caused the queen to head the nest members away to establish a new nest at different site but within the same geographical area.
Along the degree of intercolonial relatedness, intracolonial pedigree relatedness (r) was calculated based on the effective number of drones which mated with the queen of each nest (Crozier, 1970 cited in [5]). In the present study, the effective number of drones was found for each nest and locus ranged between 7.8962.1 and 14.2963.07, followed by a mean average of intracolonial relatedness of 0.28160.019. Oldroyd et al. (1996) reported that the A. dorsata queen mates with an effective number of 20.066.6 and a mean average intracolonial relatedness of 0.2960.007 [23]. They determined the (r) levels using three microsatellites and a large number of workers per nest . Despite the number of microsatellites and workers per nest used in our study and that of Oldroyd et al. (1996), a similar mean average of intracolonial relatedness was found. Our results approach was in line with the view of Takazaki and Nei (1996) [20] that using more microsatellite loci is better for obtaining genetic relationships when the closely related samples are being studied and the average heterozygosity is high.
Our finding on the nesting behaviour of the A. dorsata bees sampled along the coastal, submerged Melaleuca forest in Marang (Malaysia, Terengganu) should be taken into account by the Forestry Department of Peninsular Malaysia in their effort to formulated possible sustainable methods of nest product harvesting rather than the present practice of cutting off the entire nest.

Study Area
This study was conducted in Marang district, which is located in the state of Terengganu at the northern east of Peninsular Malaysia, between the upper left of 5u019 N, 103u119 E and the lower right of 4u509 N, 103u249 E. Mangrove, Melaleuca, Acacia, rubber and coconut trees are the dominant trees species in this area [11]. The locations of these trees are not in any way restricted information. Hence, no specific permits were required for this study and the authors also confirm that the nests sampled did not involve any endangered or protected species.
In this research, we followed Saberioon et al.'s (2010) prediction map for finding the solitary nests of A. dorsata in Marang district [11]. The samplings were done during the harvesting season of 2010. The map was prepared for the A. dorsata nesting sites [11], and the positions of the solitary nests were added to that map for the present study ( Figure 3). The nests numbered 1 to 9, 10 to 19, and 20 to 23 were sampled during the last week of June, third and fourth weeks of July, and third week of August of 2010, respectively. The distribution data of the solitary nests is shown in and Table 5. The nests numbered 1 to 20 were solitary nests, while those numbered 21 to 23 were the aggregated nests on three different bee trees (i.e. up to 20 nests on a bee tree). Therefore, in this study, only the low solitary nests were analyzed to achieve the research objectives and kept the aggregated nests for the next research. Furthermore, one-kilometre distance was estimated as the forage area of each nest [11]. Then, at least onekilometre distance was kept between the pair wise nests of each plot (Figure 3).

Sampling and DNA Extraction
Based on the effective number of mating in A. dorsata (20.066.6) [23], 20 bee workers per nest were collected. The entire nests were removed and the bee workers were taken from the honey part of each nest separately. The bee workers were preserved at 95% of ethanol and stored at 280uC for DNA extraction later.
The total genomic DNA was isolated from the thorax of each bee worker of each nest using the WizardH Genomic DNA Purification Kit (Promega, USA). The protocol of DNA extraction followed the animal tissue protocol, which is available at http:// www. promega.com/resources/protocols/. The concentrations of the extracted DNA stocks were measured using a spectrophotometer (Eppendorf Biophotometer Plus). The concentrations of DNA ranged from 200-370 ng. Then, the needed volume of DNA stocks was diluted to obtain 10 ng/ul concentration of working DNA. After testing for the purity and dilution, all the DNAs were labelled. So that each DNA sample can be traced to the worker and the nest from this was originated. The DNAs were kept at 220uC for PCR amplification later.

Analysis of Microsatellite Polymorphisms
Twenty two (22) heterospecific primer sets of DNA microsatellites (including Ad3, A007, AB024, A088, AB124, Ap036, Ap043, Ap049, Ap068, Ap085, Ap207, Ap226, Ap243, Ap273, Ap297, At003, Ag005a, A107, A76, B124 and A14) were screened and 15 polymorphic loci (Table 1) were selected for achieving the objectives of the research using the reaction conditions described in Paar et al. (2004). In the present research, however, the nonlabelled primers were used. The loci were size fractionated on 8% of denaturating polyacrylamide gels, visualized by a modified silver staining protocol [24] and then their sizes were estimated by a  comparison to a 20-bp ladder marker (Fermentas). The fragments were scored using the UVIDoc version 99.02 software. The DNA fragment sizes as the base pairs were entered into the Excel sheets as co-dominant markers and were saved for data analysis.

Statistical Analysis
In social insects that follow the haplodiploidy sex determination system, the genetic structure may be determined using the genotypes of the queen and paternal genotypes. In this study, the workers' genotypes per nest of A. dorsata were converted to COLONY [25] software format. This software determines the inferred genotypes of the queen and paternal genotypes of each nest. The inferred genotypes were subjected to genetic analysis. The inferred queen's genotype of each nest was conducted to calculate the number of alleles per locus (N a ), effective number of alleles (N e ), observed and expected heterozygosity (H o and H e ), polymorphism information content (PIC), and relatedness (R) among the queens which head the low solitary nests separately.
The genetic differentiation measure (F ST ) was estimated using the paternal inferred genotypes. According to Paar et al. (2004), the paternal genotypes were diploidized as co-dominant marker. GENALEX version 6.41 [26] was used to estimate the following measures of genetic variation N a , N e , H o , H e , and R, and genetic differentiation structure (F ST ). Hierarchical genetic structuring based on the analysis of molecular variance (AMOVA) [27] and bootstrap resampling were executed by using GENALEX version 6.41 [26]. Meanwhile, PIC values were calculated using the microsatellite toolkit package [28]. Hardy Weinberg equilibrium (HWE) tests and genotypic linkage disequilibrium were performed using GENEPOP version 3.3 [29], whereas the intercolonial relatedness between the solitary nests of A. dorsata was computed based on the queen's inferred genotype using an algorithm as described by Queller and Goodnight (1989) [30].

Conclusion
This study was conducted to show the genetic relatedness among the solitary nests of A. dorsata in the Marang district, Malaysia. The results revealed the high intracolonial relatedness within the studied nests based on the data obtained from the use of 15 single locus DNA microsatellite markers. The presence of appreciable levels of intercolonial relatedness between several single pair of nests indicated the necessity for a comprehensive study involving a large sample size to confirm whether the honey hunters in this area of Malaysia harvest the same nest at a different sites and dates during a single harvesting season or not. Such research is needed to formulate effective management strategies to replace the current common unsustainable A. dorsata nest harvesting practice of whole solitary nest removal. Genetic studies on the aggregate nests should also be done in the future.