Phylogenetic position of Bopyroides hippolytes, with comments on the rearrangement of the mitochondrial genome in isopods (Isopoda: Epicaridea: Bopyridae)

Background Classification of parasitic bopyrids has traditionally been based on morphological characteristics, but phylogenetic relationships have remained elusive due to limited information provided by morphological data and tendency for loss of morphological features as a result of parasitic lifestyle. Subfamily Argeiinae was separated from Bopyrinae based on morphological evidence, although the assignment of all genera has not been phylogenetically evaluated. Bopyroides hippolytes has been traditionally classified in Bopyrinae, but divergent morphological characters make this assignment questionable. To investigate the relationship of bopyrines, we sequenced the complete mitochondrial genome of B. hippolytes and four mitochondrial genes of two other Bopyrinae. Results The phylogenetic trees based on separate and combined cox1and 18S sequence data recovered Bopyridae as robustly monophyletic, but Bopyrinae as polyphyletic. Bopyroides hippolytes was a close sister to Argeia pugettensis, type species to Argeiinae. Mitochondrial phylogenomics also suggested that B. hippolytes was close to Argeiinae. We also found a novel gene order in B. hippolytes compared to other isopods. Conclusions Bopyroides hippolytes should be excluded from the Bopyrinae and has a close affinity with Argeia pugettensis based on molecular and morphological data. The conserved syntenic blocks of mitochondrial gene order have distinctive characteristics at a subordinal level and may be helpful for understanding the higher taxonomic level relationships of Isopoda. Supplementary Information The online version contains supplementary material available at 10.1186/s12864-022-08513-9.


Introduction
Bopyridae Rafinesque, 1815 is a parasitic family of isopod crustaceans, with 10 subfamilies, 208 genera and 634 species (Boyko et al., 2008 onwards) [1]. Bopyrids stunt the growth and reduce the reproductive abilities of their hosts [2]. Molecular phylogenetic research on Bopyridae has been scarce, partly because of the challenge of assembling sequenceable material. Molecular phylogenetic analyses began only recently. Boyko et al. [3] assessed phylogenetic relationships among Bopyridea Rafinesque, 1815 and Cryptoniscoidea Kossmann, 1880 using 18S sequence data, and recovered Bopyridae as monophyletic. Yu et al. [4] and An et al. [5] sequenced and analyzed Open Access *Correspondence: anjianmei@sxnu.edu.cn; anjianmei@hotmail.com † Ruiwen Wu and Rongxiu Guo contributed equally to this work. 1 School of Life Science, Shanxi Normal University, Taiyuan 030031, P. R. China Full list of author information is available at the end of the article the mitochondrial genome of the bopyrids Gyge ovalis (Shiino, 1939) and Argeia pugettensis Dana, 1853, respectively, but the paucity of comparative sequence data prevented resolution of bopyrid relationships at the subfamily level.
With some exceptions, members of each bopyrid subfamily are restricted to hosts from one decapod infraorder [6]. Bopyrinae, Argeiinae and Pseudioninae infest the branchial chambers of caridean shrimp. Shiino [7] grouped a series of bopyrid genera, including Argeia Dana, 1853, Stegoalpheon Chopra, 1923, Bopyrosa Nierstrasz & Brender à Brandis, 1923, Parargeia Hansen, 1897 into his Bopyrus-group, later recognized as the Bopyrinae. Markham [8] placed these four genera together with Argeiopsis, into his newly erected subfamily Argeiinae. He noted that Argeiinae can be distinguished from Bopyrinae in that the female of the latter has a head that is not oval or fusiform and is usually fused with the pereon, an oval or deltoid body outline, some or all pleomeres fused at least on one side, lateral plates and uropods that are greatly reduced or absent, and pleopods that are generally biramous.
Bopyroides Stimpson, 1864 is currently placed in the Bopyrinae but does not fit well there. Females have an oval head, distinctly separated from pereon, all pleomeres are distinct without any fusion, have prominent dorsolateral bosses and tergal projections, and pleopods are reduced or uniramous tuberculiform [9][10][11]. Thus, there remains uncertainty about the boundaries between the Bopyrinae and Argeiinae, and the validity of Argeiinae has also been questioned [12].
Key to species of Bopyroides: 1 Pleomere 6 of female with lateral, biramous exten-sion…………………B. cluthae -Pleomere 6 of female not laterally extended, round or truncate……………2 2 Female with four pairs of uniramous pleopods, head of male with straight posterior margin…B. hippolytes -Female with five pairs of uniramous pleopods, head of male with curved posterior margin…B. shiinoi We sequenced the mitogenome and 18S rRNA of three species currently assigned to the Bopyrinae: Bopyroides hippolytes, Bopyrella malensis Bourdon, 1980, and Parabopyrella cf. mortenseni (Bourdon, 1980), and analyzed their phylogenetic placements among bopyrids based on this and other available sequence data. We also compared gene arrangements across all available mitogenomes of the Isopoda and describe a novel mitochondrial gene order in bopyrids.

Table 1 Classification history and synonym for the three Bopyroides species
We addressed the following questions: (1) Is the Bopyrinae monophyletic? (2) What is the phylogenetic position of Bopyroides hippolytes? (3) How has mitochondrial gene rearrangement in the Isopoda?

Sequence alignment and data partitions
Sequence information for the four datasets after alignment and trimming are shown in Table 2. The subset partitions and best-fit models from Partition Finder and Model Finder are presented in Table S1.

Phylogenetic position of B. hippolytes and molecular phylogeny of Bopyridae
The six phylogenetic trees obtained based on separate and combined cox1 and 18S sequence data with ML and BI were congruent ( Fig. 1). All analyses recovered most species currently assigned to the Bopyrinae (Probopyrus pandalicola, P. pacificiensis, P. buitendijki, Parabopyrella mortenseni, Bopyrella malensis) as robustly monophyletic, except for Bopyroides hippolytes, which was always a close sister to Argeia pugettensis. Our Bopyroides hippolytes sequence matches others in GenBank, further confirming the identification. The relationship of Argeiinae + Hemiarthrinae had strong support.

Mitochondrial phylogenomics of Isopoda
Phylogenetic analysis based on AA sequences of the 13 PCGs recovered Phreatoicidea, then Asellota, as sister to all other sampled suborders in the Isopoda, albeit with modest support. All suborders except Cymothoida were monophyletic, although the polyphyly of Cymothoida also lacked robust support (Fig. 2). Bopyroides hippolytes was recovered as sister to Argeia (Argeiinae), distant from the other two Bopyrinae species.

Novel mitochondrial gene arrangements
Bopyroides hippolytes has a unique mitochondrial gene order that differs from that of all other known in Isopoda (Fig. 3). The circular mitochondrial genome of B. hippolytes is 15,329 bp, with C + G content of 39.9%. It possesses the standard 13 PCGs and two rRNA genes (12S and 16S), but only 19 tRNA genes, as three tRNA genes (trnK, trnW and trnI) are missing compared to the standard animal mitogenome (Fig. 3). Nineteen genes are encoded on the +strand (light strand), whereas 4 protein-coding genes, 10 tRNAs and one rRNAs are located on the -strand (heavy strand). In base composition, B. hippolytes exhibits negative AT skew (− 0.129) and positive GC skew (0.039).
Isopods have undergone numerous gene duplication and deletion events in the mitogenome among taxa sampled to date, except for the basal Phreatoicidea, which retains a standard mitochondrial genome of 37 genes. While only tRNA genes were lost or duplicated in most isopod clades, PCGs were also duplicated in Argeia pugettensis. Expect Limnoria quadripunctata, the putative CR (control region) of all available isopod mitogenomes is located between the rrnS (+strand) and cob (−strand) genes.

Phylogenetic relationships and taxonomic implications
Markham [8] separated the Argeiinae from the Bopyrinae based on several morphological features, including the shape of the head and pleopods, but considered these subfamilies to be closely related. Bopyroides hippolytes was originally described in Bopyrus, the type genus of Bopyrinae (Kröyer, 1838). While this widespread and common species has received substantial attention, and was transferred to its own genus, its subfamilial classification has not been discussed and remained accepted as Bopyrinae [6, 11, 13, 17-20, 22, 23].
The Bopyrinae was polyphyletic in all our phylogenetic analyses, with Bopyroides hippolytes separating from other bopyrine genera, but close to Argeia pugettensis. We suggest that Bopyroides hippolytes should be excluded from Bopyrinae and has a close affinity with Argeia pugettensis. This conclusion is also supported by some morphological data, such as the dorsolateral boss and tergal projection of Bopyroides hippolytes, described by Shiino [19], Allen [20], and Rybakov and Avdeev [11] that are characteristics of Argeiinae, not Bopyrinae. The two other species assigned to Bopyroides, B. cluthae and B. shiinoi, also have distinct pleomeres, prominent dorsolateral bosses and tergal projections [9,11,20] and likely are indeed congeneric. So, the boundary between Bopyrinae and Argeiinae is obscure, and the correct rank of Bopyroides need more data.
Our phylogenetic results are congruent with Boyko et al's [3] analysis of bopyrids using 18S sequence data but extends it with greater sampling of Bopyrinae and Argeiinae. These two subfamilies are well separated on the phylogeny suggesting that their similarities perceived by Markham [24] are the result of convergence. There remain limited sequence data for this large family of parasitic isopods and further work will likely lead to additional changes in their classification.

Gene rearrangement
Mitochondrial gene rearrangement provides useful information for understanding relationships at higher taxonomic levels [25,26]. Duplication and deletion of tRNA are common in the rearrangement of mitochondrial  [27][28][29], likely caused by the replication slippage mechanisms [30,31].
Most of the available isopod mitochondrial genomes in GenBank are partial or incomplete, with only 11 mitogenomes being complete (as of 30 Dec. 2020). We evaluated gene rearrangements by assembling the 12 available complete isopod mitogenomes, including Bopyroides hippolytes sequenced in this study. Among all complete isopod mitogenomes analyzed, tRNA deletions occur in all species except Limnoriidea and Phreatoicidea. These suborders are early cladogenesises among isopods and appear to retain the ancestral complement of tRNAs. Two general changes in the tRNA genes of all isopods include a reduction or loss of trnI and reduction of trnC. trnI was missing in both Epicaridea and Cymothoida (Fig. 3), and while it was retained in the other suborders, its cloverleaf structure was incomplete (loss of D-loop or TΨC). The D-loop region was lacking in the predicted secondary structures of trnC for all mitogenome analyzed. This phenomenon was also described by Kilpert and Podsiadlowski [25], who considered that this feature might be a putative autapomorphy of Isopoda. It is noticeable that an unusual deletion of trnK appears in B. hippolytes, reported for the first time in isopods. The area between trnR and trnH was considered a 'hot spot' of mitochondrial gene rearrangement in Isopoda [32,33]. As the most parsimonious explanation of gene order change in this region, Kilpert and Podsiadlowski [25] assumed multiple translocation events. Because of the mixture of inversions and genome shuffling, tandemduplication/random loss models were not a better way to explain the gene rearrangement [25,34]. Crustacean taxa usually exhibit negative overall GC skews and positive AT skews on the heavy strand, but many studies have found that the bias for this strand is inversed in isopod mitogenomes [32][33][34][35]. B. hippolytes mitogenome in this study also exhibits negative AT skews and positive GC skews. It is considered to be the result of the architectural hypervariability and frequent inversions of the origin of mitochondrial replication (ORI) located in the control region (CR), where the changed replication order of two mitochondrial DNA strands consequently resulted in an inversed strand asymmetry [33,35].
Comparison of mitochondrial gene order across Isopoda shows relative stability within suborders, but substantial differences among suborders (Fig. 4). In Epicaridea, we speculated three small conserved syntenic blocks of the mitogenomes. Both deletion and duplication of genes occurred in this suborder, particularly Argeia pugettensis. Whereas in Cymothoida and Oniscidea, we supposed two large conserved blocks. Comparison of conserved regions of different suborders can indicate that mitogenomes gene order is an especially useful tool for higher taxonomy of isopod species. Sequencing of additional isopod mitogenomes promises to be a fertile ground for research.

Conclusions
Our phylogenetic analyses based on cox1, 18S sequence and mitochondrial genome suggest that Bopyrinae is polyphyletic and Bopyroides hippolytes should be excluded from the Bopyrinae. We found a novel gene order in B. hippolytes compared to other isopods. The comparison of mitochondrial gene order shows that conserved syntenic blocks have distinctive characteristic at a subordinal level, and may be helpful for understanding the higher taxonomic level relationships of Isopoda.

Taxon sampling
Bopyroides hippolytes (Fig. 4)  Voucher specimens are deposited in the Florida Museum, University of Florida (UF).

DNA extraction, amplification, sequencing, and annotation
Total genomic DNA was extracted from eggs or the pereon of female specimens using the genomic DNA rapid extraction kit (Aidlab Biotechnologies Co., Ltd) according to the manufacturer's instructions. 18S rRNA gene region was amplified and sequenced using the same primers from Boyko et al. [3]. PCR conditions were: denaturation at 98 °C for 2 min, 40 cycles of 98 °C for 10 s, 50 °C for 15 s, and 68 °C for 1 min, and a final extension of 72 °C for 10 min. The complete mitochondrial genome was amplified by PCR using 12 pairs of primers (Table  S5). The mitochondrial genome was sequenced and annotated following our previous study [4].
The complete mitogenome was obtained from Bopyroides hippolytes (GenBank accession number:

Gene arrangement comparisons
PhyloSuite [35] was used to batch-download the 11 complete isopod mitochondrial genomes available from Gen-Bank, and to assess genomic features and gene order. Phylograms and gene orders were visualized in iTOL [36]. For the purposes of visualization, we arbitrarily designated the beginning of the cox1 gene as position 1 in each genome (pointing in the direction of cox2).

Phylogenetic analyses
Most of the available sequence data for Bopyridae are mitochondrial genes and 18S rRNA. We constructed four datasets to assess phylogenetic relationships of bopyrids (Tables S2, S3, S4): (1) cox1 dataset; (2) 18 s rRNA dataset; (3) combined cox1 + 18 s rRNA dataset; (4) mitogenome dataset. Nucleotide sequences were used in the first three datasets and amino acid sequences were used in the mitogenome analyses. Missing sequences were treated as missing data. MAFFT [37,38] was used to align sequences: nucleotide and amino acid sequences were aligned in batches (using codon and normalalignment modes, respectively) with "-auto" strategy, whereas 18S rRNA gene was aligned using Q-INS-i algorithm, which takes secondary structure information into account. We used Gblocks v0.91b (http:// molev ol. cmima. csic. es/ castr esana/ Gbloc ks_ server. html) [39] to eliminate the ambiguous sequences after alignment, as they impact phylogenetic analyses [40]. Parameters were set as follows: type of sequence was set to codons for PCGs alignments; the minimum length of a block was set to 3 bp for PCGs and 2 bp for rRNA genes, and gap positions (−b5) were allowed with half. DNAsp v5 [41] and MEGA 5.03 [42] were used to calculate sequence composition and variability.
We tested the performance of homogeneous substitution models (using Maximum Likelihood (ML), and Bayesian Inference (BI)) for the single-gene dataset and COI-18S dataset. Data were partitioned by gene and codon position for the BI and ML analyses. PartitionFinder ver. 2.1.1 [43] was used to select partition schemes for BI, and ModelFinder [44] was used to select these for ML analysis in IQ-TREE, using the corrected Akaike Information Criterion (AICc). MrBayes ver. 3.2.6 [45] was used for BI, with four simultaneous runs with four chains each run for 10 million generations, sampling every 1000 trees. The first 25% of these trees were discarded as burn-in when computing the consensus tree (50% majority rule). Sufficient mixing of the chains was considered to be reached when the average standard deviation of split frequencies was below 0.01. ML analyses were conducted in IQ-TREE [46] with 1000 ultrafast BS [47].
Bayesian Inference with the CAT-GTR model based on mitogenome amino acid (AA) data outperform partitioned homogeneous models in isopods [35]. Consequently, we tested the performance of heterogeneous CAT-GTR model for mitogenome dataset, using AA sequence of 13 protein-coding genes (13PCGs). The CAT-GTR inference was implemented in PhyloBayes-MPI 1.7a on the beta version of the Cipres server [48], with default parameters (burnin = 500, invariable sites automatically removed from the alignment, two MCMC chains), and the analysis was stopped when the conditions indicated that a good run was reached (Phy-loBayes manual: maxdiff < 0.1 and minimum effective size > 300).