Samples of Capitellidae collected from China (Shandong, Hainan, Fujian, Guangxi, and Hong Kong) are shown in Fig. 1, Table 1. All specimens were fixed and preserved in 95% ethanol. Identification was carried out following Magalhães and Blake (2022) by the number of thoracic chaetigers, the hooks and capillaries arrangement in the thoracic and abdominal regions, and the presence or absence of achaetous segment or branchiae.

Figure 1. 

Photos of capitellids newly sequenced in the present study. A. Barantolla sp.; B. Mediomastus sp.; C. Notomastus sp. A; D. Notomastus sp. B; E. Notodasus sp. A; F. Notodasus sp. B; G. Notodasus sp. C. Figures C, F were taken by Su. X.; the others were by Yang D.

In this study, the mtgenomes and nuclear genes (ncgenes, 18S, 28S, and H3) of eight species of Capitellids were obtained with NGS (Table 1). Several parapodia were used to extract the total DNA using the CTAB method. The DNA was sent to Novogene Co. Ltd. (Tianjin, China) or the Centre for Genomic Sciences (HKU, Hong Kong, China) for pair-end sequencing on the Illumina HiSeq platform with a read length of 100 bp or 150 bp (Suppl. material 1: table S3). Raw data quality was evaluated using FastQC v0.11.9 under default settings (Brown et al. 2017). Sequence assembly was performed using CLC Genomics Workbench version 7.0.3 (CLCbio, Arhus, Denmark) or NOVOPlasty v2.7.0 (Dierckxsens et al. 2017), SPAdes-3.15.4 (Prjibelski et al. 2020), and GetOrganelle (Jin et al. 2020). The fragments of the mtgenome of Capitella teleta were also obtained by PCR and sequenced by first-generation sequencing (FGS). The primers used in PCR are listed in Suppl. material 1: table S4. The fragments were assembled into one contig by CAP3 (Huang and Madan 1999).

The contig containing the mtgenome or ncgenes was selected using BLAST V2.13.0 (Camacho and Madden 2013). The mtgenome was annotated using MITOS2 (Bernt et al. 2013b). Since MITOS2 failed to detect the atp8 in Notomastus sp. B, the hmmsearch was used to detect putative atp8 following Sun et al. (2021). In detail, all possible open reading frames (ORFs) in the forward and reverse strands were predicted using ORFfinder (https://www.ncbi.nlm.nih.gov/orffinder/). A comprehensive database containing atp8 amino acid sequences from 249 species of Annelida (Suppl. material 1: table S5) was constructed by hmmbuild via HMMER v3.2.1 (Potter et al. 2018). Subsequently, the predicted ORFs (amino acids) were annotated by comparing them with the database using hmmsearch by turning all heuristic filters off (-max) to accept low-scoring matches as positive ones. The boundaries and the start or stop codons of PCGs were manually checked by sequence alignment with that of closely related species. All sequences of the mtgenomes and ncgenes obtained in this study have been deposited in GenBank (Table 1). The GC-skew and AT-skew were calculated according to Perna and Kocher (1995).

The datasets containing different genes or taxa were used for phylogenetic analysis. A total of 65 taxa were included in the five-gene dataset (i.e., 16S, 18S, 28S, H3, and cox1) (Suppl. material 1: table S6). For nine capitellids with the mtgenome, i.e., eight new in this study and Notomastus sp. LC661358.1 (Kobayashi et al. 2022b), four different datasets of nucleotide sequences were used to determine their evolutionary relationships: 37 mtgenes, 13 PCGs, ncgenes (i.e., 18S, 28S, and H3), and 37 mtgenes + ncgenes. Opheliidae and Echiuroidea were used as outgroups with Opheliidae as the root. Each gene was aligned using the MAFFT version 7 under default settings (Katoh et al. 2017). Gene concatenation was managed in PhyloSuite (Zhang et al. 2020). The Gblocks Server was used to remove poorly aligned positions or regions by less stringent selection settings (Talavera and Castresana 2007). The best model was detected by PartitionFinder2-2.1.1 in nucleotide mode (Lanfear et al. 2017), and the predicted best model has been listed in Suppl. material 1: table S7. The phylogenetic tree was constructed using both Maximum Likelihood (ML) and Bayesian Inference (BI) methods in PhyloSuite (Zhang et al. 2020), with ML analysis performed using IQ-TREE (Nguyen et al. 2015) with 10,000 standard bootstrap replicates, and BI analysis conducted using MrBayes (Ronquist et al. 2012) with a sampling frequency of 1,000 replicates.

Pairwise comparisons of mtgene order patterns of both capitellid and outgroups’ mtgenomes were performed using CREx2 (Bernt et al. 2007), and the results were used to assess mtgene rearrangements, including transposition (T), reverse transposition (RT), reversal (R), and tandem-duplication-random-loss (TDRL). Because of the high variability of tRNAs, only PCGs and rRNAs were used to compare the gene arrangement and deduce the hypothetical gene order of their common ancestor.

Base substitutions were analyzed in 13 PCGs of capitellid mtgenomes. Nucleotide sequences were aligned based on their corresponding amino acid sequences using the MUSCLE module (Edgar 2004) in MEGA11 (Tamura et al. 2021) to ensure accurate codon alignment for the Ka/Ks calculation. The aligned sequences were then used to calculate the Ka (non-synonymous substitutions) to Ks (synonymous substitutions) ratio (Ka/Ks) using KaKs_Calculator 3.0 (Zhang 2022). Ka/Ks > 1 indicates genes undergoing positive selection, Ka/Ks = 1 for neutral evolution, and Ka/Ks < 1 for purifying selection (Yang and Nielsen 2000). Nucleotide diversity of the 13 PCGs and two rRNAs from Capitellidae were detected by DnaSP 6.12.03 (Librado and Rozas 2009) and the gggenes (Wilkins 2023) for visualization.

Raw data first reported in this study are available in the NCBI with BioProject No. PRJNA1232974. BioSample accessions are SAMN47255958-SAMN47255963. The raw data of Notomastus sp. A and Notodasus sp. A are not included since they are corrupted. Please feel free to contact us if anyone needs their assembled contigs. The GenBank accessions of mtgenome, 18S, 28S, and H3 are shown below: PP133660–PP133665 and PQ010756–PQ010758 for mtgenome; PV277003–PV277010 for 18S-28S contigs; PQ365479–PQ365486 for partial 18S; PQ365487–PQ365494 for partial 28S; PQ631138–PQ631145 for H3.

The mtgenome length of nine capitellids ranged from 14,571 bp to 18,972 bp (Table 1, Suppl. material 1: table S8). Both mtgenomes and PCGs of most of the capitellids are AT-rich (mtgenome: 49.71% to 64.44%; PCGs: 50.05% to 63.96%) (Suppl. material 1: table S8). The AT-skews are -0.282 to 0.081 for mtgenomes and -0.366 to 0.056 for PCGs, respectively; the GC-skews are -0.362 to 0.298 for mtgenomes and -0.422 to 0.298 for PCGs (Suppl. material 1: table S8). For most species, full mtgenomes and PCGs have negative AT-skews and GC-skews, tRNAs have positive AT-skews and GC-skews, and rRNAs have negative GC-skews and positive AT-skews. Only C. teleta has negative AT-skews and positive GC-skews for all mtgenes (Fig. 2). ATG or ATA is the start codon for most PCGs. Both complete (TAA or TAG) and incomplete stop codons (TA- and T-) are present (Suppl. material 1: table S9).

Figure 2. 

Comparison of AT- and GC-skews of the full genomes, protein-coding genes (PCGs), tRNAs, and rRNAs for the mtgenomes of Capitellidae.

All the capitellid mtgenomes contain 13 PCGs, 2 rRNA genes, and 22 tRNAs. The mtgenome of Barantolla sp. has the same gene order as that of Mediomastus sp. except for the location of trnE. There is no conserved gene block for the others (Fig. 3). A different number of long non-coding regions (>100 bp) between genes of the capitellid mtgenome (Fig. 3, Suppl. material 1: table S10): one in Notomastus sp. B and Notodasus sp. B, two in C. teleta and Notomastus sp. A, and five in Notomastus sp. and Notodasus sp. A. Moreover, a long non-coding region (considered as Group II intron) is in cox1 of newly sequenced Notomastus sp. B and Notodasus sp. A & C. Moreover, the putative control region is in different locations of the capitellid mtgenome.

Figure 3. 

Comparison of gene orders of capitellid mtgenomes. Conserved gene clusters in Annelida are indicated by different colors. Black circles represent non-coding regions with > 100 bp sequences between genes. Triangles stand for the long non-coding region (Group II intron) within the cox1. Red circles represent the beginning and end of the contig where the putative control region might be located.

The five-gene dataset was used to reconstruct the phylogenetic tree in Capitellidae (Suppl. material 1: table S6). After alignment, the total length of the five genes was 3,257 bp, containing 453 bp for cox1, 373 bp for 16S, 1,473 bp for 18S, 640 bp for 28S, and 318 bp for H3 (Suppl. material 1: table S11). The Capitellidae clade had a high support value for the ML analysis (93%) and a low support value (#) for the BI analysis. Asterisks stand for the highest support values, and number signs stand for the support value < 50 for ML or < 0.5 for BI or not supported. The same applies hereinafter. Capitellidae primarily contained two clades, i.e., Clade A and B. They are lowly supported by the ML and BI analysis (A, ML/BI: 64/#; B, ML/BI: 54/#). Clade A contains Barantolla, Decamastus, Heteromastus, Leiochrides, Mediomastus, Mastobranchus, and Notomastus; Clade B covers species from all genera included except Mastobranchus. Capitella is monophyletic (Fig. 4). Seven genera (i.e., Barantolla, Decamastus, Heteromastus, Leiochrides, Mediomastus, Notomastus, Notodasus) are paraphyletic or polyphyletic. Notably, the support values are only moderate for most nodes in the tree.

Figure 4. 

Phylogenetic tree of Capitellidae based on concatenated data of 16S, 18S, 28S, H3, and cox1. The tree shown is from ML analysis. The species newly sequenced for this study are in red. Bootstrap support values of ML (left) and posterior probability of BI (right) are indicated above the nodes. Asterisks stand for the highest support values, and number signs stand for the support value < 50 for ML or < 0.5 for BI or not supported.

Relationships of nine capitellids with mtgenome in phylogenetic trees based on 37 mtgenes and 13 PCGs are identical, while the location of C. teleta is different in the tree based on 37 mtgenes + ncgenes (Fig. 5). In all the trees, the clade Capitellidae is well supported by both the ML and BI analyses (ML/BI: */*). Two clades are present within Capitellidae: Mediomastus and Barantolla consistently clustered together with fully supported (ML/BI: */*) and form Clade 1, and species of Notomastus, Notodasus, and Capitella form Clade 2 with a high support value (ML/BI: */* in all the trees) (Fig. 5). Notomastus sp. A + (Notomastus sp. + Notomastus sp. B) form Clade 2a, and (C. teleta + Notodasus sp. A) + (Notodasus sp. B + Notodasus sp. C) form another clade for both 37 mtgenes and 13 PCGs (Fig. 5A). When nuclear genes are included, the location of C. teleta is different, and C. teleta + (Notodasus spp. + Notomastus spp.) from Clade 2 (Fig. 5B). Moreover, the support values of clades with C. teleta are low in the 37-mtgenes tree, and the support values increase when nuclear genes are added.

Figure 5. 

The phylogenetic tree of Capitellidae based on 37 mtgenes, 13 PCGs (A), or 37 mtgenes + ncgenes (18S + 28S + H3) (B). The tree shown is from ML analysis based on 37 mtgenes + ncgenes (18S + 28S + H3) and 37 mtgenes. Asterisks stand for the highest support values, and number signs stand for the support value < 50 for ML or < 0.5 for BI or not supported.

The evolutionary patterns of PCGs and rRNA are shown in Fig. 6. A side-by-side comparison of the phylogenetic tree and gene order by TreeREx revealed the evolutionary pattern of PCGs and rRNAs (Fig. 6). The gene cluster cox1-cox2-atp8 is present in the hypothetical ancestor of Capitellidae & Urechis, but not in Capitellidae only. One transposition (and5-nad4l) happened from the capitellid ancestor to the common ancestor of Mediomastus and Barantolla, and three transpositions (nad6, cox3, nad2) from the capitellid ancestor to the common ancestor of Capitella, Notomastus, and Notodasus (Fig. 6, Suppl. material 1: fig. S1). A particularly complex process was detected from the common ancestor of Notomastus and Notodasus to the present species, and more than half of the mtgenes are involved in the rearrangement scenario (Fig. 6, Suppl. material 1: fig. S1).

Figure 6. 

Possible mechanisms of mtgene rearrangement in Capitellidae estimated using TreeREx. Conserved gene clusters in Annelida are indicated by different colors. Only PCGs and rRNAs are included. R, reversal; T, transposition; RT, reverse transposition; TDRL, tandem-duplication-random-loss.

The Ka/Ks of the 13 PCGs of Capitellidae ranged from 0.031 (cox1) to 0.385 (atp8), suggesting that they have undergone strong purification selection (Fig. 7A). When comparing the Ka/Ks of PCGs in different genera, the relatively conserved Mediomastus & Barantolla revealed much lower Ka/Ks values for each gene than Notodasus or Notomastus (Fig. 7B). Nucleotide diversity (Pi) analyses of 13 PCGs and 2 rRNAs from capitellid mtgenomes showed that the average nucleotide diversity of Barantolla & Mediomastus was lower than that of Notodasus and Notomastus (Fig. 7C). The Pi values of rrnS and atp6 are low (< 0.1). The Pi values of nad5, nad1, and nad4 are high (>0.5) (Fig. 7C).

Figure 7. 

Base substitution rate of the 13 mitochondrial PCGs in Capitellidae. A. Ka/Ks for each mtgene; B. Ka/Ks of Notodasus and Notomastus and Mediomastus & Barantolla. Ka/Ks < 1 indicates purification selection. (C) Nucleotide diversity analysis of 13 PCGs + two rRNAs based on Capitellidae. The yellow line represents Barantolla & Mediomastus, the red line represents Notodasus, and the blue line represents Notomastus.

Traditionally, mtgenomes in marine annelids were thought to be compact. However, recent studies have identified long non-coding regions in Siboglinidae (Li et al. 2015), Owenia fusiformis (Weigert et al. 2016), Syllidae (Aguado et al. 2016), Serpulidae (Seixas et al. 2017; Sun et al. 2021), Decemunciger (Bernardino et al. 2017), certain deep-sea species of the Polynoidae (Zhang et al. 2018; Hiley et al. 2024), Armandia (Kobayashi et al. 2022b), and Synelmis (Huč et al. 2024). In this study, the gene orders in Barantolla and Mediomastus are identical (except trnE), but there is substantial variability in Capitella, Notomastus, and Notodasus (Fig. 3). In contrast, long non-coding regions (>100 bp) were not found in Barantolla and Mediomastus and were found in species of Capitella, Notomastus spp., and Notodasus spp. (excluding Notodasus sp. C), suggesting that these regions might be involved in mtgene rearrangement and have once been present and removed afterward in Notodasus sp. C. This implies that different strategies or mechanisms may govern the rearrangement of capitellid mtgenomes.

The increased substitution rate may lead to a high rate of gene rearrangement (Shao et al. 2003). Our results indicate that the substitution rate of capitellid mtgenes varies. In species with relatively conserved gene orders, such as Barantolla and Mediomastus, the substitution rates for each PCG are significantly lower than those in Notomastus and Notodasus (Fig. 7B). This suggests that greater variability in mtgene order may correspond to higher substitution rates, as previous studies showed (Bernt et al. 2013a; Shao et al. 2003).

The gene order in the nine capitellid mtgenomes differs significantly from that of other annelids, and most conserved gene clusters in Annelida were split in capitellids (Jennings and Halanych 2005; Struck et al. 2023). Also, they exhibit intragenus variability. Previous studies showed that Echiura is the sister taxon to Capitellidae (Tomioka et al. 2018; Kobayashi 2023). The conserved gene cluster cox1-cox2-atp8 is present in both Urechis spp. and Notomastus sp. (Tomioka et al. 2018; Kobayashi et al. 2022b), which supports their close relationship, but the similarity of mtgene order between Echiura and Capitellidae was concluded on a limited sampling size. In our results, cox1-cox2-atp8 is not present in the hypothetical ancestral mtgene order of Capitellidae (Fig. 6), indicating the evolution of mtgene order in Capitellidae is more complex than previously thought. Further work needs to be done on more taxa covering more genera to more accurately predict the ancestral mtgene order in Capitellidae.

In recent years, phylogenetic analyses of the Capitellidae, based on limited molecular data and taxa, have produced conflicting results (Rousset et al. 2007; Goto et al. 2013). The most comprehensive phylogenetic analysis of Capitellidae utilized a concatenated gene marker dataset, including 18S, 28S, H3, and cox1 (Tomioka et al. 2018). Our findings align with this earlier work, confirming that Capitella is monophyletic while most of the other examined genera are paraphyletic or polyphyletic. Moreover, the topology of the tree changed when different outgroups were included for the five-gene tree (Fig. 4, Suppl. material 1: fig. S2). This discrepancy may arise from the under-sampling of taxonomic units or the use of limited genetic markers since there were only one or two markers for nearly half of the taxa (40% in this study for Fig. 4, Suppl. material 1: table S6).

Additionally, the accuracy of species identification by morphology is crucial, as this will affect our interpretation of molecular phylogenetic analyses. Due to its low mutation rates, the 18S is not ideal for distinguishing closely related species in Capitellidae, particularly for Dasybranchus: there is only partial 18S for D. caducus and partial 18S and 28S for Dasybranchus sp. The 18S of both D. caducus and Dasybranchus sp. have a high sequence identity compared to those of three Notodasus spp. in this study, 99.1–100% and 98.9–99.8%, respectively (Suppl. material 1: table S12). This high level of conservation in 18S suggests that this gene lacks the resolution for species-level delimitation. In contrast, for 28S, Dasybranchus sp. showed 99.6% identity with Notodasus sp. C and 93% identity with Notodasus sp. B (Suppl. material 1: table S12). More evidence is needed to determine whether the D. caducus in Fig. 4 is the same species as Notodasus sp. C.

Mtgenomes are increasingly employed in phylogenetic studies due to the wealth of genotypic information they provide (Boore 1999). They have a higher likelihood of accurately tracking short internodes compared to ncgenes (Moore 1995). The tree topology of Capitellidae presented in this paper, based on datasets that include ncgenes, is different from trees constructed with mtgenes only (Fig. 5). This indicates that mtgenomes and ncgenes might have different evolution trajectories, and they together could contribute to elucidating the complex relationships within Capitellidae. Our results with the combined mtgenes and ncgenes showed that Barantolla is the sister group to Mediomastus, with Notomastus and Notodasus being monophyletic. Notably, the phylogenetic position of Capitella is different in the tree based on ncgenes, mtgenes, and mtgenes + ncgenes (Fig. 5, Suppl. material 1: fig. S3). It is important to note that the phylogenetic placement of Capitella remains unclear. When mitochondrial rRNA and tRNA genes are included in the phylogenetic dataset, the clade that includes Capitella generally shows lower support values (Fig. 5A). In trees constructed using ncgenes only, Capitella appears as the sister taxon to all the other capitellids (Suppl. material 1: fig. S3). Notably, their phylogenetic tree support values were not high (59/0.84 for ncgenes; 69/0.76 for 18S-28S), and Capitella might be misplaced (Suppl. material 1: fig. S3). This might be due to the difficulty in the alignment of the ribosomal genes. The tree topology might be greatly affected by long-branch attraction. Phylogenetic trees based on nad2, nad4, and nad6 exhibit the same topology as that based on the mtgenome, while atp6, cob, nad5, and rrnL exhibit the same topology as that based on the mtgenes + ncgenes (Suppl. material 1: fig. S4). When Capitella was not considered, most trees based on a single gene supported a close relationship between Barantolla and Mediomastus and the monophyly of Notomastus and Notodasus. It is worth noting that our results are based on a limited taxon sampling when considering there are 41 genera in Capitellidae. Incorporating more taxa into future studies is essential to better elucidate the relationships among these genera.

Since the first report of Group II introns in Nephtys (Vallès et al. 2008), they have been found in different taxa of polychaete, such as Glyceridae (Richter et al. 2015), Decemunciger (Bernardino et al. 2017), Travisiidae (Kobayashi et al. 2022a), Notomastus and Armandia (Kobayashi et al. 2022b), and Chaetopteridae (Wu et al. 2024). The length of Group II intron varies in different species: 1,819 in Nephtys sp. (Vallès et al. 2008), 1,647 in Decemunciger sp. (Bernardino et al. 2017), 1,646 in Mesochaetopterus tingkokensis, and 365 bp in Phyllochaetopterus sp. (Wu et al. 2024). Group II intron in the cox1 gene is in the newly sequenced Notomastus sp. B, Notodasus sp. A & C. Its length in these capitellids ranges from 273 to 447 bp, and its insertion position into the cox1 genes is the same. This result indicates that the Group II intron might be present in the common ancestor of Notomastus and Notodasus. However, answering whether they are all shorter than those of some other polychaete species requires the sequencing of more capitellids.