Mitochondrial genomes of two Sinochloraspecies (Orthoptera): novel genome rearrangements and recognition sequence of replication origin

Genome organization

We sequenced the complete mitogenome of S. longifissa (18,133 bp) and the nearly complete mitogenome of S. retrolateralis (17,209 bp) with a partial CR. The mitogenomes of S. longifissa and S. retrolateralis have been deposited in the GenBank database under accession numbers of KC467055 and KC467056, respectively. They are currently the largest Orthoptera mitogenomes on GenBank. Their large sizes are due to two large non-coding regions, i.e., the CR and one intergenic spacer (IGS) located between trnS(UCN) and nad1. Both mitogenomes (Table 1) resemble the available ensiferan mitogenomes [4] and the proposed ancestral insect mitogenome [1] in gene content and orientation, tRNA anticodons, and tRNA/rRNA structures (Figure 1). Additionally, there are several gene overlaps, such as the open reading frames (ORF) of atp6-atp8 and nad4L-nad4, each of which overlaps by seven nucleotides (Table 1). A few other IGS are also present in the mitogenomes (Table 1).

Nucleotide composition

Like other ensiferans, the two mitogenomes show comparatively low A+T bias (69.05% in S. longifissa and 70.08% in S. retrolateralis) (Additional file 1), compared with most caeliferans [4]. The AT skew is 0.0017 in S. longifissa and 0.0063 in S. retrolateralis, while the GC skew is -0.2983 in the former and -0.3218 in the latter, indicating weak A skew and strong C skew in the J-strand. Through comparison of nucleotide skew in all sequenced Orthoptera, the skew divergence between the two suborders was detected. Within the Ensifera, the AT skew values of the whole J-strand range from -0.1 to 0.1 and GC skew values are lower than -0.25, showing the same skew patterns as in Sinochlora species (Figure 2). By contrast, within the Caelifera except a tridactylid Ellipes minuta[4], absolute values of AT skew and GC skew of the J-strand range from 0.1 to 0.25, indicating pronounced A skew and C skew (Figure 2). The synonymous fourfold degenerate third codon positions (P_4fd) suffer less selective constraints, and thus could indicate background mutational pressures on nucleotide skew [15]. The skew patterns at the P_4fd between the two suborders are consistent with those on the J-strand (Figure 2). It is proposed that nucleotide composition bias and skew are caused by the selection-mutation-drift equilibrium in the molecular evolution [16]. The asystematical directional mutation pressure and corresponding deaminations of cytosine and adenine in the mitogenomes are involved in the mutation processes [17].

Start and stop codons

All PCGs except cox1 and nad6 start with typical ATN codons. A previous study concerning Orthoptera mitogenomes suggested that the cox1 gene may start with an irregular tetranucleotide codon AUGA [4]. Similar irregular tetranucleotide codons were also proposed as start codons in Drosophila[18, 19]. However, there is no experimental evidence for the use of a 4-bp start codon in any creature. Recent research on characteristics of mature mRNA and rRNA genes from D. melanogaster mitochondria showed that UCG serves as the start codon of its cox1 gene [20]. The result is consistent with predictions of the start of cox1 suggested by comparison of conserved amino acid positions [21]. By comparing amino acid sequences of cox1 in all sequenced Orthoptera, we observed that the first conserved codon is the TCN serine codon, downstream of which at least three codon positions are also conserved and there are no standard ATN bases (Additional file 2). Thus, in both Sinochlora species, the conserved TCT serine codon may also serve as a start codon. Concerning the nad6 gene, it has been proposed to start with ATN codons in other sequenced orthopterans; however, there are no conserved amino acid positions within the span of over 50 amino acids downstream of trnP (Additional file 3). Then we propose that the nad6 may start with CTG in S. longifissa and TTG in S. retrolateralis, considering previous designation of such start codons [22, 23]. The two start codons, creating no IGS or overlap between trnP and nad6 genes, also appear to be more plausible in the evolutionary economic perspective [22].

Two standard stop codons TAA/TAG and two incomplete stop codons T/TA are utilized in the PCGs. For S. longifissa, six PCGs (nad2, atp8, atp6, cox3, nad3, and nad4L) terminate with TAA, one (nad1) with TAG, one (nad6) with TA, and the other four (cox1, nad5, nad4, and cob) with T. The stop codons differ between the two congeners in that TAA is utilized as the cox2 stop codon, and TA is utilized as the atp6 and cox3 stop codon in S. retrolateralis. A partial T or TA stop codon, which was proposed to create complete TAA stop codons via posttranscriptional polyadenylation [24], is also present in other metazoan mitochondrial genes.

Novel gene order rearrangement involving the control region

One of our most significant findings is the novel gene order rearrangement "rrnS-trnI-trnM-nad2-CR-trnQ-trnW" in both Sinochlora species (Figure 1). The two mitogenomes are the first representatives that have a large-scale translocation involving the CR in Orthoptera. Such gene order has not been observed in other sequenced insect mitogenomes. The gene cluster rrnS-CR-trnI-trnQ-trnM-nad2 in the proposed ancestral mitogenome has been discovered to function as a hot spot for gene rearrangements in arthropod mitogenomes. Duplicate trnI and partial trnQ genes have been discovered in the CR close to the rrnS gene in some blowflies [25, 26]. In a mantid, a complicate set of repeat units dispersing in both ends of the CR also translocated between trnM and nad2[27]. A plague thrips displayed CR duplications that are distant from the rRNA genes [28]. Various locations of the CR and neighboring tRNA genes have been exhibited in lice [29]. Location variability of the putative CR has also been reported as a product of the tRNA gene translocation, such as rrnS-trnQ-CR-trnI-trnS(UCN)-trnM-nad2 in a springtail [30], or the tRNA and rRNA gene translocation, such as rrnL-trnV-trnS(UCN)-trnC-CR-rrnS in a mite [31]. The CR and the neighboring genes have also been found to be duplicated in some ticks [32].

Organization of the control region

For both Sinochlora species, the putative CR is composed of two major sections with different nucleotide compositions: a highly A+T-biased region (ATR) adjacent to nad2 and several GC-rich TR units adjacent to trnQ. For the former, the A+T content is 75.70% in S. longifissa and 75.61% in S. retrolateralis. By contrast, the TR region has low bias toward A+T (55.52% in S. longifissa and 58.26% in S. retrolateralis). The ATR is conserved between the two congeners in length, nucleotide composition and secondary structure. Its length is 1,391 bp in S. longifissa and 1,386 bp in S. retrolateralis, consistent with the A+T-rich region of other orthopterans [4, 33]. The two species share 82.6% sequence identity in this region. A long T-stretch and several conserved potential hairpin structures are thought to play significant roles in initiating and/or regulating the transcription and replication of mtDNA [8]. They were also found in the both Sinochlora species. The long T-stretch (> 18 bp) is situated in the initial quarter of the ATR on the J-strand, and a similar T-stretch (> 8 bp) is situated in the initial third of the ATR on the minority strand (N-strand) (Figure 3A). The two species also share 97.10% sequence identity in the domain located between these two T-stretches. The domain has seven similar sites capable of forming 8- to 33-bp stem-loop structures (Figure 3B). Additionally, there are three similar sites for the formation of 11- to 25-bp stable hairpins, with 7- to 11-bp loops supported by 1–5 GC matches in the stems, adjacent to the 5’ end of the T-stretch on the J-strand (Figure 3B).

In S. longifissa, the TR region is 1,731 bp in length, comprising nine full 175-bp copies and a partial 153-bp copy. In S. retrolateralis, a total of five 175-bp repeat copies have been successfully sequenced from both ends of the TR; however, we failed to sequence through the whole TR region due to the presence of a large number of repeat copies. We estimate that there are twelve 175-bp tandem copies in S. retrolateralis based on the length (~2,100 bp) of the TR indicated by gel electrophoresis. The consensus TR motifs of the two species share 61.2% sequence identity. In S. longifissa, four of the nine TR motifs are identical to the consensus TR motif, and the other five have a few mutations and/or deletions. Among the nine TR motifs, all but the fourth possess an ORF. The consensus ORF could be translated into 58 amino acids, starting with lysine (AAA) and ending with valine (GTG) (Additional file 4). Low shared identities (< 45%) with nuclear sequences in search of the GenBank database exclude the possibility that the ORF is transferred from the nuclear genome. It has 47.1% identity with the antisense strand of nad4, and 45.9% identity with the sense strand of cox1. The shared bases are scattered in alignment, among which over 60% are A and T (Additional file 5). This suggests that the ORF might not be obtained through horizontal transfer between mitochondria. In contrast, no ORF could be found in the S. retrolateralis TR motif.

The TR motif could be duplicated through slipped-strand mispairing [34, 35]. Moreover, potential stem-loop structures in a repeated unit and its flanking part have been demonstrated to cause an increase in slipped-strand mispairing frequency [35, 36]. Stem-loop structures were detected in the TR region of the two Sinochlora species. For example, in the first TR unit and its junctions at both ends in S. longifissa, the nucleotide sequence can potentially form seven 8- to 27-bp hairpins with 5- to 10-bp loops supported by 1–5 GC matches in the stem (Figure 3C). In the corresponding regions in S. retrolateralis, the nucleotide sequence can potentially form six 10- to 29-bp hairpins with 5- to 15-bp loops supported by 1–4 GC matches in the stem (Figure 3D). Similar complicated hairpin structures in TR were also detected in the louse Bothriometopus[23] and termites [37, 38]. However, there are no conserved stem-loops among the TRs and at the joints in these insects.

Unassigned intergenic spacers and trnS(UCN)pseudogenes

Eleven IGSs with identical locations are present in S. longifissa (totalling 318 bp) and S. retrolateralis (263 bp); one additional IGS is situated between atp6 and cox3 in S. longfissa (2 bp) and between cox2 and trnK in S. retrolateralis (1 bp) (Table 1). The largest IGS (161 bp for S. longifissa and 94 bp for S. retrolateralis) lies between trnS(UCN) and nad1. By comparing the IGS of currently available Orthoptera mitogenomes, we found a 7-bp conserved motif (THYTHDA) downstream the nad1 across Orthoptera, with the only exception of Mekongiella xizangensis (Additional file 6). The conserved motif has been proposed as a binding site for mitochondrial transcription termination factor (mtTERM) in Orthoptera [4]. Similar conserved motifs, which were proposed as a binding site of the mtTERM, between trnS(UCN) and nad1 have also been found in Lepidoptera [39], and Coleoptera [40]. It has been confirmed that one of the two binding sites of mtTERM lies downstream of nad1 in Drosophila[41, 42].

In the IGS region between trnS(UCN) and nad1, we detected the vestige of one additional trnS(UCN) gene copy in both Sinochlora species. While both trnS(UCN) pseudogenes possess the same anticodon to the functional trnS(UCN), they are unable to form stable cloverleaf secondary structures (Figure 4). Relative-rate tests show that both pseudogenes evolved more than twice as fast as their functional paralogs (Table 2). This suggests that the redundant trnS(UCN) pseudogenes have experienced more relaxed selective constraints. In addition, the IGS region in S. longifissa could be folded into three 19- to 25-bp stable stem-loop structures (Figure 4C). The anticodon lies in the stem of the second hairpin and the 7-bp conserved motif aforementioned lies adjacent to the 3’ end of the third hairpin (Figure 4C). In contrast, there is no hairpin in the IGS region in S. retrolateralis.

Mechanism of genome rearrangements

The predominant mechanism of mitogenome rearrangements is considered to be partial genome duplication followed by a random [34, 43, 44] or non-random loss of the duplicated gene copies [45]. For the rearrangement involving the CR neighbourhood, we propose that the ancestor of the two Sinochlora species possessed the same gene order as ancestral insects [1]. Firstly, the gene cluster CR-trnI-trnQ-trnM-nad2 was duplicated, likely promoted by the stem-loop structures detected in the CR (Figure 3C,D). Subsequently, non-random deletions might happen to the entire subset of the duplicated genes, dependent on the gene copy’s transcriptional polarity and location in the genome, because the CR includes a transcriptional control sequence [45]. Namely, the redundant genes which possess the same polarity, such as the redundant gene cluster trnI-trnM-nad2, were successively deleted during replication (Figure 5). By contrast, persistence of the trnS(UCN) pseudogenes may have resulted from the tandem duplication of trnS(UCN) followed by a subsequent random loss of the trnS(UCN) paralogs. Apparently, more knowledge of transcription signals in insect mitogenomes and sampling additional closely related taxa will be helpful to understand the rearrangement mechanism.

Recognition sequences of replication origins

The T-stretch is involved in the recognition of the O_R of mtDNA at least among holometabolous insects, whereas similar T-stretch was not found in the upstream portion of the O_R for L. migratoria (hemimetabolous insect) [8]. However, after a thorough re-investigation, we detected long T-stretches or variants of T-stretch in the putative CR of available Orthoptera mitogenomes (Additional file 7). The long T-stretch can be observed in the available sequences of all katydids (Figure 6), except Ruspolia dubia, and varies in size from 10 (E. cheni) to 21 bp (S. retrolateralis). The T-stretches could participate in the formation of or be positioned adjacent to a possible stem-loop structure. By contrast, in crickets and grasshoppers, T-stretch variants can be detected within a potential stem-loop structure (Figure 7). In the T-stretch variants, a few transitions between thymine and cytosine occurred to form (T)nC(T)n sequences. The T-stretch variants lie proximal to the middle part of the CR in the crickets, varying in size from 13 to 19 bp (Figure 7A). The T-stretch variants vary in size from 16 to 18 bp in grasshoppers (Figure 7B), with the exception of Alulatettix yunnanensis. In A. yunnanensis, a small T-stretch could be observed in the stem portion of a hairpin structure, similar to that observed in katydids (Figure 7B).

Human mtDNA synthesis is initiated from the sites near the stem base of the secondary structure located around the light strand origin, and the small T-stretch (6–11 bp) located in the loop portion participates in the initiation process [46]. The secondary structure in human and other vertebrate mtDNA is very similar to the aforementioned stem-loop structure in the orthopterans [8]. Therefore, the long T-stretch, which participates in the formation of or is adjacent to a possible stem-loop structure in katydids, and the T-stretch variants within the stem-loop structure of crickets and most grasshoppers, could also play a crucial role in mtDNA replication initiation. Notably, the T-stretch in cockroaches and termites also participates in the formation of a certain stem-loop structure (Additional file 8). However, the stem-loop structure around the mtDNA O_R could not be detected in Drosophila[8]. Therefore, the O_R recognition sequences of mtDNA, although generally detected in Orthoptera, have diverged not only among Orthoptera but also throughout insect evolution.

Taxon sampling and mitochondria DNA extraction

Specimens of S. longifissa and S. retrolateralis were collected from Wuyi Mountain, Fujian, South China in 2005. The specimens were preserved in 95% ethanol and stored at 4°C. The mitochondria were isolated as previously described [49], and mtDNA was extracted with the DNeasy Blood & Tissue Kit (Qiagen).

Genome determination

First, short gene regions within individual genes (cox1, cox3, cytb, nad1, rrnL, and rrnS) were amplified and sequenced using listed primers (Additional file 9). Then the obtained sequences were used to design specific primers for amplifying overlapping fragments spanning the whole mitogenomes.

Fragments larger than 3 kb were amplified using TaKaRa LA Taq™ (Takara, Dalian, China), with the following cycling conditions: an initial denaturation at 94°C for 3 min, followed by 36 cycles of denaturation at 94°C for 30 s, annealing at 50–57°C for 30 s, and extension at 68°C for 3–8 min (1 kb/min), with a final elongation at 68°C for 6 min after the last cycle. 16 fragments smaller than 3 kb were performed using TaKaRa ExTaq™ or TaKaRa rTaq™ (Takara, Dalian, China), with the following cycling conditions: an initial denaturation at 94°C for 3 min, followed by 36 cycles of denaturation at 94°C for 30 s, annealing at 45–60°C for 30 s, and extension at 72°C for 1–2 min (1 kb/min), with a final elongation at 72°C for 6 min after the last cycle. After purification with AxyPrep™ DNA Gel Extraction Kit, most PCR products were directly sequenced by means of primer walking, and other fragments were cloned into the pGEM-T Easy vector (Promega, USA) prior to sequencing.

Concerning the CR, long PCR amplicons were successfully amplified, which encompassed the entire CR from nad1 to cox2 genes for both species, but gel electrophoresis showed multiple bands. The largest and brightest band was chosen to clone into the pGEM-T Easy vector for sequencing. Sequencing primers were designed from flanking regions of the whole TR region and subsequently a 400-bp sequence at each end was obtained. The obtained TR units were analyzed with the TRF4.0 software [50] in order to design suitable primers for walking. For S. longifissa, the complete TR region was sequenced using specific primers that was designed based on a mutant poly C (8 continuous C) (Additional file 4) in one of the TR units. For S. retrolateralis, the length of the PCR product indicates that the TR region is about 2,100 bp, suggesting that there are 12 complete tandem repeats; however, only 5 of these tandem repeats were sequenced.

Sequence assembly, annotation and secondary structure prediction

The complete mitogenome sequences were assembled using the SeqMan program from the Lasergene package software (DNAStar, Madison, WI). tRNA genes were identified by their cloverleaf secondary structures using tRNAscan-SE 1.21 [51]. The locations of 13 PCGs and rRNA genes were determined by comparison of homologous sequences with other sequenced orthopterans using the CLUSTAL W programs [52]. Nucleotide composition statistics, nucleotide bias and skew of the orthopterans (Additional file 10) except those at P_4fd, were retrieved from the METAMiGA database [53]. Skewness was calculated to describe strand bias [54], which measures the relative number of As to Ts (AT skew = [A - T]/[A + T]) and Gs to Cs (GC skew = [G - C]/[G + C]). We obtained codon usage from the SMS2 website [55] and computed AT and GC skew at the P_4fd on the J-strand. Potential stem-loop structures of the polyneopterans (Additional file 10) were predicted by the Mfold software [56]. Repeat sequences were identified with the TRF4.0 software. ORF was detected using the MEGA5 software [57].

Phylogenetic analysis

Divergence and substitution rates between trnS(UCN) genes and pseudogenes were investigated, using trnS(UCN) of Myrmecophilus manni[33] as outgroup. Gapped positions were eliminated from the resulted alignment. The remaining 44 sites (including 26 parsimony informative) were used to reconstruct a distance phylogeny by Minimum Evolution using JC distances [58], due to the low number of sites analyzed. Relative-rate tests [59] were used to calculate substitution rates employing RRTree version 1.1.11 [60], presuming JC distances.