The complete mitochondrial genome of Apis cerana-southern China (Hymenoptera: Apidae) and insights into the phylogenetics

The mitochondrial genome provides crucial information for understanding the evolution and phylogeny of various Apis cerana populations. Apis cerana-Southern China is a unique ecological type of Asian bee, A. cerana, primarily found in the coastal mountains of South China. Here, we used PacBio-HiFi sequencing to obtain the complete mitochondrial genome of A. cerana-Southern China and infer the phylogenetic relationships between A. cerana-Southern China and other A. cerana ecotypes. The mitochondrial genome of A. cerana-Southern China contains 16,137 bp and includes 13 protein-coding genes (PCGs), 22 tRNA genes, 2 rRNA genes, and an AT-rich region. The overall base composition is as follows: A (42.41%), C (9.59%), G (6.18%), and T (41.82%). The combined percentage of A and T (84.23%) is significantly higher than that of G and C. Among the 37 genes, 23 were located on the majority strand, while the remaining 14 were located on the minority strand. The phylogenetic tree based on the 13 PCGs showed that the genetic distance between A. cerana-Southern China, A. cerana-Central China, and A. cerana-Aba was closer. The complete mitochondrial genome sequence reported here would be useful for further phylogenetic analysis and conservation genetic studies in A. cerana-Southern China.

1 Introduction

The honey bee Apis cerana is an important pollinator as it provides not only pollination services but also ingredients for healthy human foods (Chen et al., 2017). However, the populations of A. cerana have declined significantly in recent decades due to climate change, human activities, and bee diseases (Cao et al., 2020; Li et al., 2022; Matias et al., 2017). Mapping its genetic variation is crucial for understanding population-level health, histories, and the potential capacity to respond to environmental changes. Thus, the study of genetic diversity is fundamental to the conservation and utilization of honey bee resources, and it is particularly crucial for A. cerana populations that are experiencing rapid declines. Therefore, it is important to elucidate the mechanisms of A. cerana environmental adaptation with sensitive morphological characteristics or molecular markers.

Assessments of honey bee genetic diversity have included morphological comparisons and the identification of molecular markers. Relying solely on morphological characteristics to differentiate A. cerana groups comes with limitations, but the use of molecular markers is currently recognized as a more reliable and accurate method (Li et al., 2019). Molecular genetic markers are one of the most powerful tools for genomic analysis, linking inherited traits to potential gene variants by locating and identifying heritable DNA sequences (Duran et al., 2009). The main molecular genetic markers commonly used in honey bees include restriction fragment length polymorphisms (RFLP), randomly amplified polymorphic DNA (RAPD), repeated simple sequences (SSR), single nucleotide polymorphisms (SNP), and mitochondrial DNA (mtDNA) markers (Qamer et al., 2021; Syromyatnikov et al., 2018; Wang et al., 2022).

At present, mtDNA markers are among the most frequently used molecular marker techniques. Due to their rapid rate of evolution, simple structural composition, and abundant genetic information, mitochondria have become crucial for studying the evolutionary origins and genetic diversity of organisms (Huang et al., 2021; Mandal et al., 2014; Wang et al., 2020). A. cerana exhibits a consistent matrilineal mitochondrial pattern, it is considered as an ideal organism for studying genetic diversity through the mitochondrial genome. The mitochondrial genome provides crucial information for mastering molecular evolution and phylogeny. The evolution of the A. mellifera mitochondrial genome has revealed genomic differences in A. mellifera at the subspecies level (Eimanifar et al., 2018). An evolutionary analysis of A. cerana from China and Nepal using COI-COII sequences showed that genetic differentiation was co-induced (Tan et al., 2016). Previous research has shown that among honey bees, the mitochondrial genome differs by species, subspecies, and geographical groups. The mitochondrial genome can be used to explore the genetic diversity and evolution of A. cerana populations (Tang et al., 2023).

In this study, we report the complete mitochondrial genome of A. cerana-Southern China. The analysis of the sequence and structure provides valuable insights into the possible mutations in the mitochondrial genome of honey bees to a large extent. It also sheds light on the evolution and environmental adaptation of A. cerana. This study can provide important references for the conservation and utilization of honey bee germplasm resources, as well as for the breeding and evolution of honey bees.

2 Materials and methods

2.1 Sample collection and DNA sequencing

The A. cerana-Southern China samples were obtained from the Hongchang Town, Chaonan District, Shantou, Guangdong Province, China (N 23.31, E 115.95). Fifty drone pupae were collected from one colony, then placed in absolute ethanol and stored in a −20 °C freezer. The genomic DNA was extracted and sent to Berry Genomics Co., Ltd. (Beijing, China) for PacBio HiFi sequencing.

2.2 Mitochondrial genome assembly and annotation

PacBio HiFi reads generated for A. cerana–Southern China were first used to assemble the nuclear genome, including mitochondrial sequences, using Hifiasm v0.19.5-r587 (Cheng et al., 2021; Cheng et al., 2022). The resulting contigs from the Hifiasm assembly were indexed to construct a local BLAST database. A reference A. cerana mitochondrial genome (GenBank accession NC_014295) was downloaded and used solely to identify mitochondrial contigs within the Hifiasm assembly by BLAST searches, rather than for genome assembly. The identified mitochondrial contigs were then extracted for downstream analyses. Preliminary annotation of the mitochondrial genome was performed using the MITOS Web Server (Bernt et al., 2013). The annotation results were subsequently refined and validated using GeSeq implemented in Chlorobox (Tillich et al., 2017) by comparison with previously reported protein-coding and RNA genes from related species. A circular map of the mitochondrial genome was generated using Chloroplot (Zheng et al., 2020). The secondary structures of tRNA genes were predicted using the MITOS Web Server under the invertebrate mitochondrial genetic code (Bernt et al., 2013), and the tRNA secondary structure diagrams were manually redrawn using Adobe Illustrator.

2.3 Mitochondrial genome sequence analysis

BioEdit (version 7.0.9.0) (Hall, 1999) was used to calculate the base composition and codon usage of the complete A. cerana-Southern China mitochondrial genome sequence and the protein-coding genes (PCGs), tRNA genes, rRNA genes and AT-rich regions. Codon usage bias, represented by relative synonymous codon usage (RSCU), was calculated using codonW software (retrieved from https://codonw.sourceforge.net/) and was plotted into a bar graph using JSHYCloud (http://cloud.genepioneer.com:9929). The RSCU greater than 1 indicated a preference for certain amino acids. The AT and GC biases were calculated according to the following formulas: AT skew = (A-T)/(A + T) and GC skew = (G-C)/(G + C).

2.4 Phylogenetic analyses

The mitochondrial genome of A. cerana–Southern China, assembled from PacBio HiFi sequencing data, has been deposited in the NCBI GenBank database under the accession number PP692293. In addition to this newly generated mitogenome, a total of 18 honey bee mitochondrial genomes were included in this study, comprising 15 previously published sequences retrieved from GenBank. The mitochondrial genomes of Apis mellifera sinisxinyuan and Apis mellifera capensis were selected as outgroups (Table 1). Because mitochondrial tRNA and rRNA genes are highly conserved in both nucleotide sequences and secondary structures, their inclusion provides limited phylogenetic resolution. Therefore, 13 mitochondrial protein-coding genes (PCGs) were used to infer the phylogenetic relationships among A. cerana populations from different geographic regions (Figure 1).

Table 1

Table 1. Information about the sequences of phylogenetic tree.

Figure 1

Figure 1. Geographic distribution of Apis cerana samples used for phylogenetic analysis.

The 13 PCG sequences were input into PhyloSuite (V1.2.3) software (Xiang et al., 2023; Zhang et al., 2020) and aligned using MAFFT. Then, the PCG sequences were concatenated. Modelfind was used for evolutionary model selection with a BIC value. The evolutionary tree was constructed using MrBayes and an IQ tree with default parameters.

2.5 Chromosome anchoring and mito–nuclear synteny analysis

The PacBio HiFi–assembled nuclear genome was anchored to chromosomes using RagTag (Alonge et al., 2022) against the Apis cerana chromosome-level reference (GCA_029169275.1), after which the complete mitochondrial genome served as the BLASTn query to assess mito-nuclear collinearity (synteny).

3 Results

3.1 Mitogenome organization and base composition

Using PacBio HiFi sequencing data, we assembled a high-quality mitochondrial genome of A. cerana–Southern China (GenBank accession PP692293). PacBio HiFi sequencing generated 64.44 Gb of high-quality long-read data from drone pupae, with an average read length of approximately 17.8 kb and a mapping ratio of 99.90%, providing sufficient coverage for reliable genome assembly. The Hifiasm assembly produced a high-quality nuclear genome (total length: 231.9 Mb) with 40 contigs and a contig N50 of 11.78 Mb, supported by a BUSCO completeness score of 97.90%.

We assembled the complete A. cerana–Southern China (GenBank: PP692293) mitochondrial genome (Figure 2). The complete mitochondrial genome is 16,137 bp in length and contains 37 genes, including 13 protein-coding genes (PCGs), 22 tRNA genes, two rRNA genes, and an AT-rich region (Table 2). Among these genes, 23 are located on the J strand (majority strand), whereas the remaining 14 are encoded on the N strand (minority strand) (Figure 2; Table 2). The nucleotide composition of the complete mitochondrial genome is A = 42.41%, T = 41.82%, C = 9.59%, and G = 6.18%, resulting in a high A + T content of 84.23%. Most protein-coding genes (PCGs) initiated with standard ATN start codons, while several genes used TTA as an alternative initiation codon. In several genes, including ND5, ND4, ND4L, and ND1, TTA is annotated as an alternative initiation codon, based on conserved open reading frames and comparative alignment with closely related Apis mitochondrial genomes.

Figure 2

Figure 2. Graphical maps of the mitochondrial genome of A. cerana-Southern China. The inside circles show the G + C content.

Table 2

Table 2. Annotations of the mitochondrial genome of A. cerana-Southern China.

Regarding termination codons, most PCGs terminate with the conventional mitochondrial stop codon TAA. In a small number of cases, truncated stop codons were inferred, which are presumed to be completed through post-transcriptional polyadenylation, a common mechanism in insect mitochondrial genomes. The intergenic regions of the A. cerana–Southern China mitochondrial genome range from 0 to 229 bp in length.

3.2 Protein-coding genes

The total length of the 13 PCGs in A. cerana-Southern China was 11,055 bp, accounting for 68.51% of the entire mitochondrial genome (Table 3). The length and base composition of the PCGs were varied, the 13 PCGs ranged from 162 bp (ATP8) to 1668 bp (ND5). Except for the A + T content of COX1 and COX2, the A + T content of the other PCGs was greater than 80%, indicating an AT preference. According to the AT-skew and GC-skew results, most of the PCGs were biased toward T, while all PCGs showed a preference for codon C. We conducted an analysis of the frequency and RSCU values for 13 PCGs in the A. cerana-Southern China mitochondrial genome. The RSCU values for the 13 PCGs were greater than 1, except for those of the two amino acids Met and Trp, which were each encoded by only one codon (Figure 2). This suggests that all codons were biased. The TTA codon had the highest RSCU value in A. cerana-Southern China (3.84), followed by the AGA codon (3.46). NNA-type codons accounted for 41.93% of the total number of codons, indicating that codons with an A at the third position were the most frequent. This study conducted a synteny analysis of the 13 protein-coding genes (PCGs) across 18 mitochondrial genomes (Figure 3). The results show that the order and strand orientation of the PCGs are highly conserved among genomes, with overall synteny similarity frequently reaching ≥0.95 and no large-scale rearrangements that would alter the ordering of coding regions. The differences observed are largely confined to noncoding hypervariable regions or a few tRNA clusters, manifested as changes in intergenic spacer length and minor shifts at annotation boundaries. Together with the previously noted AT bias and RSCU skew, these findings indicate that the architectural framework of the mitochondrial genome is broadly stable, whereas evolutionary changes at the sequence and codon-usage levels are more likely to drive differentiation among populations or geographic lineages. Accordingly, structural rearrangements are unlikely to be the primary mechanism underlying mitochondrial diversification in this species, and subsequent comparative and selection analyses can focus—without confounding from structural change—on nucleotide substitutions, codon-usage bias, and the ratio of nonsynonymous to synonymous substitutions (dN/dS). A comparative analysis of the lengths of the 13 mitochondrial protein-coding genes (PCGs) across 18 mitogenomes revealed highly consistent gene sizes among samples, with only minor variation detected in a few loci (notably ATP8, ND3, and ND6). ND5 was the longest gene, whereas ATP8 was the shortest, and the lengths of the COX subunits and CYTB were particularly invariant. These patterns indicate that the coding regions are under strong functional constraint and are overall conserved; the small length differences observed are likely attributable to short indels within AT-rich segments and to differences in start/stop codon usage or annotation boundaries. Consistent with our observations of AT bias and codon-usage skew, these results suggest that mitochondrial divergence in this taxon is driven more by subtle sequence and codon-level evolution than by changes in genome architecture (Figure 4).

Table 3

Table 3. Composition and skew in the mitochondrial genome of honey bee A. cerana-Southern China.

Figure 3

Figure 3. Synteny of 13 mitochondrial protein-coding genes across 18 genomes.

Figure 4

Figure 4. Length variation of the 13 mitochondrial protein-coding genes across 18 mitogenome.

3.3 rRNA and tRNA genes

The lrRNA and srRNA of A. cerana-Southern China were 1324 bp and 773 bp in length, respectively. lrRNA and srRNA were typically separated by trnV (Table 2). The A + T and G + C content of rRNA genes were 82.36% and 17.64%, and the AT skew and GC skew were −0.01 and −0.32, respectively, suggesting an apparent bias toward the use of T and C. Based on the length comparison (Figure 5), the srRNA is predominantly 772–787 bp and the lrRNA 1299–1337 bp, with minimal variation across the 18 mitogenomes, further supporting strong structural conservation of the rRNA region; the two rRNAs are typically separated by trnV.

Figure 5

Figure 5. Length comparison of 12S (srRNA) and 16S (lrRNA) across 18 mitogenomes.

In the mitochondrial genome of A. cerana (Southern China population), twenty-two tRNA genes were identified. Fourteen are encoded on the J-strand (major strand) and eight on the N-strand (minor strand). Except for trnS1, which lacks the dihydrouridine (DHU) arm, all tRNAs exhibit canonical cloverleaf secondary structures (Figures 6,7). The tRNA complement spans 1,468 bp in total, accounting for 9.10% of the mitogenome, with individual gene lengths ranging from 62 bp (trnS1) to 79 bp (trnP) (Table 3). Across genomes, tRNA lengths are tightly clustered—mostly 65–72 bp—with only minor 1-3 bp differences at a few loci, most notably in trn-Pro and trn-Thr; this pattern is consistent with small indels in AT-rich segments, subtle annotation-boundary shifts, and the missing DHU arm in trnS1 (Figure 8). These regions are strongly AT-rich (A + T = 87.35%, G + C = 12.65%; A + T is 6.91-fold higher than G + C). The AT-skew and GC-skew are 0.03 and −0.13, respectively, indicating a slight bias toward A over T and a pronounced bias toward C over G.

Figure 6

Figure 6. Length distribution of the 22 mitochondrial tRNA genes across 18 mitogenomes.

Figure 7

Figure 7. Codon bias statistics for 13 protein-coding genes of A. cerana-Southern China Ribosomal RNA genes and transfer RNA genes.

Figure 8

Figure 8. Predicted secondary structures for the tRNA genes in the A. cerana-Southern China mitochondrial genome.

3.4 AT-rich region

The AT-rich region of the A. cerana-Southern China mitochondrial genome extended over 818 bp and was located between srRNA and trnS1. The AT-rich region contained the highest A + T content (97.43%) in the entire mitogenome. Both the AT skew and GC skew for the AT-rich region were negative, indicating that T and C are more abundant than A and G, respectively (Table 3).

3.5 Phylogenetic analysis

Phylogenetic analysis was conducted using 13 PCGs with 17 closely related taxa (Figure 9). The results of the evolutionary tree analysis revealed that A. cerana was divided into four branches, one from Taiwan and one from Borneo, with high node support. A. cerana from Korea, Japan, and Russia formed a single colony, while those from mainland China formed an additional colony. A. cerana-Southern China and other A. cerana in mainland China represent sister groups. The genetic distance of A. cerana-Southern China was closer to that of A. cerana-Central China and A. cerana-Aba, and the differentiation time for A. cerana-South Yunnan was earlier.

Figure 9

Figure 9. Phylogenetic tree of honey bee, A. cerana populations. A. mellifera sinisxinyuan and A. mellifera capensis as outgroups. The numbers on the branch points of the phylogenetic tree represent bootstrap values for the tree or evolutionary branching length.

3.6 Genome-wide Distribution of Nuclear Mitochondrial DNA segments (NUMTs)

The circular Circos visualization shows mitochondrial homologous fragments mapping, via multiple links, to numerous loci on chr1–chr16, forming a multi-chromosomal, dispersed distribution rather than a continuous large syntenic block. The hits are predominantly short tracts (thin links with no extensive overlap), with relatively denser connections on chr1, chr6, and chr15–16. No long, one-to-one syntenic structure of mitochondrial coding regions was detected in the nuclear genome. This distribution pattern is consistent with nuclear mitochondrial DNA segments (NUMTs) (Figure 10).

Figure 10

Figure 10. Genome-wide distribution of nuclear mitochondrial DNA segments (NUMTs) shown by Circos.

4 Discussion

As populations of A. cerana continue to decline, it is vital to search for effective genetic markers to target and protect different types of A. cerana. Mitochondrial DNA is inherited maternally in bees, so a single bee can be used to represent an entire colony. Here, we report the complete mitochondrial genome of A. cerana-Southern China in Guangdong Province, China, which improve the data for mitochondrial molecular markers and comparative genomics of honey bee species.

Advances in sequencing technologies have enabled researchers to explore the complex genetic structure of A. cerana more accurately (Li et al., 2019; Okuyama et al., 2017). Compared with first-generation sequencing and next-generation sequencing, third-generation sequencing provides longer read lengths, which offer clear advantages for specific applications (Fan et al., 2020). In this study, we utilized PacBio HiFi sequencing to assemble the complete mitochondrial genomes of A. cerana-Southern China. Sequencing depth across the mitogenome was sufficient and broadly uniform, with an average coverage of 57.69× (range 30×–80×), supporting a high-confidence assembly (Supplementary Figure S1). Local fluctuations near AT-rich segments are expected for insect mitogenomes and did not affect gene recovery. According to our results, the full length of the A. cerana-Southern China mitochondrial genome was 16,137 bp, which is consistent with the previously reported total length of animal mitochondrial genomes ranging from 13 to 20 kb (Cameron, 2014). The mitochondrial genomes contain 37 common genes, including 13 PCGs, 22 tRNA genes, and 2 rRNA genes. All 37 mitochondrial genes exhibit a reversal of strand asymmetry on the majority strand, possibly due to the inversion of the replication origin located in the control region (Wei et al., 2010). tRNA genes and PCGs rearrangement have been previously reported in Hymenoptera (Su et al., 2023). However, this study of A. cerana-Southern China showed no rearrangement, which is considered to indicate the ground pattern of insect mitochondrial genomes. In most metazoan mitochondria, trnS1 generally has an unpaired DHU arm, while trnS2 has a standard cloverleaf structure. The A. cerana-Southern China trnS1 and trnS2 structures conform to this interpretation (Tan et al., 2011). Whether sequences tend to be A-encoded or T-encoded depends on the role they play. Strand asymmetry, also known as strand compositional bias, is typically indicated by the AT and GC skews (Perna and Kocher, 1995). We performed a codon preference analysis on 13 PCGs. The preference for the utilization of these synonymous codons is determined by several factors, including the abundance of tRNA, the mutational bias of the gene chain, gene expression level, gene length, and GC composition (Li et al., 2023). Our results show that the mitochondrial genome of A. cerana-Southern China favors NNA-type codons, which is consistent with the results of previous studies on Apis (Tang et al., 2023). Pairwise dN/dS (Ka/Ks) estimates for the 13 PCGs were all <1, indicating pervasive purifying selection (Supplementary Figure S2). Relatively elevated ratios were observed for ND1, ND3, and CYTB compared with COX genes, suggesting heterogeneous selective constraints among complexes.

Mitochondrial molecular markers provide a better understanding of the genetic structure and evolutionary history of honey bee populations, so they are widely used to study the biological origin, evolution, and genetic diversity of honey bees (Feng et al., 2024). The combination of sequencing technology and bioinformatics methods has advanced technological progress in the study of bee genome variation and environmental adaptive genome characteristics. At present, the genes commonly used for sequencing and analysis in honey bee mitochondrial genome research include tRNA^leu-COII, COII, lrRNA, CYTB and ND2. In addition, we found that the length of the AT-rich region of A. cerana is very different among different groups, and the AT-rich region may be a better subspecies division as a mitochondrial molecular marker fragment (Zhang et al., 2023). This region can also reflect the adaptive evolution of A. cerana to their living environment, and researchers are investigating the impact of this complex region on the genetic differentiation of species (Liu et al., 2014).

Phylogenetic analysis was conducted using the 13 PCGs with 17 closely related taxa. The results further suggest that the current subspecies structure is the result of multiple population divisions that occurred in a common ancestor. The results showed that all genome sequences of A. cerana were mainly divided into four branches. A. cerana from Taiwan and Borneo form their own branch. The sequences from Russia, Japan and Korea were clustered in a single group, which was consistent with a previous report (Xu et al., 2024). The genetic distance between A. cerana-Southern China, A. cerana-Central China and A. cerana-Aba close. This was followed by A. cerana-South Yunnan and A. cerana-Yungui Plateau. However, compared to these types of A. cerana, the differentiation time of A. cerana-South Yunnan was the earliest. However, there appeared to be no correlation between genetic diversity and geographical distance. Morphometric analysis revealed that A. cerana exhibits high diversity in mountainous regions and islands, a consequence of selection and genetic drift due to prolonged isolation (Ji et al., 2023; Yu et al., 2019). Previous study on honeybees in various regions, including an analysis of morphological differences between A. cerana-South Yunnan and A. cerana-Yungui Plateau, have indicated that the Tropic of Cancer serves as a geographic boundary (Zhang et al., 2014). The analysis of the correlation between the morphological indexes of A. cerana and the climate factors of their living environment found that most of the classical morphological characteristics of A. cerana, such as body size and wing length, were significantly correlated with climate characteristics (Qiu et al., 2023). This suggests that there are differences in the genetic distance of A. cerana, even when the physical distance is close.

The populations of A. cerana are prone to genetic differentiation and exhibit a high level of diversity. In contrast to previous studies that provided only limited information about A. cerana at the genetic level, the mitochondrial genome-wide analysis of A. cerana in our study offers a more detailed genetic structure of honeybees. The effective population size and genetic diversity of A. cerana in China are high, making the populations and the species within this range in general vulnerable to the effects of genetic drift and selection. Due to human and environmental damage, this could lead to the extinction of local populations. Therefore, further research is required to determine the optimal approach for conserving the genetic resources of specific A. cerana populations. In this paper, the study of the mitochondrial genome of A. cerana may be limited. Its genetic information needs further exploration, possibly in conjunction with other data, such as bee morphology in different environments. This approach would lead to a more comprehensive understanding of its biological characteristics and evolutionary history.

5 Conclusion

In conclusion, we assembled a high-quality mitochondrial genome of A. cerana-Southern China using PacBio HiFi sequencing technology. We have detailed the mitochondrial genome information of A. cerana-Southern China. Phylogenetic studies have indicated that the mitochondrial genome is a reliable molecular marker for studying the phylogeny of A. cerana. By comparing the mitochondrial genomes, we found that A. cerana exhibits high genetic diversity at the mitochondrial genomic level across different environments. The results of the phylogenetic tree indicated that the differentiation time of A. cerana-Southern China was the latest compared to other sister groups. This study provided data for the gene bank of A. cerana. It is hoped that more specimens can be collected to conduct a more comprehensive and systematic study of honeybees. The morphological data and gene fragments can be utilized for combined analyses in the future.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/genbank/, PP692293.

Ethics statement

Ethical approval was not required for the study involving animals in accordance with the local legislation and institutional requirements because these were wild honey bees; no ethical approval was required.

Author contributions

XD: Conceptualization, Data curation, Formal Analysis, Funding acquisition, Investigation, Writing – original draft, Writing – review and editing. XY: Data curation, Supervision, Validation, Writing – review and editing. JC: Data curation, Project administration, Supervision, Writing – review and editing. RS: Resources, Writing – original draft, Conceptualization, Data curation, Methodology, Project administration, Software, Validation, Writing – review and editing. JL: Resources, Writing – original draft, Formal Analysis, Funding acquisition, Visualization.

Funding

The author(s) declared that financial support was received for this work and/or its publication. This research is supported by the Yiwu City Science and Technology Bureau and Yiwu Industrial and Commercial College Program under Grant Nos. XM2025DCZ0227 and XJKJ2502YB. This work was supported by Jinhua Key Laboratory of Robot Intelligent Welding Technology and by Jinhua Public Welfare Technology Application Research Project in 2024. General project of Zhejiang Province Philosophy and Social Science Planning Project: Analysis of Online Learning Behavior and Intervention Pathways for Multimodal Data Fusion (22NDJC194YB). Yiwu Science and Technology Plan Project: Research and Monitoring of Microbial Detection Methods Based on TDLAS Technology, Project Number: 25-1-015.

Acknowledgements

Thanks to Anhui Gaohe Biotechnology Co., Ltd. for providing sequencing and analysis.

Conflict of interest

The author(s) declared that this work was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Generative AI statement

The author(s) declared that generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2025.1737945/full#supplementary-material

Footnotes

^{Abbreviations:}mt, Mitochondrial; PCGs, Protein-coding genes; tRNA, genes Transfer RNA genes; rRNA, genes Ribosomal RNA genes.

References

Alonge, M., Lebeigle, L., Kirsche, M., Jenike, K., Ou, S., Aganezov, S., et al. (2022). Automated assembly scaffolding using RagTag elevates a new tomato system for high-throughput genome editing. Genome Biol. 23, 258. doi:10.1186/s13059-022-02823-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Bernt, M., Donath, A., Jühling, F., Externbrink, F., Florentz, C., Fritzsch, G., et al. (2013). MITOS: improved de novo metazoan mitochondrial genome annotation. Mol. Phylogenetics Evol. 69 (2), 313–319. doi:10.1016/j.ympev.2012.08.023

PubMed Abstract | CrossRef Full Text | Google Scholar

Cameron, S. L. (2014). Insect mitochondrial genomics: implications for evolution and phylogeny. Annu. Rev. Entomol. 59 (1), 95–117. doi:10.1146/annurev-ento-011613-162007

PubMed Abstract | CrossRef Full Text | Google Scholar

Cao, M. F., Tang, L., Chen, J., Zhang, X. Y., Easy, R. H., and You, P. (2020). The mitogenome of freshwater loach Homatula laxiclathra (teleostei: nemacheilidae) with phylogenetic analysis of nemacheilidae. Ecol. Evol. 10 (12), 5990–6000. doi:10.1002/ece3.6338

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, C., Liu, Z. G., Luo, Y. X., Xu, Z., Wang, S. H., Zhang, X. W., et al. (2017). Managed honeybee colony losses of the eastern honeybee (apis Cerana) in China (2011-2014). Apidologie 48 (5), 692–702. doi:10.1007/s13592-017-0514-6

CrossRef Full Text | Google Scholar

Cheng, H. Y., Concepcion, G. T., Feng, X. W., Zhang, H. W., and Li, H. (2021). Haplotype-resolved de novo assembly using phased assembly graphs with hifiasm. Nat. Methods 18 (2), 170–175. doi:10.1038/s41592-020-01056-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Cheng, H. Y., Jarvis, E. D., Fedrigo, O., Koepfli, K. P., Urban, L., Gemmell, N. J., et al. (2022). Haplotype-resolved assembly of diploid genomes without parental data. Nat. Biotechnol. 40 (9), 1332–1335. doi:10.1038/s41587-022-01261-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Duran, C., Appleby, N., and David, E. (2009). Molecular genetic markers: discovery, applications, data storage and visualisation. Curr. Bioinform 4, 16–27. doi:10.2174/157489309787158198

CrossRef Full Text | Google Scholar

Eimanifar, A., Kimball, R. T., Braun, E. L., and Ellis, J. D. (2018). Mitochondrial genome diversity and population structure of two western honey bee subspecies in the Republic of South Africa. Sci. Rep. 8 (1), 1333. doi:10.1038/s41598-018-19759-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Fan, X. Y., Tang, D., Liao, Y. H., Li, P. D., Zhang, Y., Wang, M. X., et al. (2020). Single-cell RNA-seq analysis of mouse preimplantation embryos by third-generation sequencing. PLoS Biol. 18 (12), e3001017. doi:10.1371/journal.pbio.3001017

PubMed Abstract | CrossRef Full Text | Google Scholar

Feng, G. Y., Jiao, Y. J., Ma, H. Z., Bian, H. Y., Nie, G., Huang, L. K., et al. (2024). The first two whole mitochondrial genomes for the genus Dactylis species: assembly and comparative genomics analysis. BMC Genomics 25 (1), 235. doi:10.1186/s12864-024-10145-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Hall, T. A. (1999). BioEdit: a user-friendly biological sequence alignment editor and analysis program for windows 95/98/NT. Nucleic Acids Symp. Ser. 41, 95–98. doi:10.14601/Phytopathol_Mediterr-14998u1.29

CrossRef Full Text | Google Scholar

Huang, Y. X., Wang, S., Gao, Y. Q., Chen, J. H., Wang, X. L., and Li, R. J. (2021). Comparison of mitochondrial genome and development of specific PCR primers for identifying two scuticociliates, Pseudocohnilembus persalinus and Uronema marinum. Parasites Vectors 14 (1), 318. doi:10.1186/s13071-021-04821-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Ji, C. C., Shi, W., Tang, J., Gao, J. L., Liu, F., Shan, J. Q., et al. (2023). Morphometrical analyses revealed high diversity of the eastern honey bee (Apis cerana) in mountains and islands in China. J. Apic. Res. 62 (4), 647–655. doi:10.1080/00218839.2023.2205670

CrossRef Full Text | Google Scholar

Li, Y. C., Chao, T. L., Fan, Y. H., Lou, D. L., and Wang, G. Z. (2019). Population genomics and morphological features underlying the adaptive evolution of the eastern honey bee (Apis cerana). BMC Genomics 20 (1), 869. doi:10.1186/s12864-019-6246-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, G. L., Zhao, H., Guo, D. Z., Liu, Z. G., Wang, H. F., Sun, Q. H., et al. (2022). Distinct molecular impact patterns of abamectin on Apis mellifera ligustica and Apis cerana cerana. Ecotoxicol. Environ. Saf. 232, 113242. doi:10.1016/j.ecoenv.2022.113242

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, J. X., Chen, Y. L., Liu, Y. L., Wang, C., Li, L., and Chao, Y. H. (2023). Complete mitochondrial genome of agrostis stolonifera: insights into structure, codon usage, repeats, and RNA editing. BMC Genomics 24 (1), 466. doi:10.1186/s12864-023-09573-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, J., Bu, C. P., Wipfler, B., and Liang, A. P. (2014). Comparative analysis of the mitochondrial genomes of callitettixini spittlebugs (hemiptera: cercopidae) confirms the overall high evolutionary speed of the AT-rich region but reveals the presence of short conservative elements at the tribal level. PLoS One 9 (10), e109140. doi:10.1371/journal.pone.0109140

PubMed Abstract | CrossRef Full Text | Google Scholar

Mandal, S. D., Chhakchhuak, L., Gurusubramanian, G., and Kumar, N. S. (2014). Mitochondrial markers for identification and phylogenetic studies in insects-A review. DNA Barcodes Berl. 2 (1). doi:10.2478/dna-2014-0001

CrossRef Full Text | Google Scholar

Matias, D. M. S., Leventon, J., Rau, A., Borgemeister, C., and von Wehrden, H. (2017). A review of ecosystem service benefits from wild bees across social contexts. Ambio 46 (4), 456–467. doi:10.1007/s13280-016-0844-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Okuyama, H., Wakamiya, T., Fujiwara, A., Washitani, I., and Takahashi, J. I. (2017). Complete mitochondrial genome of the honeybee Apis cerana native to two remote islands in Japan. Conserv. Genet. Resourv 9, 557–560. doi:10.1007/s12686-017-0721-5

CrossRef Full Text | Google Scholar

Perna, N. T., and Kocher, T. D. (1995). Patterns of nucleotide composition at fourfold degenerate sites of animal mitochondrial genomes. J. Mol. Evol. 41 (3), 353–358. doi:10.1007/BF00186547

PubMed Abstract | CrossRef Full Text | Google Scholar

Qamer, S., Al-Abbadi, A. A., Sajid, M., Asad, F., Khan, M. F., Khan, N. A., et al. (2021). Genetic analysis of honey bee, Apis dorsata populations using random amplified polymorphic DNA (RAPD) markers. J. King Saud. Univ. Sci. 33 (1), 101218. doi:10.1016/j.jksus.2020.10.015

CrossRef Full Text | Google Scholar

Qiu, L. F., Dong, J. X., Li, X. G., Parey, S. H., Tan, K., Orr, M., et al. (2023). Defining honeybee subspecies in an evolutionary context warrants strategized conservation. Zool. Res. 44 (3), 483–493. doi:10.24272/j.issn.2095-8137.2022.414

PubMed Abstract | CrossRef Full Text | Google Scholar

Su, C. Y., Zhu, D. H., Abe, Y., Ide, T., and Liu, Z. W. (2023). The complete mitochondrial genome and gene rearrangements in a gall wasp species, Dryocosmus liui (hymenoptera: cynipoidea: cynipidae). PeerJ (San Francisco, CA) 11, e15865. doi:10.7717/peerj.15865

PubMed Abstract | CrossRef Full Text | Google Scholar

Syromyatnikov, M., Borodachev, A., Kokina, A., and Popov, V. (2018). A molecular method for the identification of honey bee subspecies used by beekeepers in Russia. Insects 9 (1), 10. doi:10.3390/insects9010010

PubMed Abstract | CrossRef Full Text | Google Scholar

Tan, H. W., Liu, G. H., Dong, X., Lin, R. Q., Song, H. Q., Huang, S. Y., et al. (2011). The complete mitochondrial genome of the asiatic cavity-nesting honeybee Apis cerana (hymenoptera: apidae). PLoS One 6 (8), e23008. doi:10.1371/journal.pone.0023008

PubMed Abstract | CrossRef Full Text | Google Scholar

Tan, K., Qu, Y. F., Wang, Z. W., Liu, Z. W., and Engel, M. S. (2016). Haplotype diversity and genetic similarity among populations of the eastern honey bee from himalaya-southwest China and Nepal (hymenoptera: apidae). Apidologie 47 (2), 197–205. doi:10.1007/s13592-015-0390-x

CrossRef Full Text | Google Scholar

Tang, X. Y., Yao, Y. X., Li, Y. H., Song, H. L., Luo, R., Shi, P., et al. (2023). Comparison of the mitochondrial genomes of three geographical strains of Apis laboriosa indicates high genetic diversity in the black giant honeybee (hymenoptera: apidae). Ecol. Evol. 13 (2), e9782. doi:10.1002/ece3.9782

PubMed Abstract | CrossRef Full Text | Google Scholar

Tillich, M., Lehwark, P., Pellizzer, T., Ulbricht-Jones, E. S., Fischer, A., Bock, R., et al. (2017). GeSeq-versatile and accurate annotation of organelle genomes. Nucleic. acids. Res. 45 (W1), W6–W11. doi:10.1093/nar/gkx391

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, T., Shao, W. Z., Wu, Z. Y., Ren, Y. L., Yang, M. F., and Wu, G. H. (2020). Complete mitochondrial genome of Bactrocera proprediaphora using next-generation sequencing from China and its phylogenetic analysis. Mitochondrial DNA Part B-Resour 5 (2), 1596–1597. doi:10.1080/23802359.2020.1742621

CrossRef Full Text | Google Scholar

Wang, Z., Chen, W. F., Wang, H. F., Xv, B. H., and Liu, Z. G. (2022). Research progress on molecular markers in honeybee. J. Bee 42 (7), 33–37.

Google Scholar

Wei, S. J., Shi, M., Chen, X. X., Sharkey, M. J., van Achterberg, C., Ye, G. Y., et al. (2010). New views on strand asymmetry in insect mitochondrial genomes. PLoS One 5 (9), e12708. doi:10.1371/journal.pone.0012708

PubMed Abstract | CrossRef Full Text | Google Scholar

Xiang, C. Y., Gao, F., Jakovlić, I., Lei, H. P., Hu, Y., Zhang, H., et al. (2023). Using PhyloSuite for molecular phylogeny and tree-based analyses. iMeta 2 (1), 2–87. doi:10.1002/imt2.87

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, L. Q., Yang, J. D., and Lai, K. (2024). The complete mitochondrial genome of the cavity-nesting honeybee, Apis cerana abansis (insecta: hymenoptera: apidae). Cytol. Genet. 58, 136–141. doi:10.3103/S0095452724020117

CrossRef Full Text | Google Scholar

Yu, Y. L., Zhou, S. J., Zhu, X. J., Xu, X. J., Wang, W. F., Zha, L., et al. (2019). Genetic differentiation of eastern honey bee (Apis cerana) populations across Qinghai-Tibet plateau-valley landforms. Front. Genet. 10, 483. doi:10.3389/fgene.2019.00483

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Z. Y., Luo, W. T., Song, W. F., Liang, C., and Zhang, X. W. (2014). Morphology of Apis cerana in the Yunnan-Guizhou plateau and bees in south of Yunnan. Apic. China 65 (1), 12–15. doi:10.1051/apido:2003049

CrossRef Full Text | Google Scholar

Zhang, D., Gao, F. L., Jakovlic, I., Zou, H., Zhang, J., Li, W. X., et al. (2020). PhyloSuite: an integrated and scalable desktop platform for streamlined molecular sequence data management and evolutionary phylogenetics studies. Mol. Ecol. Resour. 20 (1), 348–355. doi:10.1111/1755-0998.13096

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, C. H., Wang, Y. L., Chen, H. W., and Huang, J. (2023). Comparative mitochondrial genomes between the genera amiota and phortica (diptera: drosophilidae) with evolutionary insights into D-Loop sequence variability. Genes 14 (6), 1240. doi:10.3390/genes14061240

PubMed Abstract | CrossRef Full Text | Google Scholar

Zheng, S. Y., Poczai, P., Hyvönen, J., Tang, J., and Amiryousefi, A. (2020). Chloroplot: an online program for the versatile plotting of organelle genomes. Front. Genet. 11, 576124. doi:10.3389/fgene.2020.576124

PubMed Abstract | CrossRef Full Text | Google Scholar