Assessment of Genetic Variation in Apis nigrocincta (Hymenoptera: Apidae) in Sulawesi Revealed by Partial Mitochondrial Cytochrome Oxidase I Gene Sequences

Asian cavity-nesting honey bee Apis nigrocincta, a native bee species from Sulawesi and the Philippines, plays a vital role in pollinating flowering plants in local ecosystem and agriculture. In this study, we assessed the intraspecific genetic variation of A. nigrocincta using the sequence of cytochrome c oxidase subunit I (COI). Molecular phylogenetic analysis showed that there were three main clades in A. nigrocincta specimens from Sulawesi based on their respective locations (North, Central, and South Sulawesi). Genetic distance analysis using the Kimura 2-parameter (K2P) model showed that the intraspecific genetic distance in Sulawesi specimens ranged from 0.000 to 0.055. There are 26 nucleotide polymorphic sites within Sulawesi A. nigrocincta. The variation was dominated by transition T ↔ C. The molecular identification result was supported by morphological identification. The results of the two methods agree that the specimen under study was A. nigrocincta. The result of genetic distance calculation indicated that although the tested specimens were derived from remote locations, the genetic variation was still within the range of intraspecific variation.


Introduction
ere are about 30,000 bee species distributed worldwide. About 17,000 of them have been described [1]. Of them, about 20,000 species belong to the superfamily Apoidea [2]. Apidae, the largest family within this superfamily, contains at least 5,700 species of bees [3]. Honey bees are included in the genus Apis. ese insects play a significant role in the pollination of important crops. Currently, there are nine species of honey bees known to inhabit the world: Apis dorsata, A. laboriosa, A. mellifera, A. florea, A. andreniformis, A. cerana, A. koschevnikovi, A. nigrocincta (Sulawesi, Indonesia, and Mindanao, the Philippines), and A. nuluensis (Kalimantan, Indonesia), with approximately 44 subspecies [4]. e honey bees are prone to extinction due to their haplodiploid sex-determining mechanism, limited population because of anthropogenic activities [5], and low intraspecific genetic variation. e intraspecific genetic variation in A. nigrocincta has been poorly studied. In recent years, molecular methods have been applied to study intraspecific genetic variation, species identification, and phylogenetic relationship among close-related taxa [6]. One technique used for this purpose is DNA barcoding. is method has been applied to conservation biology field because it provides information relevant to wildlife conservation management [7]. For animals, the area used for DNA barcoding is the cytochrome c oxidase I (COI) gene. In this study, we examined the intraspecific variation in A. nigrocincta to facilitate species identification and analyze their genetic distances.

Sample Collection.
e adult worker bees were captured using a swingnet from an abandoned local garden in Mahakeret village, city of Manado, North Sulawesi province, Indonesia (1°28'53"N and 124°50'21"E). e location was situated at an altitude of about 40 meters above sea level. e average temperature at sample collection place (March 2019) was between 29°C and 32°C with humidity of 85-89%. e percentage of vegetation cover was 80%.

Morphometric Analysis.
e observed morphological characters included length of following parts: head (HD), antenna (AT), proboscis (PB), thorax (TO), abdomen (AB), fore-wing (FW), hind-wing (FL), midleg (ML), and hind-leg (HL) using a caliper [8].  [6]. e condition of PCR included 2 min of initial denaturation at 95°C followed by 35 cycles of denaturation at 98°C for 10 sec, annealing at 54°C for 30 sec, elongation at 68°C for 45 sec and additional extension for 5 min at 68°C. e PCR products were sequenced bidirectionally using the same PCR primer pairs at 1st BASE DNA Sequencing Services Malaysia.

Analysis of COI Data.
Chromatograms were subjected to the procedure as was done previously [9]. e clean COI sequence of the specimen was deposited in GenBank (http://www.ncbi.nlm.nih.gov). Identification was performed using BLAST identity search provided by the same platform. e clean sequences were aligned using Clustal O (1.2.1) multiple sequence alignment (http://www.ebi.ac.uk/ Tools/msa/clustalo) with other allied honey bee COI sequences from different parts in Sulawesi retrieved from GenBank. e evolutionary history was inferred by the Maximum Likelihood (ML) method based on the parameter (K2P) model [10]. Evolutionary analyses were conducted in MEGA v10.0.4 [11].

Morphometric Analysis.
e photograph of honey bee specimen obtained from Manado Mahakeret is shown in Figure 1. e specimen hereinafter referred to as CAL01. e results of the nine major morphological characteristics of the honey bee depicted by acronyms as length of head (HD), antenna (AT), proboscis (PB), thorax (TO), abdomen (AB), fore-wing (FW), hind-wing (FL), midleg (ML), and hind-leg (HL) are presented in Table 1.
ese characters were compared with the characters of A. nigrocincta provided by Hadisoesilo et al. [12]. e morphometric analysis has been used previously to study genetic variability in honey bees [13].

Molecular Analysis.
e cytochrome oxidase I (COI) sequence of A. nigrocincta CAL01 has been deposited in GenBank with accession number MK880239. e complete BLAST search is presented in Table 2. e location of all Sulawesi specimens is shown in Figure 2. Because many sequences have 72% query cover, the sequences of some specimens were cut to match the length of other specimens. Even so, percent identity ranged between 92.24% and 100%. With a 72% query cover, the specimen CAL01 had a 100% identity with the DQ020233 specimen from Bogani Nani Wartabone National Park, Gorontalo, Indonesia. With a 100% cover query, the specimen CAL01 had a 92.24% identity with the KY834222 specimen from Islamabad, Pakistan. e molecular phylogenetic analysis by the ML method based on the K2P model is shown in Figure 3. e tree reveals that there were three main clades in Sulawesi A. nigrocincta specimens. Estimation of genetic distance amongst A. nigrocincta is shown in Table 3. Analysis was performed using the K2P model [10] integrated in MEGA v10.0.4 [11]. e intraspecific genetic distance ranged from 0.000 to 0.055 (excluding specimen from Pakistan). Table 4 shows the polymorphic nucleotides of the COI. Twenty-six polymorphic sites were detected within A. nigrocincta. Variation of intraspecific COI gene showed that transition T ↔ C dominated the polymorphic pattern.
Estimation of substitution using the maximum composite likelihood method can be seen in Table 5. Different transition substitution rates are shown in bold and transversion substitution shown in italics. e sum of the values of v was equal to 100. e number of transition substitution was 86.89, while the number of substitution was 13.11. erefore, the transition/transversion ratio (ti/tv) was 6.63. However, sampling size strongly influences the ti/tv. e estimated maximum likelihood ti/tv bias (R) was 9.05 under the K2P model [10].

Discussion
Referring to the morphometric analysis conducted by Hadisoesilo et al. [12], as well as other physical characteristics, it is believed that the honey bee being studied was Apis nigrocincta (Figure 1). Nevertheless, this method still has   short comings due to lack of expertise from the researchers. us, a simpler but more accurate method was needed. One method used to identify insects that has been widely used is DNA barcoding using the COI gene. is COI gene has been widely used in evaluating inter-and intraspecific diversity in insects [14]. It is also used to complement traditional morphological-based identification to get more accurate species identification results [6]. Several studies of DNA barcoding in honey bees using the COI gene have been carried out by several researchers [5,[15][16][17][18]. erefore, molecular-based identification was still carried out in this research. After being confirmed using the COI gene, it was ascertained that the specimen being studied was A. nigrocincta. Figure 3, the specimens were grouped based on location. e first group was from North Sulawesi, the second was from Central Sulawesi, and the third was from South Sulawesi. However, the specimen DQ020219 from Kebun Kopi Central Sulawesi was clustered together with the North Sulawesi group. is can be because the location of Kebun Kopi is very close to Bogani Nani Wartabone National Park which is situated in North Sulawesi. Low nucleotide variations in the data presented in this research in Sulawesi specimens indicated that there was no geographical isolation, hence the gene flew among unrestricted populations, especially in nearby locations. e specimen KY834222 from Islamabad Pakistan was observed to be out of the Sulawesi group. In the phylogenetic tree created by involving its sister species, A. cerana, it was seen that specimen KY834222 was grouped with A. cerana japonica (data not shown). is phenomenon was supported by studies conducted by Damus and Otis [19] which stated that several subspecies of A. cerana and A. cerana japonica were confirmed as being distinct from the rest of A. cerana.

As shown in
is indicates that Apis nigrocincta found in Pakistan was most likely the part of the group A. cerana japonica. Another finding showed that A. nigrocincta from Sulawesi mainland and Sangihe Island were embedded in the A. cerana group [20]. e phylogenetic analysis was consistent with morphological and molecular evidence, indicating that A. nigrocincta had similarities with A. cerana [21,22]. e closest genetic distance was 0.000, between haplotype CAL01 and DQ020233 from Bogani Nani Wartabone National Park, Gorontalo (Indonesia). e farthest genetic distance was 0.082, between haplotype CAL01 and KY834222 from Islamabad (Pakistan), and 0.052 with DQ020223 from Bantimurung, South Sulawesi (Indonesia). e genetic distance was 0.055 between NC038114 from Sangihe Island North Sulawesi (Indonesia) and DQ020223. With this value, Takahasi et al. [23] stated that A. nigrocincta maintained a high specific genetic diversity on Sulawesi Island. e result of genetic distance calculation indicated that although the tested specimens were derived from remote locations, it was still within the range of intraspecific variation, except with specimen from Pakistan. Maximum pairwise distance of A. mellifera in 10 locations across China and Pakistan was 0.039 [18]. e result amongst Pakistan's haplotypes was 0.027 and China's haplotypes was 0.010. e maximum pairwise distance 0.053 among A. cerana was detected between Flores' and Taiwan's haplotypes [24]. e closer the genetic distance, the closer the kinship among the organisms being compared [25]. Different geographical conditions and remote locations can cause fairly high genetic  Figure 3: Molecular phylogenetic analysis by maximum likelihood method based on the Kimura 2-parameter model [5].    variation in a species and are usually characterized by morphological differences [26,27]. e lack of data collection on various species and the dependence of the morphological identification process have led to debate for species that have similar morphological structures but are classified as different species due to genetic variation [28,29]. e greater the value of genetic distance between populations or individuals, the more isolated they are from one another. Genetic distance indicates the possibility of the influence of geographical isolation on a population [30,31].
No information was provided earlier on the polymorphic sites of A. nigrocincta. However, six polymorphic sites were detected in Bombus koreanus in China [32] and Apis mellifera in United States [33]. e transition T ↔ C was found dominating the variation in A. cerana in India [34]. Due to molecular mechanism, transitions are more often found in a higher frequency than transversions. Furthermore, the transition tends to produce less amino acid substitution. erefore, the transition tends to be more stable, and is a silent substitution as a single nucleotide polymorphism in one population [25]. e ability of COI gene to differentiate between A. nigcrocinta group from different locations in Sulawesi has been demonstrated in this study. is molecular marker was also able to identify A. nigrocincta successfully based on the similarity sequence with other A. nigrocincta specimens.
is finding provides valuable information on the effective use of COI as molecular marker in research on phylogenetic studies.

Conclusions
e current study assessed the intraspecific genetic variation in A. nigrocincta using COI gene sequences. e results showed that there were three main clades of Sulawesi specimen, namely, North, Central, and South Sulawesi. is finding suggested that genetic variation in A. nigrocincta in Sulawesi is still in the category of intraspecific variation.

Data Availability
e data related to this article are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.