The complete mitochondrial genome of the house dust mite Dermatophagoides pteronyssinus (Trouessart): a novel gene arrangement among arthropods

Background The apparent scarcity of available sequence data has greatly impeded evolutionary studies in Acari (mites and ticks). This subclass encompasses over 48,000 species and forms the largest group within the Arachnida. Although mitochondrial genomes are widely utilised for phylogenetic and population genetic studies, only 20 mitochondrial genomes of Acari have been determined, of which only one belongs to the diverse order of the Sarcoptiformes. In this study, we describe the mitochondrial genome of the European house dust mite Dermatophagoides pteronyssinus, the most important member of this largely neglected group. Results The mitochondrial genome of D. pteronyssinus is a circular DNA molecule of 14,203 bp. It contains the complete set of 37 genes (13 protein coding genes, 2 rRNA genes and 22 tRNA genes), usually present in metazoan mitochondrial genomes. The mitochondrial gene order differs considerably from that of other Acari mitochondrial genomes. Compared to the mitochondrial genome of Limulus polyphemus, considered as the ancestral arthropod pattern, only 11 of the 38 gene boundaries are conserved. The majority strand has a 72.6% AT-content but a GC-skew of 0.194. This skew is the reverse of that normally observed for typical animal mitochondrial genomes. A microsatellite was detected in a large non-coding region (286 bp), which probably functions as the control region. Almost all tRNA genes lack a T-arm, provoking the formation of canonical cloverleaf tRNA-structures, and both rRNA genes are considerably reduced in size. Finally, the genomic sequence was used to perform a phylogenetic study. Both maximum likelihood and Bayesian inference analysis clustered D. pteronyssinus with Steganacarus magnus, forming a sistergroup of the Trombidiformes. Conclusion Although the mitochondrial genome of D. pteronyssinus shares different features with previously characterised Acari mitochondrial genomes, it is unique in many ways. Gene order is extremely rearranged and represents a new pattern within the Acari. Both tRNAs and rRNAs are truncated, corroborating the theory of the functional co-evolution of these molecules. Furthermore, the strong and reversed GC- and AT-skews suggest the inversion of the control region as an evolutionary event. Finally, phylogenetic analysis using concatenated mt gene sequences succeeded in recovering Acari relationships concordant with traditional views of phylogeny of Acari.


Background
Approximately 48,000 mite and tick species (Arthropoda: Chelicerata: Arachnida: Acari) [1,2] have been described to this day. Since the number of undescribed species is thought to be twenty-fold higher, this subclass is by far the most species-rich group among the Arachnida [2]. Acari diversified 400 million years ago and currently three major lineages are recognised: Opilioacariformes, Acariformes and Parasitiformes [2,3]. The Acariformes comprise two major groups, the Trombidiformes and Sarcoptiformes [2,4]. Two of the most prominent members of the Sarcoptiformes are the European house dust mite Dermatophagoides pteronyssinus (Trouessart, 1897) and the American house dust mite Dermatophagoides farinae (Hughes, 1961), both belonging to the family of the Pyroglyphidae (cohort Astigmata). Pyroglyphid mites are typical inhabitants of animal nests. In the human environment, they are mainly found in upholstery, textile floor covers and beddings, where they primarily feed on the skin scale fraction in house dust [5]. About 40 years ago, house dust mites were first recognized as one of the major sources of allergens in house dust [6]. The allergenic proteins are found in high concentrations in mite faeces, which, after drying and pulverizing, become airborne and can be inhaled. The presence of these allergens in sensitive persons is able to cause diseases like asthma, dermatitis and rhinitis [7,8]. In countries with a temperate climate, 6 to 35 per cent of the population is sensitive to house dust mite-derived allergens [9].
Complete mitochondrial (mt) genome sequences are becoming increasingly important for effective evolutionary and population studies. Mt genome sequences are not only more informative than shorter sequences of individual genes, they also provide sets of genome-level characters, such as the relative position of different genes, RNA secondary structures and modes of control of replication and transcription [10][11][12][13]. However, the applicability of mt genomes as a marker of highly divergent lineages is still controversial [14,15] and remains to be elucidated [16]. In addition, unravelling mt genomes can be of economic importance as well, since several chemical classes of pesticides target mt proteins. Well-known acaricides like acequinocyl and fluacrypyrim affect mt electron transport through the inhibition of the mt encoded cytochrome b in complex III [17]. Also, the economically important class of METI (Mitochondrial Electron Transfer Inhibitors)-acaricides target the mt complex I, although their exact molecular target has not yet been elucidated. Recently, resistance to the acaricide bifenazate was shown to be caused by mutations in the mt encoded cytochrome b and to have evolved rapidly through a short stage of mt heteroplasmy [18].
At present, the mt genomes of 20 species belonging to the Acari are available at NCBI ( [19], status January 10, 2009). Most of the submitted sequences have the typical features of metazoan mt genomes. They are circular, between 13 and 20 kb in length, contain a coding region with 37 genes (22 tRNAs, 2 rRNAs and 13 protein coding genes) and a relatively small non-coding region. The latter is mostly AT-rich and fulfils a role in the initiation of replication and transcription [20,21]. Compared to this typical configuration, the mt genomes of Steganacarus magnus, Metaseiulus occidentalis and Leptotrombidium pallidum show some abnormal features. S. magnus lacks 16 of the 22 tRNAs normally present in mt genomes [22]. M. occidentalis has a unusually large mt genome (24.9 kb) resulting from a duplication event of a large fragment of the codon region. Despite its large size, genes coding for nad6 and nad3 were not found during the initial annotation process [23]. L. pallidum on the other hand has 38 mt genes due to a duplication of the 16S-rRNA [24].
In this study, we analyse the complete mt genome of a member of the Sarcoptiformes, the European house dust mite D. pteronyssinus, after obtaining the complete sequence using a long PCR approach.

Genome organisation
The mt genome of D. pteronyssinus was amplified, using long PCR, in three overlapping fragments. The final assembled sequence was 14,203 bp [GenBank: EU884425; Fig. 1], making it the fifth smallest sequenced genome within the Acari. Only the mt genomes of Tetranychus urticae (13,103 bp), Leptotrombidium akamushi (13,698 bp), Leptotrombidium deliense (13,731 bp) and S. magnus (13,818 bp) are smaller (Table 1 [ [25][26][27]). As nonspecific amplification artefacts and incomplete coverage of genes are well-known drawbacks of a PCR approach [28], we checked the genome size by restriction digest on rolling circle amplified mtDNA (Fig. 2). This approach confirmed the sequence size, considering that the relative mobility of mtDNA restriction fragments can show slight (5-12%) deviations compared to their sequence length [29]. The mt genome of D. pteronyssinus is the first mt sequence of a mite belonging to the Astigmata and is together with the mt genome of S. magnus the only representative from the order of the Sarcoptiformes. Adding this genome to the database resulted in 21 publicly available Acari mtDNA sequences. Twelve belong to species in the superorder of the Parasitiformes whereas nineamong which D. pteronyssinus -belong to species in the superorder of the Acariformes.
The mt genome of the horseshoe crab Limulus polyphemus is considered to represent the ground pattern for arthropod mt genomes [30,31]. Comparing the D. pteronyssinus genome to this sequence revealed that only 11 of the 38 gene boundaries in L. polyphemus are conserved in D. pteronyssinus ( Fig. 3 [32]). Moreover, by making use of the pattern search function in the Mitome-database ( [33], status January 10, 2009), the mt gene order of D. pteronyssinus appeared to be unique among arthropods.
Remarkably, the relative position of trnL 2 (between nad1 and 16S-rRNA), which differentiates the Chelicerata, Myriapoda and Onychophora from the Insecta and Crustacea according to Boore [12,34], is not conserved. However, Boore's hypothesis was based on mt genome data from only 2 Chelicerata that were available in 1998. At present, 41 complete chelicerate mt genomes are available in the NCBI-database ( [19], status January 10, 2009). Out of these, only 29 depict the specific arrangement of trnL 2 between nad1 and 16S-rRNA (see additional file 1 for an overview of gene arrangements of chelicerate mt genomes). This illustrates that care should be taken when general rules are deduced from limited datasets.
Mt gene arrangements have already provided strong support toward the resolution of several long-standing controversial phylogenetic relationships [12]. Surprisingly, the mt gene order of D. pteronyssinus differs considerably from that of other mites (see additional file 1). Comparing the D. pteronyssinus mt genome to the mt sequence of the oribatid S. magnus [22], the closest relative of D. pteronyssinus, revealed that only 6 of the 22 gene boundaries in S. magnus are conserved in D. pteronyssinus (Fig. 3).
Extending this analysis to the other Acari mt genomes showed that in several cases the set of neighboured genes that were not separated during the evolution [35] was greater between members of different superorders (e.g. D. pteronyssinus (Acariformes) and Rhipicephalus sanguineus (Parasitiformes)) than between members of the same superorder (e.g. D. pteronyssinus (Acariformes) and T. urticae (Acariformes)) (Table2). Exclusion of tRNAs in our analysis showed a similar trend, suggesting that protein coding genes were also involved in mt gene rearrangements. These results indicate that mt gene orders seem less useful for deduction of phylogenetic relationships between superorders within the Acari. However, comparing gene order might be more powerful to establish phylogenetic relations within families, as was previously proposed [14,36]. In the case of the Ixodidae family, it was shown that the division of Prostriata (Ixodes sp.) and Metastriata (R. sanguineus, Amblyomma triguttatum, Haemaphysalis flava) could be linked to mt gene arrangements [37,38].

Base composition and codon usage
The overall AT-content of the mt genome of D. pteronyssinus is 72.6% (Table 1). This is within the range of the average AT-content of Acari mt genomes (74.6 +/-4.0%). The high AT-content is reflected in the codon usage ( Table  3 [ Metazoan mt genomes usually present a clear strand bias in nucleotide composition [40]. This is probably due to asymmetric patterns of mutations during transcription and replication when one strand remains transiently in a single-stranded state, making it more vulnerable to DNA damage [41]. However, in the case of mtDNA-replication, this hypothesis is not without controversy [42][43][44][45]. The strand bias in nucleotide composition can be measured as GC-and AT-skews ((G%-C%)/(G%+C%) and (A%-T%)/ (A%+T%), respectively) [46]. The overall GC-and ATskews of the J-strand of the D. pteronyssinus mt genome are 0.194 and -0.199, respectively. These are the most extreme values encountered within mite mt genomes up till now (Table 1) and they are reversed compared to the usual strand biases of metazoan mtDNA (negative GC-skew and positive AT-skew for the J-strand). Moreover, a positive GC-skew for mite mt genomes seems to be rare since at present, it was only encountered in Varroa destructor. Although hypothetical, it could be the result of a strand swap of the control region [40]. This region contains all initiation sites for transcription [47] and an inversion of the control region is expected to produce a global reversal of asymmetric mutational constraints in the mtDNA, resulting with time in a complete reversal of strand compositional bias [40]. The asymmetrical directional mutation pressure is also reflected in the codon usage of genes oriented in opposite directions [48]. Whereas NNG and NNU codons are preferred over NNA and NNC codons on the J-strand, genes on the N-strand show the exact opposite trend (see additional file 2 for an across-strand (N and J) comparison of frequencies of codons ending with the same nucleotide).
Schematic representation of the mt genome of D. pteronyssinus  (Table 1). Compared to other mite mt proteins, cox1, cox2 and cytB are best conserved. On the other hand, atp8, nad6 and nad4L showed lowest similarity values (see additional file 3 for the average identity and similarity % of mt proteins of D. pteronyssinus).
Start and stop codons were determined based on alignments with the corresponding genes and proteins of other mite species. In the case of stop codons, we could also benefit from available expressed sequence tags (ESTs) of D. pteronyssinus (n = 1797) and D. farinae (n = 1735) ( Fig.  1) [49]. As for other metazoan mt proteins, unorthodox initiation codons are used [20] (see additional file 4 for start and stop codons of protein coding genes of Acari mt genomes). Eight genes (cox2, atp6, cox3, nad3, nad6, nad4L, nad4, cytB) use the standard ATG start codon, 3 genes (cox1, nad1, nad2) start with ATA and nad5 initiates with ATT. atp8 most likely starts with codon TTG.
Eleven genes employ a complete translation termination codon, either TAG (cox1, cox3) or TAA (cox2, atp8, atp6, nad1, nad3, nad6, nad4L, nad5, cytB). With the exception of nad3, atp8 and nad4L, D. pteronyssinus ESTs for all these genes confirmed the position of the stop codon ( Fig. 1, see additional file 5 for sequences of Dermatophagoides ESTs covering the mt genome of D. pteronyssinus). Berthier et al. [50] showed that the adjacent genes, nad4L/nad4 and atp8/atp6, were transcribed and translated as a bicistronic mRNA in the model organism Drosophila melanogaster. However, as no ESTs were found that aligned with the nad4L/nad4 and atp8/atp6 gene boundaries, it could not be confirmed whether this was also the case for D. pteronyssinus. Despite its efficiency, the use of sequence alignments to determine the position of stop codons resulted in several cases in overlapping genes. For example, based on a highly conserved tryptophan at the C-terminal end of Acari nad3 proteins, a stop codon was positioned despite the resulting 17 bp overlap with trnaR. The two remaining genes (nad2 and nad4) are likely equipped with a truncated stop codon (T). Polyadenylation of the mRNA is needed in these cases to form a fully functional TAA stop codon [51]. Although speculative, ESTs of D. farinae confirm the truncated stop of nad4 (Fig. 1, see additional file 5).

Transfer RNAs
Fourteen tRNAs are encoded on the J-strand and 8 on the N-strand (Fig. 1). Secondary structures were predicted for all tRNAs (Fig. 4). With the exception of trnS 1 (UCU instead of GCU in L. pallidum) and trnP (UGG instead of AGG in S. magnus), all anticodon sequences were identical to those of L. pallidum and S. magnus, the only acariform mites for which tRNA secondary structures have been reported [22,24]. Usually, T is in the first anticodon position for tRNAs that recognise either four-fold degenerate codon families or NNR codons. G is usually in this position only to specifically recognize NNY-codons [52]. Except for trnM, all of the D. pteronyssinus mt tRNAs follow this pattern. trnM has the anticodon CAT (to recognise both ATG and ATA), which is the case for almost all animal mt systems [52] (Fig. 4).
Only one tRNA lacks the D-arm: trnS 1 , as is common for most metazoans. With the exception of trnC, trnV and trnS 1 , all tRNAs have T-arm variable loops (TV replacement loops) instead of the T-arm. Similar structures were found for tRNAs of L. pallidum [24] and S. magnus [22]. The absence of the T-arm is a typical feature for tRNAs of Chelicerata belonging to the orders of the Araneae, Scorpiones and Thelyphonida. However, other taxa within the Chelicerata (Amblypygi, Opiliones, Ricinulei, Solifugae and ticks) possess typical metazoan cloverleaf tRNAs [13]. Masta and Boore [13] suggested a multi-step evolutionary process in an attempt to understand how so many tRNAs in these chelicerate groups could lose their T-arm. According to this speculative theory, changes in mt ribosomes, resulting in the fact that the loss of arms from tRNAs was tolerated [53,54], and/or changes in specific elongation factors [55][56][57] are considered as a first step in this process. Figure 3 Mitochondrial gene arrangement of Limulus polyphemus, Dermatophagoides pteronyssinus and Steganacarus magnus. Graphical linearisation of mt genomes is presented according to [32]. Gene sizes are not drawn to scale. J stands for majority and N for minority strand. Protein coding and rRNA genes are abbreviated as in the abbreviations section. tRNA genes are abbreviated using the one-letter amino acid code, with L 1 = CUN; L 2 = UUR; S 1 = AGN; S 2 = UCN. White boxes represent genes with the same relative position as in the arthropod ground pattern, L. polyphemus. Light-gray boxes represent genes that changed positions relative to L. polyphemus; dark-gray boxes represent genes that changed both position and orientation. Circular dots between the genes of D. pteronyssinus represent conserved gene boundaries compared to L. polyphemus. Square dots between the genes of S. magnus represent conserved gene boundaries compared to D. pteronyssinus.  *bold numbers represent pairwise common interval distances between mt gene orders of Acari (37 genes in total), while italic numbers represent pairwise common interval distances between mt gene orders of Acari without tRNAs (15 genes in total) (L. pallidum and M. occidentalis were excluded from the dataset as the genomes contain duplicated genes whereas other genes are absent). Accession numbers of Acari mt genomes are listed in Table 1. a similar gene order as L. polyphemus, Ornithodoros sp. and C. capensis b similar gene order as H. flava and A. triguttatum c similar gene order as L. deliense d due to lack of tRNAs only protein coding gene and rRNA order comparison was possible * RSCU is the number of times a particular codon is observed relative to the number of times a codon would be observed in the absence of any codon usage bias [39].

Mitochondrial gene arrangement of Limulus polyphemus, Dermatophagoides pteronyssinus and Steganacarus magnus
Inferred secondary structures of the 22 mitochondrial tRNAs from D. pteronyssinus

). A maximum of 3 mismatches in this stem is found in trnR.
In contrast, almost all tRNAs (18) possess a completely matched 5 bp anticodon stem. trnC, trnS 1 and trnN have a single mismatch whereas trnY has two mismatches in this stem. All tRNAs, except trnL 2 , have a symmetric anticodon loop consisting of 2 bp up-and 2 bp downstream of the 3 bp anticodon. The anticodon loop of trnL 2 consists of 2 nucleotides preceding the anticodon and 3 nucleotides immediately following it. This kind of aberrant anticodon loops have also been reported for the twohumped camel Camelus bactrianus ferus (trnS 2 ) [58] and the scorpion Mesobuthus gibbosus (trnH and trnN) [59]. As mentioned before, sequences of some tRNAs overlap with neighbouring genes. The extreme examples are trnR, trnS 2 and trnV. trnR overlaps with the adjacent gene nad3 on the same strand for 17 bp at its 3'-end whereas trnS 2 overlaps with the adjacent gene trnM on the same strand for 12 bp at its 3'-end. trnV overlaps with the adjacent gene trnP on the opposite strand for 11 bp at its 3'-end and with trnK on the opposite strand for 7 bp at its 5'-start. Despite these overlaps, we consider these genes not likely to be pseudogenes. First of all, their sequence is relatively well conserved when compared to corresponding genes of other Acari. Secondly, besides sequence conservation they depict a conserved secondary structure. Thirdly, an EST [GenBank: CB284825] of the related species D. farinae was found corresponding to the region covering trnR, trnM and trnS 2 of D. pteronyssinus indicating that the genes are expressed (see additional file 6 for an alignment of trnR, trnM and trnS 2 of D. pteronyssinus with an EST of D. farinae). Finally, and most importantly, stem mismatches and sequence overlap are not uncommon for mt tRNAs of arachnids [13,60], and are probably repaired by a posttranscriptional editing process [54,61].

Non-coding regions
The largest non-coding region (286 bp) is flanked by trnF and trnS 1 . It is highly enriched in AT (91.61%) and can form stable stem-loop secondary structures. Based on these features, it possibly functions as a control region [20,62]. With the exception of T. urticae (95.45%), it has the highest AT-content of all Acari mt control regions ( Table 1). The position of the non-coding region differs from most insect and arachnid mt genomes, where the region is mostly located in close proximity to 12S-rRNA ( [62], see additional file 1).
Based on the sequence pattern, the control region can be subdivided in a repeat region and a stem-loop region. The first region (11,528 bp) contains several ATrepeats. In order to verify the exact number of repeats we resequenced this region. For this purpose, two flanking primers, Dp-Ms-F and Dp-Ms-R, were synthesised span-ning approximately 700 bp. The PCR product was cloned and ten independent clones were sequenced. This revealed that the number of AT-repeats varied between 7 to 28, suggesting that this domain can be considered as a microsatellite [63]. This is remarkable as a mt microsatellite was never reported before for species belonging to the Chelicerata. Also in metazoan mtDNA such microsatellites are rare and have, to our knowledge, only been reported for butterflies [64], a dragonfish, Scleropages formosus [65], a bat, Myotis bechsteinii [66], a turtle, Pelomedusa subrufa [67] and several seal species [68][69][70].
The second region (11,768 bp) holds two short palindromic sequences, TACAT and ATGTA, which are conserved in mt genomes of mammals [71] and fishes [65,72]. They can form a stable stem-loop structure (Fig.  5-A2), which might be involved as a recognition site for the arrest of J-strand synthesis [71]. Near this region other stem-loop structures could be folded (Fig. 5-A) but none of them had flanking sequences similar to those that are conserved in the control region of the mt genome of insects [62] and metastriate ticks [37].
As described before, four other stretches of non-coding nucleotides were found outside the control region. These short sequences can fold into stable stem-loop structures ( Fig. 5-B, C, D) which may function as splicing recognition sites during processing of the transcripts [73].  Table 1).

Ribosomal RNAs
The 12S-rRNA and 16S-rRNA genes of Leptotrombidium species (Acariformes: Trombidiformes: Trombiculidae) are 23.4% and 23.5% shorter than their counterparts in Drosophila yakuba. This substantial reduction is mainly caused by the loss of stem-loop structures at the 5'-end of the rRNA genes [74]. To identify whether similar domains are absent in the rRNAs of D. pteronyssinus, we constructed their secondary structures (Fig. 6 [75,76]). This revealed that the D. pteronyssinus 12S-rRNA indeed lacks similar stem-loops as L. pallidum, compared to D. yakuba. The structure also revealed 1 additional stem-loop (stem-loop 1) not present in 12S-rRNA of L. pallidum. Like in L. palli-dum, one stem-loop replaces three stem-loops (24, 25 and 26) whereas another replaces a region of four stem-loops (39, 40, 41 and 42) of the D. yakuba 12S-rRNA [74]. Based on the modelled structure in combination with an alignment of other acariform 12S-rRNAs, the greatest sequence conservation was found in the loop region of stem-loops 21 and 27 and the region between stem-loops 48 and 50. 16S-rRNA and 12S-rRNA secondary structures of the mitochondrial genome of D. pteronyssinus Figure 6  16S-rRNA and 12S-rRNA secondary structures of the mitochondrial genome of D. pteronyssinus. The numbering of the stem-loops is after de Rijk et al. [75] for 16S-rRNA and after van de Peer et al. [76] for 12S-rRNA. Blue coloured nucleotides show 100% identity when aligned to 12S-rRNA and 16S-rRNA genes from other Acariformes (as listed in Table 1). Inferred Watson-Crick bonds are illustrated by lines, whereas GU bonds are illustrated by dots.
It is still an open question how relatively well-conserved structures such as rRNAs can dramatically decrease in size while remaining functional. Wolstenholme et al. [53] and Masta [54] suggested a correlation between the occurrence of truncated rRNAs (compared to Drosophila) and the loss of the T-arm in tRNAs. The coincidence of short rRNAs and missing T-arms in tRNAs was also observed in S. magnus, L. pallidum and D. pteronyssinus. Other acariform mites like T. urticae, Ascoschoengastia sp. and Walchia hayashii also exhibit short rRNAs (Table 1) and the prediction of their tRNA secondary structures could further support this hypothesis. However, examples contradicting this hypothesis also occur e.g. pulmonate gastropods with tRNAs lacking T-arms have no truncated rRNAs. Therefore, it remains possible that truncation of both tRNAs and rRNA genes only reflects an independent trend towards minimisation of the mt genome as suggested by Yamazaki et al. [77].

Phylogenetic analysis
A phylogenetic tree was constructed based on nucleotide and amino acid sequences from all mt protein coding genes of Acari. The ILD-test [78] indicated a significant incongruence (P = 0.01) among data set partitions for nucleotide alignments and low congruence (P = 0.07) among data set partitions for amino acid alignments. A considerable debate exists on the utility of this test [79][80][81][82][83][84]. However, the principle of Kluge [85] implies that all data should always be included in a combined analysis for any phylogenetic problem and therefore we combined data partitions for both amino acid and nucleotide alignments for phylogenetic analysis. A maximum parsimony (MP) analysis based on nucleotide alignments (data not shown) grouped V. destructor (Parasitiformes) within the Acariformes, close to D. pteronyssinus. This is in contrast with the generally accepted view on the phylogeny of the Acariformes and Parasitiformes [2,3]. As mentioned before, V. destructor and D. pteronyssinus both have a reversal of asymmetrical mutation pattern. When such reversals occurred independently, D. pteronyssinus and V. destructor could have acquired a similar base composition and as a consequence group together due to the longbranch attraction (LBA) phenomenon [40,86]. Modelbased methods such as maximum likelihood (ML) and Bayesian inference (BI) are less sensitive to LBA [40,87] and were for this reason considered for phylogenetic analysis.
ML and BI analysis performed on the amino acid data set resolved trees with an identical topology (Fig. 7-A) in which D. pteronyssinus clusters with S. magnus, forming a sistergroup of the Trombidiformes. This is in agreement with the most recent views on the classification of the Acariformes [2][3][4]. The nucleotide data set resulted in similar trees, confirming the evolutionary position of D. pteronyssinus (Fig. 7-B). The only major inconsistency over the trees was the position of T. urticae. Although this species is generally considered as a member of the Trombidiformes [2][3][4], it was clustered with the sarcoptiform mites D. pteronyssinus and S. magnus in the trees based on the nucleotide dataset. (Fig. 7-B). However, the position in the different trees is questionable as it is supported by low bootstrap values/Bayesian posterior probabilities (Fig. 7-A/B). Adding additional mt genome data from closely related taxa of T. urticae and from taxa located between T. urticae and Trombiculidae would probably position T. urticae with higher support values within the Trombidiformes.
In the trees based on the nucleotide dataset, H. flava is, compared to A. triguttatum, evolutionary closer related to R. sanguineus while in the trees based on the amino acid dataset this is the opposite. However, as the clustering of H. flava and R. sanguineus is in agreement with the most recent views on the classification of the Ixodida [88,89], we consider the nucleotide topology as the most correct one. Murrell et al. [90] considers the Parasitiformes to be paraphyletic with respect to the Opilioacariformes, but as there are no complete mt genomes of Opilioacariformes available, we were not able to verify this hypothesis.

Conclusion
This is the first description of a complete mt genome of a species belonging to the Astigmata, a cohort within the Sarcoptiformes. Although the length, gene and AT-content are similar to other Acari mtDNA, the mt genome of D. pteronyssinus exhibits some interesting features. The gene order of D. pteronyssinus is completely different from that of other Acari mt genomes. Gene order comparison indicated that mt gene orders seem less useful for deduction of phylogenetic relationships between superorders within the Acari. GC-and AT-skews of the J-strand were very large and reversed as compared to those found in most metazoan mtDNA.
Compared to parasitiform mites, both D. pteronyssinus rRNAs were considerably shorter and almost all transfer RNAs lacked the T-arm. It would be interesting to investigate whether the occurrence of truncated rRNAs and the loss of the T-arm in tRNAs are correlated or just a trend Phylogenetic trees of Acari relationships Figure 7 Phylogenetic trees of Acari relationships. Trees were inferred from amino acid (A) and nucleotide (B) datasets. All protein coding gene sequences were aligned and concatenated; ambiguously aligned regions were omitted by Gblocks 0.91b [105]. Trees were rooted with two outgroup taxa (L. polyphemus and L. migratoria). Numbers behind the branching points are percentages from Bayesian posterior probabilities (left) and ML bootstrapping (right). Accession numbers for the different Acari mt genomes are listed in Table 1. toward minimisation of the mt genome. Finally, phylogenetic analysis using concatenated mt gene sequences succeeded in recovering Acari relationships concordant with traditional views of phylogeny of Acari.

Mite identification
Upon arrival in the laboratory, mites were identified as D. pteronyssinus by J. Witters (ILVO, Belgium) and F. Th. M. Spieksma (Laboratory of Aerobiology, LUMC, The Netherlands) using morphological characteristics. To back up this identification, molecular techniques were applied. For this purpose DNA was extracted and used as a template for PCR. Primers 12SID-F and 12SID-R (see additional file 7 for primer sequences) successfully amplified a 316 bp fragment. BLASTn searches against non-redundant nucleotide sequences using the amplified fragment as query resulted in a perfect match with a mt 12S-rRNA sequence of D. pteronyssinus [GenBank: AF529911].

Mite strain, mass rearing and isolation
The initial D. pteronyssinus culture was provided by D. Bylemans (Janssen Pharmaceutica, Belgium). Mites were cultured on a 1:1 mixture of Premium Gold (Vitacraft, Germany) and beard shavings at 75% R.H., 25°C and permanent dark conditions [91,92]. Mites were isolated from the colony using a modified heat-escape technique [93,94]. Briefly, mite cultures were transferred to small plastic petri dishes (75 mm in diameter, 28 mm high) with a lid on top. These dishes were placed in the dark on a hot plate set at 45°C (Bekso, Belgium). After 15-20 minutes the mites moved away from the heat source, formed groups on the lid of the petri dish and could be collected using a fine hair brush.

DNA extraction
Approximately 1000 D. pteronyssinus mites were collected in an Eppendorf tube and were ground in 800 μl SDS-lysis buffer (400 mM NaCl, 200 mM TRIS, 10 mM EDTA, 2% SDS) using a small sterile plastic pestle (Eppendorf, Germany). After incubation for 30 min at 60°C under continuous rotation, a standard phenol-chloroform extraction was performed [95]. Total genomic DNA was precipitated with 0.7 volumes of isopropanol at 4°C for 1 hour, centrifuged for 45 minutes at 21,000 × g and washed with 70% ethanol. Precipitated DNA was resolved in 50 μl 0.1 M Tris pH 8.2.

Amplification of the mt genome
Primers COXI-F and 12S-R, based on partial D. pteronyssinus cox1 and 12S-rRNA sequences [GenBank: AY525570 and AF529911, respectively] (see additional file 7 for primer sequences), successfully amplified a 4.6 kb sequence of the mt genome of D. pteronyssinus. Degenerate primers CYTB-F-Deg and CYTB-R-Deg (see additional file 7 for primer sequences), designed on conserved regions of Acari cytB, amplified a partial cytB sequence from D. pteronyssinus. A specific primer COXI-R, designed from the 3' end of the 4.5 kb sequence in combination with the primer CYTB-F, designed from the partial cytB sequence, successfully amplified a 2.2 kb sequence. Another primer CYTB-R, designed from the 5'-end of this 2.2 kb sequence, in combination with the primer 12S-F successfully amplified a 8.6 kb sequence, making the mt genome sequence complete.

Annotation and bioinformatics analysis
The complete genomic sequence was assembled and annotated using VectorNTI (Invitrogen, Belgium) according to Masta and Boore [60]. Open reading frames (ORFs) were identified with the program Getorf from the EMBOSS-package [96]. The obtained ORFs were used as query in BLASTp [97] searches against the non-redundant protein database at NCBI. Two large non-protein-coding regions were candidates for the rRNAs (16S and 12S respectively). The boundaries were identified based on alignments and secondary structures of rRNA genes of other mite species. Sixteen of the 22 tRNAs were identified by tRNA-scan SE [98] with a cove cutoff score of 0.1 and the tRNA-model set to "nematode mito". The remaining tRNAs (trnM, trnV, trnY, trnS 1 , trnI, trnC) were determined in the unannotated regions by sequence similarity to tRNAs of other mite species. In order to obtain additional information on mt gene boundaries, BLASTn [97] searches of D. pteronyssinus tRNA, rRNA and protein encoding nucleotide sequences were carried out against ESTs [49] restricted to Dermatophagoides sequences (n = 3532). ESTs with statistically significant matches (E-value cutoff: 0.1) were collected, checked for vector contamination and aligned by Clustal W [99] as implemented in BioEdit 7.0.1 [100] against the appropriate nucleotide sequence of D. pteronyssinus. MatGAT 2.02 was used to calculate similarity and identity values [101] of mt proteins. The identification of gene subsets that appear consecutively in different genomes was performed by common interval distance analysis using CREx [102] (see additional file 8 for input data of the CREx program).

Construction of secondary structures of RNAs and noncoding regions
Secondary structures of tRNAs were determined following the method of Masta and Boore [60]. Secondary structures of tRNAs were drawn with CorelDraw 12.0 (Corel Corporation, Canada). The rRNA genes of D. pteronyssinus were aligned with those of other Acariformes and conserved areas were identified. These regions were mapped on the published structures of L. pallidum rRNA [74]. Regions lacking significant homology were folded using Mfold [103]. Secondary structures of rRNAs were drawn using the RnaViz2 program [104] and afterwards modified with CorelDraw 12.0 (Corel Corporation, Canada). Secondary structures of non-coding regions were folded using Mfold [103]. When multiple secondary structures were possible, the most stable (lowest free energy (-ΔG)) one was preferred. Drawing and editing of these structures was done in a similar way as for rRNA secondary structures.

Rolling circle amplification and restriction enzyme digestion
Extraction and rolling circle amplification of the mtDNA of D. pteronyssinus was done according to Van Leeuwen et al. [18]. Rolling circle amplified mtDNA was digested with two enzymes (XmnI and EcoRI; New England Biolabs) following the manufacturer's instructions. Restriction digests were fractionated by agarose gel electrophoresis as described before.

Phylogenetic analysis
Sequence data were obtained from 21 Acari species (for GenBank accession numbers see Table 1) and two outgroup taxa (Limulus polyphemus [GenBank: NC_003057] and Locusta migratoria [Genbank: NC_001712]). Only mite species with a completely sequenced mt genome were selected. Alignments from all mt protein-coding genes were used in phylogenetic analysis. Amino acid sequences and nucleotide sequences were aligned by Clustal W [99] as implemented in BioEdit 7.0.1 [100]. The nucleotide alignment was generated based on the protein alignment using codon alignment. Ambiguously aligned parts were omitted from the analysis by making use of Gblocks 0.91b [105], with default block parameters except for changing "allowed gap positions" to "with half". Abascal et al. [106] recently presented evidence that some insects and ticks use a modified mitochondrial code, with AGG coding for lysine rather than serine as in the standard invertebrate mitochondrial code. As 10 out of 20 Acari species in our dataset are ticks all positions aligning to AGG codons in the final amino acid alignment were removed.
For the nucleotide alignments the "codons" option was used in Gblocks 0.91b [105]. Due to the results of a saturation analysis [107] on single codon positions, implemented in DAMBE 4.2.13 [108], third codon positions were eliminated from the nucleotide alignment. An incongruence length difference test (ILD-test) [78] as implemented in PAUP* (version 4.0b10; [109]) was used to assess congruence among gene partitions.
Model selection was done with ProtTest 1.4 [110] for amino acid sequences and with Modeltest 3.7 [111] for nucleotide sequences. According to the Akaike information criterion, the mtART+G+I+F model was optimum for phylogenetic analysis with amino acid alignments and the GTR+I+G model was optimal for analysis with nucleotide alignments.
Two different analyses were performed. (1) Maximum likelihood (ML) analysis was performed using Treefinder [112], bootstrapping with 1000 pseudoreplicates (2) Bayesian inference (BI) was done with MrBayes 3.1.2 [113]. As the mtART model is not implemented in the current version of MrBayes, the mtREV+G+I model was used for phylogenetic analysis with the amino acid alignment. Four chains ran for 1,000,000 generations, while tree sampling was done every 100 generations. Burnin was calculated when the average standard deviation of split frequencies had declined to < 0.01. The remaining trees were used to calculate Bayesian posterior probabilities (BPP).