Complete Mitochondrial Genome Sequencing of the Calanoid Copepod Bestiolina similis (Sewell, 1914) from the Chennai Coast

1. Introduction

Copepods, the most abundant metazoans found on Earth, represent a diverse group of organisms consisting of thousands of species distributed worldwide [1,2]. Among them, the order Calanoida stands out as a particularly rich and significant component of marine ecosystems, encompassing a large number of species that exhibit subtle variations in their morphological characteristics [3,4]. Accurate identification of these species poses a challenge due to their cryptic nature, requiring careful examination and differentiation. To overcome this challenge, innovative approaches incorporating molecular markers have been developed, enabling more reliable and precise species identification [5].

Mitochondrial genes have proven to be highly variable and valuable tools in revealing taxonomically significant variation among copepod species. This variation enables clear discrimination between closely related species and facilitates the recognition of genetic divergence within conspecific groups, often associated with geographical isolation and cryptic speciation [6]. Consequently, mitochondrial (mt) genes have gained significant attention in recent years for investigating the taxonomic position of copepods. Among these genes, mitochondrial cytochrome oxidase I (mtCoI) has emerged as a common marker for rapid species identification of copepods [7]. However, while mtCoI is useful for species identification, it does not provide comprehensive insights into the phylogeny of copepod species, as it is known to be a fast-evolving gene [8]. In contrast, recent studies utilizing complete mitochondrial genome-based phylogenetic analyses have generated considerable interest in exploring the deep branching patterns of arthropods. Mitochondrial genome-based phylogenies have frequently revealed intriguing and sometimes striking topologies compared to traditional phylogenetic approaches [9,10].

The family Paracalanidae comprises dominant calanoid copepods in tropical oceans, playing a crucial role in the grazing behavior of ocean ecosystems [11]. Among them, Bestiolina similis is a prominent calanoid species found abundant in some parts of coast of India [12]. It belongs to the genus Bestolina (Andronov, 1991) within the Paracalanidae family, consisting of 11 species [13]. Bestiolina is a new genus which was previously considered along with genus Acrocalanus. However, Acrocalanus species are characterized by their absence of a fifth leg in female [14]. On the contrary, Bestiolina species has a small rudimentary fifth leg in female. In the past two decades, a significant number of Bestiolina species have been documented [15,16,17,18]. However, their small size (>1mm), inadequate sampling techniques (using mesh size above 200 µm) and misidentification with copepodite stages of other paracalanid species makes them challenging to study [19]. Due to the following reasons, accurate taxonomic identification can be challenging, prompting the use of molecular approaches for taxonomy. Therefore, the objective of this study is to extract and sequence the complete mitochondrial genome of the marine copepod B. similis for the first time and deposit it in GenBank for further investigations. Currently, only a limited number of copepod mitochondrial genomes have been sequenced, with only a few belonging to the calanoid group. Hence, this study will construct a phylogenetic tree based on the extracted mitochondrial genome, incorporating the available data, to enhance our understanding of the evolutionary relationships among calanoid copepods. This research contributes to the expanding knowledge of copepod biodiversity and provides valuable insights into their taxonomic and phylogenetic relationships in marine ecosystems.

2. Materials and Methods

2.1. Copepods collection and identification

Zooplankton samples were collected from the Marina beach (13°03'00'' N, 80°16'56.64'' E, India) using a standard Bolton silk zooplankton net with a mesh size of 150μ five nautical miles away from the coast during the month of August, 2022. The collected samples were preserved in 95% ethanol and transported to laboratory. Copepods were identified up to species level by carefully dissecting and studying the key characters [20,21].

2.2. DNA Extraction and quality Control

DNA from copepods was extracted using the Qiagen DNeasy Blood and Tissue Kit (Cat No.69506). The copepods were pelleted through centrifugation, and the supernatant was carefully removed, leaving behind the pellet, which was then air-dried. Subsequently, each microcentrifuge tube containing the pellet received ATL buffer and Proteinase K, followed by a brief vortexing and incubation at 56⁰C for 2 hours. RNase A treatment (MP Biomedicals; Cat. No.210107683) was carried out at 65⁰C for 20 minutes. The resulting lysate was thoroughly mixed with half the volume of absolute alcohol and loaded into a DNeasy mini spin column, placed within a 2 ml collection tube. The tubes were centrifuged at 8000 rpm for 1 minute, and the flow-through was discarded. The remaining column wash steps were conducted as per the manufacturer's instructions. Finally, DNA was eluted from the column using 10 mM Tris HCl, pH 8.0.

2.3. Library Preparation for mitochondrial genome

The DNA library preparation was carried out following the QIASeq FX DNA Library Preparation protocol (Cat#180475) in accordance with the manufacturer's instructions. The process involved enzymatic fragmentation, end-repair, and A-tailing of 1 ng to 10 ng of Qubit quantified DNA in a one-tube reaction using the FX Enzyme Mix provided in the QIASeq FX DNA kit. Subsequently, the end-repaired and adenylated fragments underwent adapter ligation, where an index-incorporated Illumina adapter was ligated to generate the sequencing library. The library was then subjected to 10 to 12 cycles of Indexing-PCR, including initial denaturation at 98˚C for 2 minutes, cycling (98˚C for 20 seconds, 60˚C for 30 seconds, 72˚C for 30 seconds), and final extension at 72˚C for 1 minute, to enrich the adapter-tagged fragments. Lastly, the amplified library was purified using Sera-MagTM Select beads (Cytiva, # 29343057), followed by a thorough quality control check of the library.

2.4. Illumina sequencing and data analysis

Illumina-compatible sequencing library was quantified by Qubit fluorometer (Thermo Fisher Scientific, MA, USA) and its fragment size distribution was analyzed on Agilent 2200 TapeStation. The libraries were subsequently subjected to paired-end sequencing on an Illumina NovoSeq 6000 sequencer for 150 cycles, following the manufacturer's instructions (Illumina, San Diego, USA). The sequencing generated approximately 26.8 to 34.3 million Illumina short reads data. To ensure data quality, the raw reads were trimmed using Trimgalore-v0.4.01 to remove adapter sequences and low-quality reads. The processed high-quality reads were then utilized for assembly, gene prediction, and annotation using the MitoZ-v3.4 tool, a specialized toolkit for animal mitochondrial genome assembly, annotation, and visualization.

2.5. Phylogenetic tree construction

The phylogenetic trees of mtDNA were constructed using the maximum likelihood (ML) method through MEGA-X software [22]. Existing mtDNA data from the NCBI database of relevant species were utilized for comparison and to produce the phylogenetic tree. Only data representing the genus level were considered for the analysis. Optimal trees were generated through heuristic searches with 1000 bootstraps [23]. The evolutionary history was inferred from the constructed trees using the Jukes-Cantor model [24]. A bootstrap consensus tree was obtained from 1000 replicates, representing the evolutionary history of the analyzed taxa [25]. Branches corresponding to partitions reproduced in less than 50% of bootstrap replicates were collapsed. The initial tree(s) for the heuristic search were automatically obtained using Neighbor-Joining and BioNJ algorithms applied to a matrix of pairwise distances estimated through the Maximum Composite Likelihood (MCL) approach. The topology with the most favorable log likelihood value was selected. Species distances were compared at the phylogenetic levels of calanoid copepods to assess their evolutionary relationships.

3. Results

3.1. Overview of the mitochondrial genome of B. similis

The complete mitochondrial genome sequence of B. similis found to be 23,704 base pair in length, including 13 PCGs and 20 tRNA genes (Figure 1)(Genbank: SUB12990372). It comprises a two longest non-coding region from the beginning of the genome (0 bp to ∼6.5k bp) to the end (∼21.9k bp to 23.7k bp). Other non-coding regions were below 1000 bp in length. Overall nucleotide composition of the complete mitochondrial genome is of 43.42% Adenine, 43.22% Thymine, 6.48% Guanine and 6.88% Cytosine.

3.2. Protein coding regions of the genome

The mtDNA of B. similis revealed 13 protein-coding genes. Table 1 shows the list of protein coding genes observed from the mt DNA. ND5, COX1, ND4 and CYB were observed to be the bigger genes in the genome with 1.713, 1572, 1299 and 1137 bp. ATP8 was observed to be the smallest gene with 156 bp. There was a big non-coding region with ∼1k bp observed between ATP8 and COX1.

3.3. Regions of tRNA genes

The mtDNA revealed 20 tRNAs coding genes each with ∼60 bp. There was ∼2k bp gap observed between the trnD(guc) and trnR(ucg) which was occupied by COX2 and ND4. Two gene coding for tRNA-Leu were observed - trnL(uaa) and trnL(uag). Table 2 gives an account of tRNA genes of the mtDNA of B. similis.

3.4. Phylogenetic tree based on mitochondrial genome of B. similis

The phylogenetic tree, constructed based on complete mtDNA, reveals a close relationship between the observed Bestiolina species and Paracyclopina and Calanus species (Figure 2). However, other calanoid species, such as Labidocera and Temora, appear to be distantly grouped from the Bestiolina and Calanus species. This finding indicates a substantial disorientation within the order Calanoida when considering the families based on complete mtDNA.

4. Discussion

The genus Bestiolina (Andronov, 1991), previously known as Bestiola [26], was distinguished from the genus Acrocalanus by the presence of 6 setae in P3 endopodite, which was 7 in the latter [27]. Since then, taxonomic and phylogenetic investigations on Bestiolina species have garnered considerable interest. Cornils and Blanco-Bercial [27] explored the evolution of Paracalanidae species through both morphological and molecular datasets, revealing that Bestiolina species displayed a close relationship with Acrocalanus species in molecular datasets. However, their relationship differed significantly in morphological datasets.

Species belonging to the genus Bestiolina are commonly found in tropical regions across different oceans, with some speculation about their origin in the Indo-Malayan region [14]. Previous studies by Sivakumar et al. [28], Nawaz et al. [4] Umer et al. [29], Dishad Begum [30], Shanthi and Ranamibai [31] and Rajthilak et al. [32] investigated the copepod diversity in and around the Chennai coast of India; however, none of them reported the occurrence of Bestiolina species. Therefore, this study marks the first report of Bestiolina species in Chennai waters. Nawaz et al. [10] observed a overall abundance of Paracalanidae copepods in Chennai waters. Despite the dominance of paracalanid copepods, Bestiolina species were not observed throughout the study. However, B. similis has been documented as dominant in other parts of the Indian coasts [11,33]. The occurrence of Bestiolina species in the present study suggests their potential migratory abilities.

The size range of metazoan mitochondrial (mt) genomes varies from 14 to 48 kb [34], encompassing genes encoding at least 13 proteins, 22 transfer RNAs, and 2 ribosomal RNAs (rRNA). Notably, the mitogenome typically includes a substantial non-coding region that plays a role in transcription initiation and gene replication [35]. In recent times, there has been an increasing trend in complete mitochondrial genome extraction for species identification and the study of phylogenetic and molecular evolution among closely related species [36]. The protein-coding genes (PCGs) of B. similis exhibit the same genes present in other calanoid mitochondrial genes, such as C. hyperbolus [37] and E. affinis [38]. Unlike Penella and Paracyclopina species, B. similis does possess the ATP8 gene [39], indicating a genomic structure similar to other calanoid species. However, the presence of rRNA genes was not observed in the present study. While vertebrate mitogenomes typically contain 2 rRNA genes, this feature is also common among most invertebrates. However, some cnidarins have been noted to exhibit reduced RNA genes in their genome [40], which may also apply to B. similis.

The phylogenetic tree constructed based on whole genome sequences revealed that E. affinis and L. rotunda were grouped together, forming a distant cluster from B. similis and C. simillimus. Consequently, B. similis exhibited a close proximity to P. nana, which belongs to the order Cyclopoida. This grouping indicates a lack of closely related species from both calanoid and cyclopoid orders in the phylogenetic construction. On the other hand, the positioning of Eutemora and Labidocera far from calanoid species suggests different evolutionary patterns within these calanoid species. Notably, only species from 4 genera of calanoid mitochondrial genomes have been extracted and studied thus far (including this study), which accounts for the displacement of these species in the tree. Furthermore, P. nana is the only cyclopoid species for which the complete mitochondrial genome has been extracted. The difficulty in extracting complete mitochondrial genomes is attributed to the challenges in long PCR amplification. B. similis was observed to possess <80% of A-T content, which may contribute to the possibility of long PCR amplification issues, as indicated by Ki et al [41]. . However, given the limited availability of mitogenomes, it is currently not feasible to draw conclusive phylogenetic relationships between the species. Therefore, future studies are needed to further investigate the molecular evolution of calanoid species based on whole mitogenomes.