Phylogeny of Terrestrial Isopods Based on the Complete Mitochondrial Genomes, Subvert the Monophyly of Oniscidea and Ligiidae up to New Subfamily Ligiaidea of Isopoda

Background: Oniscidea is the only truly terrestrial taxon within the Crustacea, and vital to soil formation. However, the monophyly of suborder Oniscidea has been in dispute since 1995, with different studies disagreeing on whether the coastal Ligiidae are included within the suborder. To clarify the phylogenetic hypothesis of suborder Oniscidea, we sequenced the complete mitochondrial genomes of Ligia exotica (Roux, 1828) and Mongoloniscus sinensis (Dollfus, 1901). Results: Like most metazoan, the complete mitogenomes of two species with circular double strands. The structure and characters of mitogenomes of these two species are analyzed. The constructed phylogenetic analyses show that Oniscidea is polyphyletic group, with Ligia being more closely related to marine isopods (Valvifera + Cymothoida + Sphaeromatidea). Conclusions: We elevate the taxonomic status of the family Ligiidae to the suborder Ligiaidea which are with parallel rank with Oniscidea. Ligiaidea is much primitive than other exact terrestrial isopods. Crinocheta are strongly monophyly, family is more closely related to Porcellionidae rather than Armadillididae.

most basal split within a monophyletic Oniscidea clade, thus casting doubt not only in the placement of Ligiidae within the Oniscidea but also in our historical understanding of the evolutionary history of the terrestrial isopods. Considering additional support is needed for these ndings, in this study, we analyze the mitochondrial genomes (hereafter mitogenomes) of two Ligia species, as well as an additional 22 other isopod species, to explore the evolutionary path of the Oniscidea. Our goal is to clarify the phylogenetic relationships amongst the Oniscidea.

Genome organization and base composition
The mitochondrial genome of M. sinensis and L. exotica are circular, double-stranded DNA moledules, 16,018 and 14,978 bp in length respectively (Figs. 1 & 2). The complete mitogenome of M. sinensis contains 34 mt genes: 13 protein-coding genes (PCGs), 19 tRNA genes, and two rRNA genes (Fig. 1). It lacks three tRNA genes: TrnA, trnE, trnL1. The average A + T content of the M. sinensis mt genome is approximately 75.32%, which is higher than other isopods (typical range: 54.4%-71.2 ). 18 overlapping regions and 13 intergenic regions are found in the genome. (Table 1). Nucleotide frequencies of all mt genes of M. sinensis are listed ( Table 2).  The circular mitogenome of Ligia exotica is composed of 13 PCGs, 21 tRNA genes, two rRNA genes, one non-coding region, while only lacking the trnG gene (Fig. 2). The average A + T content of the L. exotica mt genome is approximately 59.13%. 17 overlapping regions and 13 intergenic regions are found in the genome (Table 3). Nucleotide frequencies of all mt genes of M. sinensis are listed (Table 4).   PCGs of M. sinensis are 11,088 bp in size, with its A + T content reaching 74.18%, the highest among all known Oniscidea. The ATG codon is the most commonly found start codon (six PCGs), with ATT next (found in nad3, nad4I, nad5, nad6, and atp8), and ATA only found in nad1. As for the terminal codons, TAA is found in 12 genes, and TAG only found in nad2 ( Table 1).
The PCGs of L. exotica are 11113 bp in size, with its A + T content being 58.06%, much lower than M. sinensis. ATA is the start codon for four PCGs (cob, nad4I, cox2, nad8), with ATG found in four PCGs (nad2, nad4, atp6, cox3), ATT in three (nad3, nad5, nad6), and ATC and CGA in one each (nad1 and cox1 respectively). As for terminal codons, TAA is found in 7 genes, TAG in ve (cox2, nad2, nad3, nad5, nad6) and T found only in Nad4 (Table 3). The un nished T codon is not counted separately, as we presumed that it would be completed ( The two species exbibit the same protein-coding genes location as 9 genes (cox1-3, atp8, atp6, nad2-3, and nad5-6) are encoded by the plus strand and four genes (cob, nad1, nad4, and nad4L) by the minus strand ( Table 2 and Table 4 Transfer and ribosomal RNA genes The two rRNAs, rrnL and rrnS, of M. sinensis are 1143 and 717 bp in size, with 79.79% and 68.48% A + T content respectively. All 19 commonly found tRNAs are present in the mitochondrial genome of M. sinensis, ranging from 51 bp (trnS2) to 73 bp in size (trnS1), and adding up to 1158 bp in total combined length. As for L. exotica, the two rRNAs, rrnL and rrnS, are 1184 and 699 bp in size, with 62.58% and 57.94% A + T content, respectively. All 21 commonly found tRNAs are present in the mitochondrial genome of L. exotica, ranging from 53 bp (trnC) to 64 bp in size (trnI, trnM, trnW), and adding up to 1278 bp in combined length. tRNA genes are distributed throughout the mitogenome and are found on both strands. The putative secondary structures of all identi ed tRNAs are shown in Fig. 3 and Fig. 4. The majority of tRNAs have a common t-shaped or clover-leaf secondary structure. Exceptions include the trnC of M. sinensis, where the DHU-arm is absent, and the trnG, trnK, trnP, trnS2, trnW and trnV, where the TΨC-arm is absent. In L. exotica, all of the secondary structures (predicted by MITOS and ARWEN) exhibit the conventional cloverleaf structure, except for trnC and trnV, which lack DHU arms and trnF, trnP, trnM, trnL2, trnK, trnD, trnR which lack the TΨC-arm.

Non-coding regions
In M. sinensis. there is one biggest non-coding region of 879 bp length located between trnT and trnS1. And other three non-coding regions are between trn-W and cob the length is 292 bp, 485 bp and 341 bp.
We have not detected a hairpin structure in the mt control region of M. sinensis (Fig. 1). There are six discontinuous repeats of length 53 bp. As for L. exotica, one non-coding regions (NCR), 662 bp in size, are located between trnE and trnW. There are no tandem repeats of sections in this region and a GC-rich region containing the putative hairpin structure (Fig. 5).

Phylogeny and gene order
The phylogenies produced using BI and ML methods show concordant topologies; however, support values differ between approaches. BI analyses produced very high statistical support while the ML topology exhibited a mixture of mostly high but several lower support values ( Fig. 7 and Fig. 8). Neither analyses recovered the monophyly of the Oniscidea. Though M. sinensis formed a clade with three species of Porcellionidae Brandt Ratzeburg,1831 (Oniscidea) which was placed within a larger clade including all other Oniscidean species included in this study except Ligia species. But L. exotica formed a clade with L. oceanica, and the clade of two species of family Ligidae was not placed into the main clade consisting subfamily Oniscidea.
The gene order of the mt monomer of M. sinensis and L. exotica are shown in Fig. 6. Comparison with the putative isopoda ground pattern mt genome [11] revealed four gene rearrangements in M. sinensis, three gene rearrangements in L. exotica. As for M. sinensis, the rst rearrangement trnR is between the cox3 and nad3, and it is between nad3 and nad1 of isopoda ground pattern. The second rearrangement, trnV, is between nad3 and nad1, but it is between 16S rRNA and nad3 in the isopoda ground pattern. The third rearrangement of trnT is between the 12S rRNA and cob, and it is between cob and nad5 in isopoda ground pattern. Another rearrangement of trnT was translocated near trnS1 and cob.
As for L. exotica, the rst rearrangement trnR is between the cox3 and nad3, and it is between nad3 and nad1 of isopoda ground pattern; The second and third rearrangement trnW and trnE are interchanging. In a word, comparing with isopoda ground pattern, two isopods have gene translocation as other known mitochondrial genomes in isopod. Missing tRNA is found in all oniscideas including the present two species, and it is universal phenomenon. So, the order of mitochondrial genes in the Oniscidea is weakly conserved.

Discussion
The object of this work was to reconstruct the evolutionary relationships of the Oniscidea based on complete mitogenomes, with a particular focus on the Ligiidae as Ligia have been proposed to represent the intermediate forms in the evolution of terrestriality in isopods [7].
Phylogeny of Oniscidea based on the morphological characters or molecular data are debated [12,18,19,20,21]. The monophyly of Oniscidea has been well supported by morphological characters such as the complex water-conducting system and reduced rst antenna with only three articles [6,8,9,12,22,23]. The rst antenna with rudimentary three articles has been regarded as a prominent synapomorphy for Oniscidea [6,9], but some have suggested it to be a possible homoplasy [6]. Hornung (2011) found the water-conducting systems within Oniscidea are changeable [24]. So, the two main morphological characters to support the monophyly of Oniscidea actually be convergent related to the terrestrial environment challenges [12]. Indeed, monophyly of Oniscidea have been denied by various molecular data [15,21,25,26,27]. But further phylogenetic analyses of suborder Oniscidea and taxonomic rank and Our phylogenetic tree does indicate that the other ten species of oniscidea are indeed closely grouped together. Three species of Armadillidiidae Brandt,1833 and Porcellionidae Brandt & Ratzeburg,1831 are also closely clustered. Monophyly of Agnaridae is supported by two mitochondrial and three nuclear genetic markers [28], and it is closely related to Porcellionidae, rather than Armadillidiidae. M. sinensis formed a clade with three species of Porcellionidae (Oniscidea), and three species of Armadillidiidae formed another clade. So, our ndings are in agreement with above research, but add further evidence to de nite the relationship between families of Oniscidea. Oniscus asellus, with its longer horizontal branches, morphologically with much atter body and more varied body than other isopods. With regards to M. sinensis, they formed a clade with three species of Porcellionidae Brandt & Ratzeburg,1831 and were indeed much closer in morphology and life environment. So, Trachelipus rathkii clustered with Cylisticus convexus, and two different heterogeneous nodes are found, by changing the anti-codon of tRNA to produce additional double tRNA genes, the results are consistent in both species. This may be the reason why two species have always been one clade. So, monophyly of Crinocheta is strongly supported, and the close relatives of Ligiaidea are Valvifera + Cymothoida + Sphaeromatidea.
The sequence of mitochondrial genomes of M. sinensis and L. exotica are highly derived compared with other known isopod species pattern. All 24 oniscideans studied to date have missing tRNAs (see Suppl. Materials 1), possibly due to the short length of tRNA. We report the complete absence of the DHU arm for the TrnC of M. sinensis and trnC, trnV of L. exotica. Complete absence of tRNA DHU arm is also found in other isopods, for instance the trnC and trnV in L. oceanica. Doublet et al., (2012) and Chandler et al., (2015) pointed out mitogenomes of isopods especially oniscideans with unique linear/circular organization [27,29]. It may be associated with an ancient and conserved constitute heteroplasmic sites, or faciliate the evolution of stable mitochondrial heteroplasmies, and then constrain the multimeric structure of mitogenomes. The heteroplasmic site alters the anticodon in tRNA genes, allow it to code for different tRNA. As for M. sinensis and L. exotica, some tRNA genes appear to be missing not by annotation software or RNA editing, and it is usual in isopods (Suppl. Materials. 1) [10,30].
In Ligia exotica, Ligia oceanica and Idotea baltica, there is a positive GC skew in the plus strand coding gene and a negative GC skew in the minus strand. Most of the genes in M. sinensis follow the above rules, and only some of the genes are reversed, such as ad4l and trn-l are in the negative chain but their GC is positively skew, and atp8 is in the positive chain but its GC is negatively skew. This is common in most other crustaceans, probably because of the reversal of the mitochondrial gene replication region.

Conclusions
The family Ligiidae with its own suborder Ligiaidea is proposed. Monophyly of Agnaridae is strong supported by the complete mitochondrial, and it is closely related to Porcellionidae rather than Armadillidiidae. Overall, the conquest of the land is one of the major events in the history of the life on Earth, so the key groups from the ocean to land should attract much more attention. In this study we add evidence to the growing consensus behind the polyphyly of Oniscidea. Thus, taxonomic revision of this group should be urgently considered. Mitogenomes and nuclear genes of more species such as genera including Caucasoligidium, Ligidoides, Ligidium, Tauroligidium, Typhloligidium and other current and previous members of the Ligiidae family should be sequenced. Morphological and ecological characters of isopods should be analyzed with multi-evidences.

Sample and DNA extraction
As only eight isopods have their complete mitogenomes and other 14 isopods with incomplete mitogenomes reported in public databases (Table 1)

Genome determination, gene annotation and sequence analysis
Libraries were constructed using Whole Genome Shotgun (WGS) strategies with libraries sequenced using paired-end (PE) approaches on the Illumina MiSeq sequencing platform. Raw data was quality tested using the Read Quality Inspection Tool FastQC (http://www.bioinformatics.ac.uk/projects/fastqc). The data quality assessment was carried out by inspecting the single base quality, base content distribution, GC content distribution, and sequence base quality. High-quality second-generation sequencing data was assembled by using A5-miseq v20150522 and SPAdesv3.9.0 to construct contig and scaffold sequences.
Newly produced sequences were extracted according to the sequenced depth of splicing-sequence using blastn (BLAST v2.2.31+) in nt bank of NCBI to the mitochondrial sequence. The mitochondrial spliced results obtained by the above different software were co-linearly analyzed by mummer v3.1 software, and the positional relationship of contigs was determined, with gaps between lled.
Results were corrected with the pilon v1.18 software to obtain the nal mitochondrial sequence. Pilon signi cantly improves draft genome assemblies by correcting bases, xing mis-assemblies and lling gaps. Mitogenomes were annotated using the MITOS webserver [32] with tRNAs re-detected using two  Map of the Mongoloniscus sinensis. mitogenome. From inside to outside, the rst circle represents the scale; the second GC skew (-); the third GC content; the fourth the protein-coding gene, tRNA, rRNA and its arrangement on the genome. Putative secondary structures for the 21 transfer RNAs of the Ligia exotica. mitogenome. All 21structures were generated using the MITOS webserver and veri ed using tRNAscan-SE and ARWEN. Most tRNAs feature a standard clover-leaf structure. Exceptions: The TΨC arm is absent from trnF, trnP, trnM, trnL2, trnK, trnD,trnR and the DHU arm is absent from trnC,trnV Hairpin structures in mitochondrial control regions of Ligia exotica (Isoppda).

Figure 6
Gene order in isopod mitogenomes.

Supplementary Files
This is a list of supplementary les associated with this preprint. Click to download.