New insight into the taxonomic resolution of the genus Bythotrephes Leydig (Crustacea: Cladocera) based on molecular data from Central Europe

The taxonomic status of the genus Bythotrephes Leydig (Crustacea: Cladocera) has been debated since the second half of the XIX century. The most widespread view of recent decades has been that Bythotrephes is a monotypic genus, which was support by preliminary molecular data. However, the recent detailed morphological revision of this genus clearly distinguishes at least seven species. Therefore, we performed a multi-lake survey in Central Europe to give new insight into the taxonomic status of Bythotrephes by combining genetic analysis with traditional morphology-based taxonomy. Based on the morphology we identified two species in Central Europe, B. brevimanus and B. lilljeborgi, as well as hybrid forms. For the genetic analysis, we used newly obtained 113 sequences of mtDNA COI gene of the 535-bp length Bythotrephes from Central Europe and sequences downloaded from GenBank. There were no significant differences between all analyzed sequences, which supports the hypothesis that Bythotrephes is a monotypic genus, with only one highly polymorphic species. On the other hand, the results of our work could point out that the COI gene is insufficient to evaluate the taxonomic status of Bythotrephes. Nonetheless, we have identified 29 new haplotypes of mtDNA COI, and one which was the same as the haplotype found in North America and Finland. Furthermore, this haplotype was the source variant from which most other haplotypes were derived.

www.nature.com/scientificreports/ is not widely distributed there and did not penetrate North America from the west 10 . However, recent research indicated the presence of Bythotrephes in lakes sediments in North America dated by at least the early 1900s 16 . The morphological traits are well described for each species, however, the genetic data are fragmentary and ambiguous. Allozyme studies from Russia indicate that the genus consists of several geographical groups, which could be separated by both morphological and genetic characters 17,18 . Based on these results, it was proposed that there were few species in the genus. This point of view was recently accepted 10,11,19 . On the other hand, allozyme studies from North America indicate a low genetic differentiation between B. longimanus and B. cederströmii [20][21][22] . Furthermore, Berg and Garton (1994) found that allozyme variation among sympatric populations of B. longimanus and B. cederströmii was even larger than that among allopatric populations, and they first suggested the existence of only one species, B. longimanus 21 . However, allozyme studies are based on protein sequences, and owing to multiple coding patterns, species-based variation in DNA sequences might not be detected 1 . Therriault et al. (2002) used sequencing of the mitochondrial COI gene to explore the genetic nature of these species. They analyzed B. cederströmii (2 lakes in Canada), B. cf. brevimanus (1 lake in Finland, 2 lakes in Germany) and B. cederströmii connectens from the Volga River in Russia, which seems to be a hybrid form. These results revealed low genetic differentiation between these species and they suggest that Bythotrephes is a monotypic genus with a single species, Bythotrephes longimanus Leydig, 1860. Furthermore, the same haplotype was identified in Lake Ontario and Lake Puruvesi (Finland) suggesting a potential source of invasion from the Baltic Sea region 1 .
Nevertheless, the records of Bythotrephes are scattered and incomplete, especially in Europe where most species or forms occur. Therefore, we performed morphological and genetic analyses of Bythotrephes populations in 20 lakes in Central Europe (NW and NE Poland, Lithuania). Our study aims to present the intra-and interpopulation genetic variability of Bythotrephes in Central Europe. We also present morphological differences between B. brevimanus and B. lilljeborgi. Furthermore, photography documentation of each individual together with its COI sequence was uploaded to the BOLD System. We believed that our molecular data combined with available sequences will allow a better understanding of phylogeographic patterns of this genus and may provide new insight into its history of invasion in North America.

Results
Morphology analysis. Based on the morphological features we identified in Central Europe two species B. brevimanus and B. lilljeborgi (Fig. 1) following keys of Korovchinsky 10,11,23 , and possible hybrid between these species. B. lilljeborgi was found in three lakes in NE Poland (Gaładuś, Białe Wigierskie and Hańcza) and one in Lithuania (Galstas). B. brevimanus was found in three lakes of NW Poland (Drawsko, Żerdno, Lubie) and in Lake Studzieniczne (NE Poland). The possible hybrids between these species were found in other lakes, while in eastern part were more like B. lilljeborgi and in west Poland were more like B. brevimanus.
The morphological differences between both species included few morphological traits, like claws of postabdomen and caudal process, thoracic limbs of first pair (tl I) with their armament, and apical end of antennal branches. The claws of postabdomen and caudal process in B. lilljeborgi are comparatively long (4.8-11.6% of body length) and directed downwards, while in B. brevimanus they are short (2.0-6.9% of body length) and directed backward (Fig. 1B, G). The B. brevimanus have generally shorter tl I with a shorter distal segment (Fig. 1H), while in B. lilljeborgi the tl I and distal segment are longer (Fig. 1C). The apical setea of the second segment of tl I in B. lilljeborgi are reduced, while in the B. brevimanus are better developed (Fig. 1C, H). On the protopodite of tl I, pseudognathobasic process armed lateraly and distaly with two outgrowths with apical setae in both species, but in B. lilljeborgi distal outgrowth with numerous spinules (Fig. 1D, I). The distal segment of upper antennal branches with small apical spine in both species, but in B. lilljeborgi also with denticles ( Fig. 1F, J). The hybrid was characterized with varied morphological traits for both species. Some of them have smaller claws of postabdomen directed downward and backward (like B. brevimanus), but the tl I and swimming antennal like B. lilljeborgi. Others have larger claws of postabdomen directed downward, but tl I was similar to B. brevimanus.
However, morphological data was inconsistent with genetic (COI) data. We observed different morphological traits for populations with the same haplotype (H7). Population from Lake Białe Wigierskie (NE Poland) had morphological traits similar to B. lilljeborgi like comparatively long and massive claws of postabdomen and caudal process inserted densely (Fig. S1B), comparatively longer tl I and distal segment of tl I (Fig. S1C), and smaller body size (Fig. S1A). While the same haplotype (H7) from Lake Słowa (NW Poland) had morphological traits similar to B. brevimanus like smaller claws of postabdomen and caudal process inserted comparatively sparsely (Fig. S1F), shorter tl I and distal segment of tl I (Fig. S1E), and larger body size (Fig. S1D).
Genetic analysis. Our analysis of a 535-bp mitochondrial DNA cytochrome c oxidase subunit 1 (COI) gene fragment amplified from 113 Bythotrephes individuals from 17 studied lakes in Poland and 3 lakes in Lithuania yielded 45 polymorphic sites with a total of 37 transitions and 8 transversions, defining 30 mtDNA COI haplotypes ( Fig. 2; GenBank accession numbers for obtained haplotypes: MZ196458-MZ196487). Twenty-nine mtDNA COI gene haplotypes were detected for the first time, while only one haplotype from our study (H8) was identical as the haplotype previously described by Therriault 1,24 , and two unpublished sequences described as B. longimanus (GenBank acc. num. GU689226 and HQ966461). In total, the 20 studied lakes haplotypes and nucleotide diversity were very high (h = 0.94, π = 1.30%; Table 1) and the number of mtDNA COI haplotypes ranged from 1 (Lake Serwy in NW Poland, and lakes Galstas and Aviris in Lithuania) to 4 (lakes Lubie and Żerdno in NW Poland, and Lake Hańcza in NE Poland) ( Fig. 2; Table 1). The highest genetic diversity was found in the lakes Lubie and Żerdno in NW Poland (h = 1.00, 0.90; π = 1.59, 1.05, respectively), whereas the corresponding values (h and π = 0) were the lowest in lakes Serwy (NE Poland), Galstas and Aviris (Lithuania). Although the number of haplotypes per lake did not www.nature.com/scientificreports/ exceed five, the number of segregating sites within studied lakes ranged from 0 to 15 (Lake Lubie). Values of the mean number of pairwise differences differed from 0 to 8.50 (Lake Lubie; Table 1). The mtDNA COI haplotypes of genus Bythotrephes were used to infer the phylogenetic tree. The ML phylogenetic analyses revealed the presence of the two evolutionary branches: East and Central (Fig. 3), and the net divergence (d A ) between branches was 0.7% (SE ± 0.3%). The majority of haplotypes found in NE Poland created the East branch, together with haplotypes H14 and H22 found in the Lithuania lakes (Figs. 3 and 4), which seems to be B. lilljeborgi, or hybrid B. lilljeborgi × B. brevimanus. Only haplotypes H1, H2, and H3 from this branch were found in two lakes from NW Poland (Lubie and Żerdno). Haplotype H36 (AF435129) recorded in the Volgograd Reservoir by Therriault et al. (2002) was grouped also in the East branch and was described as  In this phylogenetic branch, the position of the H8 haplotype indicates that it was the sequence variant source from which a significant number of haplotypes were derived (Fig. 4). Interestingly, this haplotype, apart from the Polish lakes Ińsko (NW Poland) and Gaładuś (NE Poland), was previously described in Canada Lake Simcoe and Ontario, and in Finland Lake Puruvesi (Therriault et al. 2002; AF435122) and it was identified as B. cederströmii. This phylogenetic branch also grouped haplotypes described by Therriault et al. (2002) in the European part of Bythotrephes range: Finland (H31-AF435123), Germany (H32-AF435124, H33-AF435125, H34-AF435126), Netherlands (H35-AF435128) which seems to be B. brevimanus. Therefore, our Central branch groups with other lakes in Central Europe with the domination of B. brevimanus, as well as with B. cederströmii from Canada and Finland. The East branch seems to be B. lilljeborgi and is grouped also with Volga River which was described as B. cederströmii connectens (Fig. 4). However, the photography provided by Therriault  www.nature.com/scientificreports/ The haplotype network based on the mtDNA COI sequences revealed the presence of the same Bythotrephes phylogenetic branches (Fig. 4). At least 6 mutations distinguished haplotype H8 from the Central branch from haplotypes (H22 and H36) belonging to the East branch. The East branch seems to group B. lilljeborgi from NE Poland and Lithuania, with Volga River (Fig. 4). However, the haplotype 37 from the Volga River is distinguished by at least 8 mutations (Fig. 4).
The most common mtDNA COI haplotypes H1 and H3 occurred with a frequency of 15% and 12.4%, respectively. The H1 haplotypes was found in Lake Lubie in NW Poland and lakes Garbaś, Hańcza, Wigry in NE Poland. In Lake Żerdno from NW Poland and lakes Hańcza, Studzienniczne, Serwy, Gaładuś in NE Poland the haplotype H3 was detected. Thirteen haplotypes (5 from NW Poland and 8 from NE Poland) were found only in one individual (Table S1). While haplotypes H13 and H14 occurred in all analyzed Bythotrephes taken from lakes Aviris and Galstas in Lithuania. Haplotypes 9, 18, 21, 22, 25 and 28 were found in 2-5 individuals, but each of them occurred in only one lake (Ińsko, Garbaś, Busznica, Ancia, Krzemno, and Trzesiecko; respectively) ( Fig. 2; Table S1).
Pairwise genetic differentiation values (mtDNA COI) between studied lakes ranged from 0 to 1.00 (Serwy-Galstas, Serwy-Aviris, Galstas-Aviris). The majority of comparison values of Φ ST were significant (Table S2). The first and second axes of the PCA (PC1 and PC2) performed on the data set containing 19 lakes (Lake Pile was excluded from the analysis due to being represented by only one individual) explained 29.98 and 20.82% of the total variability (Fig. 5). The analysis based on Φ ST -matrix separated almost all NE Poland lakes (except Szelment Wielki and Białe Wigierskie) from lakes localized in NW Poland and Lake Galstas in Lithuania. Studied lakes were in general separated into two groups along PC1 according to their possessing mtDNA COI representing the Central and East phylogenetic branches. The intermediate form was found in lake Lubie, where haplotypes belonging to both ML tree branches were described. With respect to the PC2 the most divergent were NE Poland lakes Garbaś, Wigry, Studzienniczne and Serwy (Fig. 5).

Discussion
The taxonomic status of the genus Bythotrephes has been debated for a long time and now is intensified due to contrary data from genetic and morphological analyses. The morphology data generally indicated the presence of few species within the Bythotrephes genus. By combining molecular data (COI gene) with the morphological description we want to give new insight into the taxonomic status of Bythotrephes. At first in the second half of the XIX century, two species B. longimanus and B. cederströmii were described, and this point of view is the most Table 1. Mitochondrial DNA cytochrome c oxidase subunit 1 (COI) gene (535-bp) diversity indices for 20 Bythotrephes lakes in Poland and Lithuania. N-sample size; Nh-number of haplotypes; h-haplotype diversity; π-nucleotide diversity (%); S-number of segregating sites; PD-mean number of pairwise differences; SE-standard error; *-due to only one individual, Pile Lake was not included in statistical comparisons.    Berg and Garton (1994) found that allozyme variation among sympatric populations of B. longimanus and B. cederströmii was larger than that among allopatric populations, and they also suggested the existence of only one species (B. longimanus) 21 . We also have found large interpopulation genetic diversity (COI), and in one lake even four haplotypes were found. The first results of sequencing COI gene by Therriault     www.nature.com/scientificreports/ North American and Europe and give better proof for Bythotrephes as a monotypic genus 1 . However, their study includes only a few sequences from North America and Europe. Therefore, we performed a multi-lake survey in Central Europe to give new insight into the taxonomic status of Bythotrephes by combining genetic analysis with traditional taxonomy. We identified two species in Central Europe, B. brevimanus and B. lilljeborgi, and also a hybrid between these species. For phylogenetic studies, we used our 113 sequences from Central Europe (Poland, Lithuania) with sequences from Finland, Germany, Netherland (seems to be B. brevimanus), and Canada (seems to be B. cederströmii). There were no significant differences between all analyzed sequences of COI gene, even though we analyzed three species (B. brevimanus, B. lilljeborgi, B. cederströmii) identified based on the morphology. Therefore, our data also support the hypothesis on the monotypic genus of Bythotrephes, with only one highly polymorphic species. The fragment of the mitochondrial COI gene was found to be an informative molecular marker for such studies, and COI data were successfully used for taxon delimitation in cladocerans 31,32 . At the same time, mitochondrial uniformity could be a consequence of hybridization and old introgression. Moreover, COI is a tricky gene-sometimes it has no resolution to discriminate two apparent species (Kotov, 2021, pers. comm., 8 May). The results of our work pointed out that the COI gene is insufficient to evaluate the taxonomic status of Bythotrephes. Most morphological characters are related to the activity of nuclear genes, therefore further analysis of nuclear genes could be more informative for the Bythotrephes genus. Nevertheless, according to morphology, one species seems quite improbable, because morphological differences, at least between some species, are very conspicuous (Korovchinsky 2021, pers. comm., 8 May). The level of gene flow among populations may be a function of dispersal ability 33 . Research suggests that widely dispersing species have high intrapopulation variation and low interpopulation variation 34,35 , while poorly dispersing species show the opposite pattern 33 . Therefore, large intrapopulation variability in our study indicated the wide dispersion abilities of Bythotrephes. Nowadays this species is considered one of the most invasive zooplankton member since they invaded North American 3,4 . Bythotrephes was first observed in North America in Lake Ontario in 1982 36 and was probably introduced via ballast water in ships that travel between St. Petersburg in Russia, and Great Lakes ports 37  We have identified the same haplotype in two Polish lakes, Ińsko (NW Poland) and Gaładuś (NE Poland). Furthermore, this haplotype (H8) in our median-joining network (Fig. 3) was the source variant from which a significant number of haplotypes were derived. The other haplotypes in Central Europe have high similarity to H8, while north-eastern Poland and Lithuania lakes were more similar with Volga River and formed East Branch. Therefore, our results pointed that central Baltic countries (Germany, Poland) could be a more accurate source of the invasion, than the northeast Baltic region. This pattern of gene flow corresponding with major sea routes connecting the Great Lakes with the North Sea and Baltic Sea countries 38 .

Methods
Field work. We collected the Bythotrephes samples in the middle of summer (2018 and 2020) from 20 lakes in Central Europe (NW Poland-8, NE Poland-9, Lithuania-3; Fig. 2). The samples from each lake were collected in triplicate by vertical hauls with a 100-µm plankton net close to the deepest point of lakes or at a location with a depth exceeding 15 m in large and deep lakes. Concentrated samples were added into a 110-ml tube and fixed in 96% ethanol. Individuals of Bythotrephes were separated from other plankton using light microscopy. Photographic documentation was made for each individual for genetic analysis and it is available at BOLDSystem (project code: BYTHO; BIN: BOLD:AAB9254; sample IDs: UwB_1_Bl1, UwB_2_Bl2, …, UwB_113_Bl113). Bythotrephes individuals for genetic analysis were selected and stored in 70% ethanol at 4 °C before DNA extraction. The map of the studied lakes in Central Europe with Bythotrephes haplotypes distribution and frequencies was created in QGIS software (version 2.18.24; www. qgis. org).
Morphological analysis. The species identification was based on the recent review of this genus by Korovchinsky 10,11,18,23 . We present the most important morpho-functional traits which differ between identified species, like claws of postabdomen, thoracic limbs of first pair (tl I) with their ornamentation 10 , pseudognathobasic process on the protopodite of tl I, and apical end of upper antennal branches. Additionally, we performed detailed measurements of body and body parts measurements according to the protocol proposed by Korovichinsky 11,23 for two populations from Lake Drawsko and Galstas. Lake Galstas represent B. cf. lilljeborgi and 'East branch' of our molecular analysis with a single haplotype (H14), which was found only in this lake (Table S1). While Lake Drawsko represents B. cf. brevimanus and 'Central branch' of our molecular analysis with haplotypes H5, H6, and H23. We also compared morphological differences in one haplotype (H7) from Lake Białe Wigierskie (NE Poland) and Lake Słowa (NW Poland). The morphological analysis was performed using an Olympus BX53 microscope with cellSens imaging software.

Molecular analyses.
We analyzed 113 individuals of Bythotrephes from 20 different lakes in Central Europe (Fig. 2). The sample sizes from given lakes ranged from 1 (Lake Pile) to 10 individuals (Lake Hańcza; Table 1). The DNA extraction was performed with the DNeasy Blood and Tissue Kit (Qiagen, Hilden, Germany) according to the manufacturer's instructions. PCR amplification of mitochondrial DNA cytochrome c oxidase subunit 1 (COI) gene was carried out with a GeneAmp PCR System 9700 in 5 μL volumes, and the reaction mixtures consisted of 2 μL of extracted genomic DNA as a template, 1.7 μL of Qiagen Multiplex PCR Master Mix (1x), 0.3 μL mix of newly designed primers (BytCOI_F: 5′-TCG GAA TTT GAG CTG GGA TAGTA-3′; BytCOI_R: 5′-TGT AAG GAG TAT AGT AAT AGC TCC -3′) and 1 μL of Qiagen nuclease-free water. The primers used in this study Statistical analyses. We used DnaSP v.5.10 41 to calculate the measures of sequence variation, including the number (Nh) and frequency of haplotypes, haplotype diversity (h), nucleotide diversity (π, in %), the number of segregating sites (S), and the mean number of pairwise differences (PD). The pairwise Φ ST values for the 19 lakes were calculated in Arlequin v.3.5.1.2. 42 , and their significance was tested using 10,000 permutations that were corrected for multiple tests by Bonferroni correction. Principal coordinate analysis (PCA) was performed on mtDNA COI Φ ST data in Genalex v.6.0 43 . To test the phylogenetic relationships among obtained in the study mtDNA cytochrome oxidase subunit 1 haplotypes, we constructed phylogenetic trees with Mega v.6.06 44 with 1,000 bootstrap replicates used to assess support for tree nodes. The optimal model of substitution 45 for our mtDNA COI sequences was determined using Akaike's information criterion 46 with jModelTest 47 and used to calculate a maximum-likelihood (ML) algorithm. The GTR + I + G model was selected as the best-fitting model for the ML tree, which was additionally rooted using mtDNA COI haplotypes from Cercopagis pengoi (AF320013 from Black Sea and AF320014 from Caspian Sea) downloaded from GenBank. We estimated net pairwise divergence (d A ) between mtDNA COI branches Central and East using the program Mega v.6.06 44