Mitochondrial genomes of blister beetles (Coleoptera, Meloidae) and two large intergenic spacers in Hycleus genera
BMC Genomics volume 18, Article number: 698 (2017)
Insect mitochondrial genomes (mitogenomes) exhibit high diversity in some lineages. The gene rearrangement and large intergenic spacer (IGS) have been reported in several Coleopteran species, although very little is known about mitogenomes of Meloidae.
We determined complete or nearly complete mitogenomes of seven meloid species. The circular genomes encode 13 protein-coding genes (PCGs), 22 transfer RNAs (tRNAs) and two ribosomal RNAs (rRNAs), and contain a control region, with gene arrangement identical to the ancestral type for insects. The evolutionary rates of all PCGs indicate that their evolution is based on purifying selection. The comparison of tRNA secondary structures indicates diverse substitution patterns in Meloidae. Remarkably, all mitogenomes of the three studied Hycleus species contain two large intergenic spacers (IGSs). IGS1 is located between trnW and trnC, including a 9 bp consensus motif. IGS2 is located between trnS2 (UCN) and nad1, containing discontinuous repeats of a pentanucleotide motif and two 18-bp repeat units in both ends. To date, IGS2 is found only in genera Hycleus across all published Coleopteran mitogenomes. The duplication/random loss model and slipped-strand mispairing are proposed as evolutionary mechanisms for the two IGSs (IGS1, IGS2). The phylogenetic analyses using MrBayes, RAxML, and PhyloBayes methods based on nucleotide and amino acid datasets of 13 PCGs from all published mitogenomes of Tenebrionoids, consistently recover the monophylies of Meloidae and Tenebrionidae. Within Meloidae, the genus Lytta clusters with Epicauta rather than with Mylabris. Although data collected thus far could not resolve the phylogenetic relationships within Meloidae, this study will assist in future mapping of the Meloidae phylogeny.
This study presents mitogenomes of seven meloid beetles. New mitogenomes retain the genomic architecture of the Coleopteran ancestor, but contain two IGSs in the three studied Hycleus species. Comparative analyses of two IGSs suggest that their evolutionary mechanisms are duplication/random loss model and slipped-strand mispairing.
The animal mitochondrial genome (mitogenome) is an informative model for phylogenetics, molecular evolution, and comparative genomic research, due to its simple genetic structure, maternal inheritance, and high evolutionary rate properties [1,2,3]. The mitogenome has been increasingly used to analyze the phylogeny and evolution of the highly diverse and rapidly radiating insects [4,5,6,7,8]. The typical mitogenome of metazoans is a circular molecule encoding a conserved set of 37 genes for 13 protein-coding genes (PCGs), two ribosomal RNA (rRNA) genes, and 22 transfer RNA (tRNA) genes, and comprising a non-coding control region . The mitogenome commonly displays exceptional economy of organization evidenced by lacking introns, few intergenic spacers, incomplete stop codons, and even overlapping genes . However, the extremely diverse Insecta also exhibits high diversity in their mitochondrial genomes, such as gene rearrangements and/or long non-coding regions (except the control region) in some lineages within Hymenoptera [10,11,12,13,14], Hemiptera [15,16,17], Dictyoptera [18, 19], Diptera , Orthoptera , Thysanoptera , Psocoptera , and Phthiraptera . Within the order Coleoptera , most beetles retain the same gene content and gene organization as the hypothesized ancestral mitogenome for Insecta , while a few possess gene rearrangements of tRNAs [7, 27, 28] and PCGs  in their mitogenomes. For example, the mitogenome of the firefly, Pyrocoelia rufa contains a 1724-bp large intergenic spacer (IGS) composed of tandem repeat units . Non-coding intergenic spacers were also found in the mitogenomes of other beetles (Trachypachus holmbergi, Priasilpha obscura, and Cyphon sp.) . Although a recent study reported that Hycleus chodschenticus has a large intergenic spacer , very little is known about the rearrangement and large non-coding regions in Meloidae mitogenomes.
Meloidae is a medium sized family within Coleoptera Tenebrionoidea, containing more than 3000 species within approximately 125 genera . Meloids are commonly referred to as blister beetles due to a defensive secretion, cantharidin. Cantharidin is an intoxicant that can be used for the removal of warts and may be effective in the treatment of primary liver cancer, leucocytopenia, chronic liver disease, neurodermatitis and other major illnesses [33, 34]. The medicinal properties of cantharidin, hypermetamorphosis development, and parasitoid habit of Meloidae species have been extensively researched [35, 36]. However, the mitogenomes of Meloidae species is less well researched. Two of the few studies of Meloidae mitogenomes found that Epicauta chinensis and Lytta caraganae retained the ancestral model of the insect mitogenome, without any gene rearrangement or long non-coding regions [31, 37]. Of the approximately 3000 species within Meloidae, only three complete mitogenomes (of Epicauta chinensis, Lytta caraganae, and Hycleus chodschenticus) have been described (excluding another without any description; Mylabris sp., JX412732.1). The lack of research considerably limits the genomic comparisons and molecular phylogenetic studies of Meloidae. Thus we believed there was an urgent need to explore the mitogenome evolution in the diverse family of Meloidae.
Consequently, we determined the mitogenomes of seven blister beetle species, representing four meloid genera. These species were E. gorhami, E. tibialis, L. caraganae, Mylabris aulica, H. phaleratus, H.marcipoli, and H. cichorii. We described the general features of the newly sequenced mitogenomes from the seven species and analyzed two IGSs in all Hycleus species to explore their evolutionary mechanisms. In addition, we attempted to assess the possibility of the IGS2 to be a molecular marker and indicators of phylogenetic relationships within Meloidea, based on mitogenomic datasets. The mitogenomes of the seven meloids will significantly add to the knowledge of Meloidae taxonomy, phylogeny, and evolution.
Genome content and gene organization
This study presents six complete mitogenomes and one nearly complete mitogenome (H. cichorii) with the absence of the control region and three tRNAs (trnI, trnQ and trnM). The total lengths of complete mitogenomes range from 15,633 to 16,003 bp. All seven new sequences were submitted to GenBank under the accession numbers listed in Table 1. Every mitogenome of the seven meloid species is a circular DNA molecule, encoding the typical 37 genes including 13 PCGs, 22 tRNAs, two rRNAs, and a putative control region. The major strand (J strand) carries most of the genes (9 PCGs and 14 tRNAs), while the remaining genes are encoded on the minor strand (N strand) (Fig. 1).
All PCGs of the seven mitogenomes use typical ATN start codons. Conventional stop codons TAA/TAG are assigned to most of the PCGs, but the cox1, cox2, nad5 and nad4 genes of all meloids terminate with the incomplete stop codon T. This terminator is adopted by cox3 genes of Epicauta and Lytta species (Additional files 1, 2, 3, 4, 5, 6 and 7: Tables S1–S7). The A + T - rich regions of meloid mitogenomes range from 1015 to 1201, with the location between rrnS and trnI (Table 2). The poly-T stretch (15 bp) was detected in control regions of all meloids, but without tandem repeats.
A + T content and codon usage
The overall A + T contents of the seven meloid mitogenomes range from 67.53% to 71.94% (Table 1), and such an A + T bias is reflected in the codon frequencies of these mitogenomes (Fig. 2). Relatively synonymous codon usages (RSCU) were calculated over all seven meloid species, excluding stop codons (Additional file 8: Table S8). The RSCU demonstrate that codons with A or T in the third position are always overused as compared to other synonymous codons. Additionally, codons TTT (Phe), TTA (Leu), ATT (Ile), and ATA (Met) are the four most frequently used codons in these meloid mitogenomes. These codons are all comprised of A or T nucleotides, which is indicative of the biased usage of A and T nucleotides in the meloid mitochondrial PCGs.
Evolutionary rates of PCGs
All available mitogenomes were used to assess the evolutionary rate of PCGs for Meloidae. The variable sites, nucleotide diversity (π), and the ratio of non-synonymous substitution (Ka) to synonymous substitution (Ks) were calculated for each PCG (Table 3). Nad6 and atp8 genes exhibit the highest level of nucleotide diversities, whereas cox1 gene is the most conserved. The Ka/Ks value is correlated with the nucleotide diversity, with the highest level in nad6 and atp8 and the lowest in cox1. Notably, the Ka/Ks ratio for every PCG is lower than 1, indicating that all PCGs are evolving under the purifying selection.
Comparison of tRNA secondary structures
All seven meloid mitogenomes encode 14 tRNAs on the J strand and 8 on the N strand (Fig. 1, Table 2). The comparative results of tRNA secondary structures are provided in Figs. 3 and 4. All tRNAs could be folded into the typical clover-leaf structure, except trnS1 (AGN) lacks a dihydrouridine (DHU) arm, which is replaced by a simple loop (Fig. 4).
The nucleotide conservation of tRNAs is markedly J strand-biased, trnA, trnE, trnG and trnK, which have the highest percentage of identical nucleotides, are all located on the J strand. The tRNAs on the N strand (i.e., trnQ, trnH, trnF) and close to the control regions (i.e., trnI and trnQ) have high level of nucleotide variation (Additional file 9: Table S9).
The mitogenomes of E. gorhami, E. tibialis, and L. caraganae contain four intergenic spacers with a total length of 23, 25, 26 bp, respectively, while the M. aulica has six intergenic spacers with a total length of 28 bp (Table 2). Unexpectedly, H. phaleratus, H. marcipoli, and H. cichorii contain 7, 8, 8 intergenic spacers with a total length of 266, 203, 175 bp, respectively (Table 2), including two large intergenic spacers (>30 bp). IGS1 is located between trnW and trnC with a length of 62, 56, 39 bp in H. phaleratus, H. marcipoli, and H. cichorii, respectively (Table 2). The IGS1 sequences exhibit relatively high similarity among the three Hycleus species, and a 9 bp long congruent motif AAATTATGG was detected in the three Hycleus species (Fig. 5a). IGS2 is located between trnS2 and nad1 with a length of 181, 123, and 110 bp in H. phaleratus, H. marcipoli, and H. cichorii, respectively (Table 2). The alignment of IGS2s among all sequenced Hycleus mitogenomes, including the recently published sequence of H. chodschenticus , shows that a pentanucleotide motif (TACTA) exists in these Hycleus species. Furthermore, H. chodschenticus, H. phaleratus, H. marcipoli, and H. cichorii include five, four, three, and two repeats (respectively) of this motif (Fig. 5b). The organization of IGS2 indicates that these four Hycleus species contain two copies of an 18 bp conserved sequence (ATACTAAAYTTTRTTAAC) in both ends of IGS2 (Fig. 5b), but other meloid beetles have only one (Fig. 6a).
We carried out MrBayes, RAxML, and PhyloBayes analyses based on nucleotide and amino acid datasets to determine the influence of different datasets and analytical methods on tree topology and node reliability. Bayesian Inference (BI) and Maximum Likelihood (ML) analyses used the same datasets to generate congruent tree topologies. BI trees had higher node support values than ML trees (Fig. 7). PhyloBayes analyses generated different tree topologies with polytomies (Additional file 10: Figure S1). Four tree topologies derived from our six phylogenetic trees are consistent with the monophylies of Meloidae and Tenebrionidae, and the basal position of Mordellidae. Differences are present in the inter-family relationships of Aderidae, Ciidae, Oedemeridae, and Prostomidae.
Within Meloidae, L. caraganae is sister taxon to the species belonging to the genus Epicauta. The Meloidae family results monophyletic and receives maximum statistical support. In Tenebrionidae, the Tenebrioninae and Diaperinae are never recovered as monophyletic groups (Fig. 7, Additional file 10: Figure S1).
All seven mitochondrial genomic arrangements share the ancestral type for insects [26, 38], as is reported in the published meloid mitogenomes [31, 37]. The nucleotide composition of mitochondrial genomes for meloids also corresponds well to the A + T bias generally observed in insect mitogenomes.
All mitogenomes of the seven meloids have incomplete stop codons, which have been described in many other insect species [39, 40]. It has been demonstrated that incomplete stop codons can produce functional stop codons in polycistronic transcription cleavage and polyadenylation processes . Remarkably, the cox3 genes of Lytta and Epicauta species possess the same incomplete stop codon, while Hycleus and Mylabris beetles utilize complete terminators (Additional files 1, 2, 3, 4, 5, 6 and 7: Tables S1–S7). The similar preference for the adoption of stop codons seems to suggest that the genus Lytta is more closely related to Epicauta than the two other genera, and this relationship was confirmed by phylogenetic results.
The evolutionary rates of all mitochondrial PCGs indicate that their evolution is based on purifying selection (Table 3), as is reported in other insects [41, 42]. However, the cytochrome oxidase subunits (cox1, cox2, and cox3) and cytochrome b (cob) have lower Ka/Ks ratios than ATPase subunits (atp8 and atp6) and NADH dehydrogenase subunits (nad1–6 and 4 L). The nucleotide diversity also shows cox and cob genes are obviously more conserved than atp and nad genes. This phenomenon indicates that the various functional genes in the mitochondria of meloids underwent different selection pressures during evolution. Furthermore, the cox1 has the slowest evolutionary rate, demonstrating that functional constraints are more powerful for this gene than positive selection.
The absence of a DHU arm of trnS1 commonly exists in many metazoan mitogenomes, including insects [30, 43, 44]. However, this tRNA (missing DHU arm) was often suggested by authors to be functional [45, 46]. Another unusual feature is the use of TCT as the trnS1 anticodon in meloids, although most arthropods adopt a GCT anticodon in trnS1. This exceptional trnS1 anticodon was also found in many other beetles [5, 30, 47].
Some mismatched pairs in stems of tRNAs (e.g., T-T in the DHU stem of trnQ and in anticodon stem of trnK; C-T in the TΨC stem of trnI; A-C in anticodon stem of trnL2 (UUR); A-G in acceptor stem of trnW), are common in insect mitogenomes and can be corrected through editing processes or may represent unusual pairings . It was not possible to model the substitution pattern of the TΨC loop in trnH (Fig. 3) because of the high variation among orthologous sequences. The increasing variation usually accompanies more compensatory base changes in stems, resulting in the tRNA more or less not conserved (Additional file 9: Table S9).
The ends of rRNA genes of meloid mitogenomes were assumed to extend to the boundaries of flanking genes because it is impossible to accurately determine by DNA sequence alone . Consequently, rrnL was assumed to fill up the blank between trnV and trnL1 (CUN) (Fig. 1), but the boundary between the rrnS and the putative control region was defined by the alignment with homologous sequences .
The control region in the insect mitogenome is equivalent to the control region of vertebrate mitogenomes, which contains the origin sites for transcription and replication [50, 51]. The six complete mitogenomes include a poly-T stretch (15 bp) that was suggested to function as a possible recognition site for the initiation of replication of the mitochondrial DNA N strand . Like other Coleopteran mitogenomes, the control regions of meloids also exhibit the highest A + T content in the whole mitogenome. This region is unlikely to be more variable than protein-coding genes due to such high A + T content and consequently limits its usefulness as a molecular marker .
All newly sequenced mitochondrial genomic arrangements share the ancestral type for insects without rearrangement, but possess large non-coding regions (except the control region) in some lineages. The intergenic spacers in the mitogenomes of E. tibialis, E. gorhami and L. caraganae are similar to those in E. chinensis mitogenome . The total length of M. aulica’s intergenic spacers is not significantly different from the former four meloids, but its genome does contain more intergenic spacers (Table 2). Unexpectedly, the whole lengths of intergenic spacers in the three Hycleus mitogenomes are much longer than those of other meloids. The most remarkable feature is the presence of two IGSs in the mitogenomes of three Hycleus species. A 494-bp long intergenic spacer was also reported in the recently published mitogenome of H. chodschenticus . Consequently, the total lengths of known Hycleus mitogenomes are longer than those of other meloid mitogenomes, but predominantly due to the presence of IGSs rather than the lengths of genes or control regions.
The mitochondrial genome typically displays exceptional economy of organization, evidenced by lack of introns, few intergenic spacers, incomplete stop codons and even overlapping genes . However, the large IGSs in mitochondrial genomes were observed in some Hymenopteran [11, 13], Hemipteran [15, 17] Dictyopteran  and Coleopteran [5, 30] insects. Previously reported IGSs contain tandem repeat units (in Pyrocoelia rufa and Evania appendigaster) [5, 11], or additional origin of replication  and similar to it . The IGS in Reduviidae bugs have unidentified open reading frames encoding 103 or 104 amino acids but without blast similarity [15, 17]. In contrast, the two IGSs in Hycleus species have no significant similarity with other genes of their mitogenomes and lack open reading frames or tandem repeats. Nevertheless, the IGS2 of all studied Hycleus species contain discontinuous repeats of a 5 bp consensus motif (TACTA) (Fig. 5b). This motif was also found in many other Coleopteran insects , similar to the 7 bp motif ‘ATACTAA’ conserved in Lepidoptera  and the hexanucleotide motif ‘THACWW’ in Hymenoptera . The pentanucleotide motif was suggested as the possible binding site of a transcription termination peptide (mtTERM), as its position signifies the end of the J strand coding region in the circular mitochondrial DNA . However, we do not know the function of this discontinuous repeat.
The discontinuous repeats and the 18 bp long consensus sequence in both ends of IGS2 within all studied Hycleus species (Fig. 6a) suggest that the slipped-strand mispairing  may be the evolutionary mechanism of this IGS. According to this theory, mispairing involves dissociation of replicating DNA strands and then misaligned reassociation (Fig. 6b), following replication or repair lead to insertions of several repeat units. Formed tandem repeat experiences random loss and/or point mutation, only the repeat units in both ends are completely retained and the residues form the IGS2 (Fig. 6c). Although we are not absolutely certain about our assumption due to the highly divergent region, the fragmented repeat units in highly divergent region (Fig. 5b) and complete repeat units at both ends (Fig. 6a) suggest that the slipped-strand mispairing is the most convincing mechanism for IGS2.
The IGS between trnW and trnC was only previously reported in Trachypachus holmbergi , while no IGS (>30 bp) at this position has yet been found in the mitogenomes for Meloidae. Although H. chodschenticus has 3-bp intergenic spacer at the same position , it is too usual as many short intergenic spacers to be considered as IGS1. Therefore, the evolutionary mechanism of IGS1 in the three Hycleus species in present study is different from H. chodschenticus. The 9-bp consistent motif and the relatively high similarity among the three Hycleus species in the present study (Fig. 5) suggest that they have the mutual mechanism of IGS1. The duplication/random loss model [1, 54] may account for the IGS1. We speculate that tandem duplication of trnW-trnC-trnY is caused by some error in DNA replication, followed by random loss of partial duplicated genes, and then the residues form the IGS1 (Fig. 8). This model was also proposed as the mechanism for rearrangements in some Hemipteran  and Mantodea mitogenomes , and the IGSs in Blaptica dubia . It is possible that the duplication/random loss event of Hycleus beetles occurred relatively early and many nucleotides were deleted during the random loss period. Consequently, the residual IGS1 has low similarity with the original tRNAs, and the IGS1 is somewhat not conserved in Hycleus.
Similar to the three Hycleus species in the present study, H. chodschenticus also possesses the IGS between trnS2 and nad1 . Although intergenic spacers with approximately 20 bp at this position are common in Coleopteran mitogenomes, no large intergenic spacers have been found within other Coleopteran lineages to date. Based on the current knowledge of IGS2 across Coleoptera, it is reasonable to assume that the emergence of IGS2 in Hycleus mitogenomes might be a special evolutionary mutation. It is difficult to identify the closely related genera Hycleus and Mylabris since they share very similar morphological characteristics. However, IGS2 might be a potential marker to distinguish Hycleus from its closely related and indistinguishable genera with respect to the sizeable intergenic spacer that exists in all studied Hycleus species but is absent in other genera,. Although we are unable to adequately confirm that IGS2 exists in all Hycleus species due to limited samples, it provides a new candidate for molecular identification of this genus. Variations in the quantity, location and sequence of intergenic spacers might also be a valuable marker for phylogenetic and evolutionary studies at lower taxonomic levels, if these details of more taxa were obtained in future.
Phylogenetic studies indicated that different datasets and inference methods influence the tree topology of Coleoptera [31, 55]. Our phylogenetic results showed that tree topologies are sensitive to datasets rather than inference methods, since the different inference methods with the same datasets generated consistent tree topologies. The heterogeneous-site model in PhyloBayes was suggested as being more reliable for phylogenetic inferences within Coleoptera [31, 55]. Our Bayesian analyses under the heterogeneous-site model are unable to resolve phylogenetic relationships within Tenebrionoidea, but perhaps this is due to insufficient taxa. Although the nucleotide dataset of PCGs was proposed as better for phylogenetic analyses at superfamily level , we could not assess the quality of different datasets because some key lineages share the same tree topologies. However, nodal support values are sensitive to inference methods. For the same datasets, ML trees had significantly lower support values than BI trees, especially at several nodes that involved the family-level relationships of Aderidae, Ciidae, Scraptiidae (Fig. 7). This is consistent with previous phylogenetic studies using mitogenomes [14, 31, 37].
Phylogenetic relationships within Tenebrionoidea are ambiguous. The inter-family relationships are also uncertain, especially for Aderidae, Ciidae, Oedemerdae, and Prostomedae, which are respectively represented by only one taxon. However, all tree topologies well recover the monophyly of Tenebrionidae, Meloidae, and Mordellidae (Fig. 7, Additional file 10: Figure S1). Tenebrioninae and Diaperinae are never recovered as monophyletic groups, as is reported by Gunter et al. . To date, no phylogeny has successfully resolved the interfamilial relationships between tenebrionoids, either by using morphological or molecular characteristics. The comprehensive phylogeny of Coleoptera [57, 58] and the largest molecular phylogeny of Tenebrionoidea  were unable to recover strong support or definitively infer the phylogenetic relationships within the superfamilies. In contrast to these phylogenetic studies based on several genes, our phylogenetic inferences may be hindered by lack of across taxon sampling rather than dataset validity.
The genus Lytta is the sister group to Epicauta rather than grouped with Mylabris within Meloidae (Fig. 7). The placement of Lytta and whether Lytta is more closely related to Mylabris or Epicauta could not be inferred from previous molecular phylogeny of the family Meloidae based on 16S rRNA and ITS2 . We confirmed this relationship by multiple inference approaches based on both nucleotide and amino acid sequences of 13 mitochondrial PCGs. Considering the high diversity and rapid radiations of insects [59, 60], mitochondrial genomes could be a better approach to resolve intractable phylogenetic relationships due to its relatively rapid rate of mutation and purely maternal inheritance [3, 61]. Consequently, we believe our conclusion that the genus Lytta is more closely related to Epicauta than Mylabris or Hycleus is reliable because it is based on complete mitogenomes and the preference of stop codon usage.
The data collected thus far regarding meloid mitogenomes could not resolve the phylogenetic relationships within Meloidae. In fact, no phylogeny of Meloidae based on either morphological or molecular characteristics has been able to successfully resolve the relationships at genera and species levels. Taxon sampling is known to be one of the most significant determinants of accurate phylogenetic inferences, particularly in species rich lineages [62, 63]. Considering the diversity of the family Meloidae and the limitation of the present molecular information, more conclusive phylogenetic results will be achievable as bio-information becomes increasingly available. This study will assist with these more conclusive phylogenetic results and future studies on taxonomy, phylogeny and systematics of Meloidae insects.
Our study presents the mitochondrial genomes of seven meloid beetles. All complete mitogenomes of meloids retain the typical gene content and organization of the ancestor. The evolutionary rates of all PCGs in the studied Meloidae indicate that their evolution is according to purifying selection. The comparison of tRNA secondary structures exhibit diverse substitution patterns in Meloidae. Two large intergenic spacers exist in the three studied Hycleus mitogenomes, and the sequence and structure of the two IGSs contributed to our conclusion regarding their possible evolutionary mechanisms. The phylogenetic results inferred from mitochondrial genomes support that the genus Lytta is more closely related to Epicauta than to Mylabris. Although data collected thus far could not resolve the phylogenetic relationships within Meloidae, this study will assist in future mapping of the Meloidae phylogeny.
Samples collection and DNA extraction
The specimens of the seven meloid species used for this study are listed in Table 1, with locality data. The fresh materials were immediately preserved in 100% ethanol and stored in a − 20 °C refrigerator. Total genomic DNA was extracted from a frozen adult using Tianamp Genomic DNA kit.
PCR amplification and sequencing
All mitochondrial genomes of these collected meloid species were generated by amplifications of nine overlapping PCR fragments. Eight fragments were amplified using common primers for all seven species designed from the aligned E. chinensis mitogenome (GenBank accession number KP692789) , only the primers of the control region were specifically designed for each species. Details of primers are given in Additional file 11: Table S10. The PCR was performed with Vazyme Taq DNA Polymerase (Mg2+ Plus Buffer) and carried out on a PTC-100 thermal cycler (BioRad, Hercules, CA). PCR conditions used were: 5 min denaturation at 94 °C; 35 cycles of 30s at 94 °C, 30s at 49–56 °C and 1–3 min (1 min/Kbp) at 72 °C; followed by 10 min extension hold at 72 °C. The PCR products were sequenced on an ABI PRISM 3730 DNA sequencer by Tsingke Biotechnology Company with primers walking on both strands.
Genome assembly and annotation
Sequences from overlapping fragments were assembled with the neighboring fragments using SeqMan program included in the Lasergene software package (DNAStar Inc., Madison, Wisc.). Protein-coding regions were identified by ORF Finder from the NCBI website with invertebrate mitochondrial genetic codes, and compared with published mitochondrial sequences by using MEGA 6.0 . Most of the tRNA genes were identified using tRNAscan-SE 1.21  with invertebrate genetic codon predictors; however, the trnS1 was predicted by alignments with other homologous genes because of its lack of dihydrouridine (DHU) arm. The rRNA gene boundaries were interpreted as the end of a bounding tRNA gene and alignment of sequences with homologous regions of known Coleopteran mitogenomes. The control regions were assumed between rrnS and trnI within all meloid mitochondrial genomes. The A + T contents, relative synonymous codon usage values, and evolutionary rates (number of variable sites, nucleotide diversity, and Ka/Ks ratios) for each PCG were calculated via MEGA 6.0 .
Phylogenetic analyses were assessed using 29 Tenebrionoidea species representing 8 families, with the Chrysomelid Diabrotica barberi (GenBank accession number NC_022935.1)  employed as the outgroup. Species’ PCGs were extracted according to GenBank annotations by using GenScalpel . All these nucleotide and amino acid sequences were aligned using MUSCLE  with the default setting. Gaps and ambiguous sites were removed from the protein alignment to generate a 10,356-bp nucleotide dataset and a corresponding amino acid dataset (3452 amino acids). The best partitioning schemes and corresponding evolutionary models were selected by PartitionFinder 1.1.1  with 13 partitions defined by genes. We set the model selection as Bayesian information criterion (BIC), unlinked branch lengths, and greedy search algorithm to estimate the best fitting schemes. The best-fit partitioning schemes and corresponding models are shown in Additional file 12: Table S11.
Phylogenetic analyses with site-homogeneous model were conducted to reconstruct the phylogenetic trees of superfamily Tenebrionoidea by using Bayesian inference (BI) and maximum likelihood (ML) methods. Bayesian phylogenetic analysis was implemented using MrBayes 3.2.2 , and ran for 10,000,000 generations sampling per 1000 generations. Bayesian posterior probabilities were estimated using the Markov chain Monte Carlo (MCMC) sampling approach. ML analysis was carried out with RAxML-HPC2 on XSEDE 8.0.24  using 1000 bootstraps to estimate the node support. Bayesian analyses with a site-heterogeneous model were performed using PhyloBayes MPI 1.5a  with two MCMC chains run after the removal of constant sites under the CAT-GTR model.
- atp8 and apt6 :
ATPase subunits 8 and 6
- Cob :
- cox1-cox3 :
Cytochrome C oxidase subunit I-III
- nad1-nad6 and nad4L :
NADH dehydrogenase subunits 1–6 and 4 L
- rrnL :
large (16S) rRNA subunit
- rrnS :
small (12S) rRNA subunit
Relatively synonymous codon usages
Moritz C, Dowling TE, Brown WM. Evolution of animal mitochondrial DNA: relevance for population biology and systematics. Annu Rev Ecol Syst. 1987;18:269–92.
Avise JC. Molecular markers, natural history and evolution. New York: Champman & Hall; 1994.
Boore JL. Survey and summary: animal mitochondrial genomes. Nucleic Acids Res. 1999;27:1767–80.
Foottit RG, Adler PH. Insect biodiversity, science and society. Oxford: Wiley-Blackwell; 2009.
Bae JS, Kim I, Sohn HD, Jin BR. The mitochondrial genome of the firefly, Pyrocoelia rufa: complete DNA sequence, genome organization, and phylogenetic analysis with other insects. Mol Phylogenet Evol. 2004;32:978–85.
Chai HN, Du YZ, Zhai BP. Characterization of the complete mitochondrial genomes of Cnaphalocrocis medinalis and Chilo suppressalis (Lepidoptera: Pyralidae). Int J Biol Sci. 2012;8:561–79.
Timmermans MJ, Vogler AP. Phylogenetically informative rearrangements in mitochondrial genomes of Coleoptera, and monophyly of aquatic elateriform beetles (Dryopoidea). Mol Phylogen Evol. 2012;63:299–304.
Coates BS. Assembly and annotation of full mitochondrial genomes for the corn rootworm species, Diabrotica virgifera and Diabrotica barberi (Insecta: Coleoptera: Chrysomelidae), using next generation sequence data. Gene. 2014;542:190–7.
Ojala D, Montoya J, Attardi G. tRNA punctuation model of RNA processing in human mitochondria. Nature. 1981;290:470–4.
Crozier RH, Crozier YC. The mitochondrial genome of the honeybee Apis mellifera: complete sequence and genome organization. Genetics. 1993;133:97–117.
Wei SJ, Tang P, Zheng LH, Shi M, Chen XX. The complete mitochondrial genome of Evania appendigaster (hymenoptera: Evaniidae) has low a+T content and a long intergenic spacer between atp8 and atp6. Mol Biol Rep. 2010;37:1931–42.
Babbucci M, Basso A, Scupola A, Patarnello T, Negrisolo E. Is it an ant or a butterfly? Convergent evolution in the mitochondrial gene order of hymenoptera and Lepidoptera. Genome Biol Evol. 2014;6:3326–43.
Rodovalho CM, Lyra ML, Ferro M, Jr MB. The mitochondrial genome of the leaf-cutter ant Atta laevigata: a mitogenome with a large number of intergenic spacers. PLoS One. 2014;9:e97117.
Wei SJ, Li Q, van Achterberg K, Chen XX. Two mitochondrial genomes from the families Bethylidae and Mutillidae: independent rearrangement of protein-coding genes and higher-level phylogeny of the hymenoptera. Mol Phylogenet Evol. 2014;77:1–10.
Dotson EM, Beard CB. Sequence and organization of the mitochondrial genome of the Chagas disease vector, Triatoma dimidiata. Insect Mol Biol. 2001;10:205–15.
Shao R, Campbell NJ, Schmidt ER, Barker SC. Increased rate of gene rearrangement in the mitochondrial genomes of three orders of Hemipteroid insects. Mol Biol Evol. 2001;18:1828–32.
Hua J, Li M, Dong P, Cui Y, Xie Q, Bu W. Comparative and phylogenomic studies on the mitochondrial genomes of Pentatomomorpha (Insecta: Hemiptera: Heteroptera). BMC Genomics. 2008;9:1–15.
Cheng XF, Zhang LP, Yu DN, Storey KB, Zhang JY. The complete mitochondrial genomes of four cockroaches (Insecta: Blattodea) and phylogenetic analyses within cockroaches. Gene. 2016;586:115–22.
Ye F, Lan XE, Zhu WB, You P. Mitochondrial genomes of praying mantises (Dictyoptera, Mantodea): rearrangement, duplication, and reassignment of tRNA genes. Sci Rep. 2016;6 doi:10.1038/srep25634.
Beckenbach AT. Mitochondrial genome sequences of Nematocera (lower Diptera): evidence of rearrangement following a complete genome duplication in a winter crane fly. Genome Biol Evol. 2012;4:89–101.
Liu C, Jia C, Ma C, Li L, Zhou S. Mitochondrial genomes of two sinochlora species (Orthoptera): novel genome rearrangements and recognition sequence of replication origin. BMC Genomics. 2013;14 303–303
Shao R, Barker SC. The highly rearranged mitochondrial genome of the plague thrips, Thrips imaginis (Insecta: Thysanoptera): convergence of two novel gene boundaries and an extraordinary arrangement of rRNA genes. Mol Biol Evol. 2003;20:362–70.
Li H, Shao R, Song F, Zhou X, Yang Q, Li Z, Cai W. Mitochondrial genomes of two Barklice, Psococerastis albimaculata and Longivalvus hyalospilus (Psocoptera: Psocomorpha): contrasting rates in mitochondrial gene rearrangement between major lineages of Psocodea. PLoS One. 2013;8:e61685.
Cameron SL, Johnson KP, Whiting MF. The mitochondrial genome of the screamer louse Bothriometopus (Phthiraptera: Ischnocera): effects of extensive gene rearrangements on the evolution of the genome. J Mol Evol. 2007;65:589–604.
Crowson RA. The biology of Coleoptera. London: Academic Press; 1981.
Boore JL, Lavrov DV, Brown WM. Gene translocation links insects and crustaceans. Nature. 1998;392:667–8.
Timmermans MJ, Barton C, Haran J, Ahrens D, Culverwell CL, Ollikainen A, Dodsworth S, Foster PG, Bocak L, Vogler AP. Family-level sampling of mitochondrial genomes in Coleoptera: compositional heterogeneity and phylogenetics. Genome Biol Evol. 8:161–75.
Wei SJ, Chen XX. Progress in research on the comparative mitogenomics of insects. Chin J Appl Entomol. 2011;2015(48):1573–85.
Andujar C, Arribas P, Linard B, Kundrata R, Bocak L, Vogler AP. The mitochondrial genome of Iberobaenia (Coleoptera: Iberobaeniidae): first rearrangement of protein-coding genes in the beetles. Mitochondr DNA. 2016. doi:10.3109/19401736.2015.1115488.
Sheffield NC, Song H, Cameron SL, Whiting MF. A comparative analysis of mitochondrial genomes in Coleoptera (Arthropoda: Insecta) and genome descriptions of six new beetles. Mol Biol Evol. 2008;25:2499–509.
Yuan ML, Zhang QL, Zhang L, Guo ZL, Liu YJ, Shen YY, Shao R. High-level phylogeny of the coleoptera inferred with mitochondrial genome sequences. Mol Phylogenet Evol. 2016;104:99–111.
Bologna MA, Oliverio M, Pitzalis M, Mariottini P. Phylogeny and evolutionary history of the blister beetles (Coleoptera, Meloidae). Mol Phylogenet Evol. 2008;48:679–93.
Wang GS. Medical uses of mylabris in ancient China and recent studies. J Ethnopharmacol. 1989;26:147–62.
Rauh R, Kahl S, Boechzelt H, Bauer R, Kaina B, Efferth T. Molecular biology of cantharidin in cancer cells. Chin Med. 2007;2:8.
Bologna MA. Coleoptera Meloidae. Calderini, Bologna: Fauna d’Italia. XXVIII; 1991.
Ghoneim KS. Enhancement of research interests in physiology and biochemistry of blister beetles (Coleoptera: Meloidae): a review. Int Res J Biochem Bioinform. 2013;3:75–92.
Du C, He S, Song X, Liao Q, Zhang X, Yue B. The complete mitochondrial genome of Epicauta chinensis (Coleoptera: Meloidae) and phylogenetic analysis among coleopteran insects. Gene. 2016;578:274–80.
Taanman JW. The mitochondrial genome: structure, transcription, translation and replication. BBA-Bioenergetics. 1999;1410:103–23.
Fenn JD, Cameron SL, Whiting MF. The complete mitochondrial genome sequence of the Mormon cricket (Anabrus Simplex: Tettigoniidae: Orthoptera) and an analysis of control region variability. Insect Mol Biol. 2007;16:239–52.
Nardi F, Carapelli A, Fanciulli PP, Dallai R, Frati F. The complete mitochondrial DNA sequence of the basal hexapod Tetrodontophora bielanensis: evidence for heteroplasmy and tRNA translocations. Mol Biol Evol. 2001;18:1293–304.
Hua J, Li M, Dong P, Cui Y, Xie Q, Bu W. Comparative and phylogenomic studies on the mitochondrial genomes of Pentatomomorpha (Insecta: Hemiptera: Heteroptera). BMC Genomics. 2008;9:610.
Zhang K, Zhu W, Rong X, Liu J, Ding X, Hong X. The complete mitochondrial genome sequence of Sogatella furcifera (Horváth) and a comparative mitogenomic analysis of three predominant rice planthoppers. Gene. 2014;533:100–9.
Salvato P, Simonato M, Battisti A, Negrisolo E. The complete mitochondrial genome of the bag-shelter moth Ochrogaster lunifer (Lepidoptera: Notodonidae). BMC Genomics. 2008;9:331.
Negrisolo E, Babbucci M, Patarnello T. The mitochondrial genome of the ascalaphid owlfly libelloides macaronius and comparative evolutionary mitochondriomics of neuropterid insects. Bmc Genomics. 2011;12:221.
Hanada T, Suzuki T, Yokogawa T, Takemoto-Hori C, Sprinzl M, Watanabe K. Translation ability of mitochondrial tRNAs-Ser with unusual secondary structures in an in vitro translation system of bovine mitochondria. Genes Cells. 2001;6:1019–30.
Stewart JB, Beckenbach AT. Phylogenetic and genomic analysis of the complete mitochondrial DNA sequence of the spotted asparagus beetle Crioceris duodecimpunctata. Mol Phylogenet Evol. 2003;26:513–26.
Friedrich M, Muqim N. Sequence and phylogenetic analysis of the complete mitochondrial genome of the flour beetle Tribolium castanaeum. Mol Phylogenet Evol. 2003;26:502–12.
Boore JL. Complete mitochondrial genome sequence of the polychaete annelid Platynereis dumerilii. Mol Biol Evol. 2001;18:1413–6.
Boore JL, Brown WM. Mitochondrial genomes of Galathealinum, Helobdella, and Platynereis: sequence and gene arrangement comparisons indicate that Pogonophora is not a phylum and Annelida and Arthropoda are not sister taxa. Mol Biol Evol. 2000;17:87–106.
Andrews RM, Kubacka I, Chinnery PF, Lightowlers RN, Turnbull DM, Howell N. Reanalysis and revision of the Cambridge reference sequence for human mitochondrial DNA. Nat Genet. 1999;23:147.
Yukuhiro K, Sezutsu H, Itoh M, Shimizu K, Banno Y. Significant levels of sequence divergence and gene rearrangements have occurred between the mitochondrial genomes of the wild mulberry silkmoth, Bombyx mandarina, and its close relative, the domesticated silkmoth, Bombyx mori. Mol Biol Evol. 2002;19:1385–9.
Zhang DX, Hewitt GM. Insect mitochondrial control region: a review of its structure, evolution and usefulness in evolutionary studies. Biochem Syst Ecol. 1997;25:99–120.
Levinson G, Gutman GA. Slipped-strand mispairing: a major mechanism for DNA sequence evolution. Mol Biol Evol. 1987;4:203–21.
Boore JL. The duplication/random loss model for gene rearrangement exemplified by mitochondrial genomes of Deuterostome animals. In: Sankoff D, Nadeau J, editors. Comparative genomics. Dordrecht: Kluwer Academic Publishers; 2000. p. 133–47.
Timmermans MJ, Barton C, Haran J, Ahrens D, Culverwell CL, Ollikainen A, Dodsworth S, Foster PG, Bocak L, Vogler AP. Family-level sampling of mitochondrial genomes in Coleoptera: compositional heterogeneity and phylogenetics. Genome Biol Evol. 2015;8:161–75.
Gunter NL, Levkaničová Z, Weir TH, Ślipiński A, Cameron SL, Bocak L. Towards a phylogeny of the Tenebrionoidea (Coleoptera). Mol Phylogenet Evol. 2014;79:305–12.
Hunt T, Bergsten J, Levkanicova Z, Papadopoulou A, John OS, Wild R, Hammond PM, Ahrens D, Balke M, Cateriho MS, Gómez-Zurita J, Ribera I, Barraclough TG, Bocakova M, Bocak L, Vogler AP. A comprehensive phylogeny of beetles reveals the evolutionary origins of a superradiation. Science. 2007;318:1913–6.
Bocak L, Barton C, Crampton-Platt A, Chesters D, Ahrens D, Vogler AP. Building the Coleoptera tree-of-life for >8000 species: composition of public DNA data and fit with Linnaean classification. Syst Entomol. 2014;39:97–110.
Caterino MS, Shull VL, Hammond PM, Vogler AP. Basal relationships of Coleoptera inferred from 18S rDNA sequences. Zool Scr. 2002;31:41–9.
Whitfield JB, Kjer KM. Ancient rapid radiations of insects: challenges for phylogenetic analysis. Annu Rev Entomol. 2008;53:449–72.
Zardoya R, Meyer A. Phylogenetic performance of mitochondrial protein-coding genes in resolving relationships among vertebrates. Mol Biol Evol. 1996;13:933–42.
Zwickl DJ, Hillis DM. Increased taxon sampling greatly reduces phylogenetic error. Syst Biol. 2002;51:588–98.
Rosenberg NA. Statistical tests for taxonomic distinctiveness from observations of monophyly. Evolution. 2007;61:317–23.
Tamura K, Stecher G, Peterson D, Filipski A, Kumar S. Mega6: Molecular evolutionary genetics analysis version 6.0. Mol Biol Evol. 2013;30:2725–9.
Lowe TM, Eddy SR. tRNAscan-SE: a program for improved detection of transfer RNA genes in genomic sequence. Nucleic Acids Res. 1997;25:955–64.
Yin YH, Du LM, Yue BS. GenScalpel: an application for sequence retrieval and extraction from the GenBank flatfile. J Hered. 2012;103:908–11.
Edgar RC. MUSCLE: multiple sequence alignment with high accuracy and high throughput. Nucleic Acids Res. 2004;32:1792–7.
Lanfear R, Calcott B, Ho SY, Guindon S. Partitionfinder: combined selection of partitioning schemes and substitution models for phylogenetic analyses. Mol Biol Evol. 2012;29:1695–701.
Ronquist F, Teslenko M, van der Mark P, Ayres DL, Darling A, Höhna S, Larget B, Liu L, Suchard MA, Huelsenbeck JP. MrBayes 3.2: Efficient Bayesian phylogenetic inference and model choice across a large model space. Syst Biol. 2012;61:539–42.
Stamatakis A. RAxML version 8: a tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics. 2014;30:1312–3.
Miller MA, Pfeiffer W, Schwartz T. Creating the CIPRES Science Gateway for inference of large phylogenetic trees. Gateway Computing Environments Workshop (GCE). Gateway Computing Environments Workshop (GCE). 2010;14:1-8.
We sincerely appreciate Zhao Pan (Hebei University) for helping to identify some samples and Megan Price (Sichuan University) for English language editing. The research was funded by the Science and technology project of Sichuan Province, 2015SZ0036 and the state standardization project of traditional Chinese medicine, ZYBZH-C-SC-53.
Science and technology project of Sichuan Province, 2015SZ0036.
State standardization project of traditional Chinese medicine, ZYBZH-C-SC-53.
Availability of data and materials
All the data supporting our findings is contained within the manuscript.
Ethics approval and consent to participate
Consent for publication
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Annotation of the Mylabris aulica mitogenome. (DOCX 21 kb)
Annotation of the Hycleus phaleratus mitogenome. (DOCX 21 kb)
Annotation of the Hycleus marcipoli mitogenome. (DOCX 21 kb)
Annotation of the Hycleus cichorii mitogenome. (DOCX 22 kb)
Annotation of the Epicauta gorhami mitogenome. (DOCX 21 kb)
Annotation of the Epicauta tibialis mitogenome. (DOCX 21 kb)
Annotation of the Lytta caraganae mitogenome. (DOCX 21 kb)
Codon usage in mitochondrial genomes of seven meloids. (DOCX 27 kb)
Summary of multiple alignments of tRNAs in meloid mitogenomes. (DOCX 15 kb)
The phylogenetic tree of 16 species from superfamily Tenebrionoidea based on the nucleotide dataset and the amino acids dataset of 13 mitochondrial protein-coding genes, inferred from PhyloBayes. The numbers abutting branches refer to Bayesian posterior probabilities. The Diabrotica barberi (Coleoptera: Chrysomelidae) was employed to root the trees as outgroup. (TIFF 1285 kb)
PCR primers used to amplify the mitochondrial genomes of seven meloids. (DOCX 19 kb)
The best-fit schemes and evolutionary models for two datasets. (DOCX 16 kb)
About this article
Cite this article
Du, C., Zhang, L., Lu, T. et al. Mitochondrial genomes of blister beetles (Coleoptera, Meloidae) and two large intergenic spacers in Hycleus genera. BMC Genomics 18, 698 (2017). https://doi.org/10.1186/s12864-017-4102-y