1 INTRODUCTION

The distribution of fish is strictly constrained by the drainage basin because of geographic isolation. In association with the development of river systems and climatic fluctuations, freshwater fish may undergo vicariance, resulting in genetic diversification and even allopatric speciation (Mayr, 1942; Rundle and Nosil, 2005; Tedesco et al., 2012; Yu et al., 2014). Freshwater fish are more susceptive to the physicochemical environment than marine fish, so their distribution patterns are thought to be more related to climate change (Comte et al., 2013). Phylogeography has provided new insights to evaluate the influence of historical events on the genetic structure and composition of modern populations, which has revolutionized population genetics and historical ecology (Avise et al., 1987). Moreover, the phylogenetic study of widespread freshwater fish is an excellent indicator for geographic events.

Gene flow among freshwater fish populations within the same catchment is common, but dispersal among different catchments and drainage basins is typically restricted because of the lack of hydrological connectivity (Mesquita et al., 2005; Huey et al., 2008). However, freshwater fish often inhabit multiple isolated drainage basins, implying dispersal among drainages has been possible in the past (Huey et al., 2011). In China, the colonization of drainages by freshwater fish can be explained by a variety of mechanisms. These include catchment rearrangements through erosion and uplift, climatic-induced fluctuations in sea level (Li, 1981), and anthropogenic translocation (Xiong et al., 2015).

The glacial periods of the Quaternary were a major factor in shaping contemporary distributions of extant flora and fauna (Provan and Bennett, 2008). Fossil and pollen record evidence also suggests that most of the biota, particularly temperate species, survive glacial maxima in lower latitude refugia (Bennett et al., 1991). Many studies have indicated that Quaternary glacial climatic oscillations are too close to explain lots of deep genetic divergences observed nowadays (Huang, 1982; Miralles and Carranza, 2010; Yu et al., 2014). Glacial impacts were greatest in Eurasia, North America, and Antarctica, where glacial ice covered major portions of their surfaces (Hewitt, 2003). China was also affected by the glacial cycle during the Quaternary Period. In the cold periods, mainland China was mainly influenced by alpine glaciers and cold air. During the warm periods, the mainland climate was dominated by warm humid air masses, which led to the retreat of alpine glaciers (Yu et al., 2014).

Cyprinids are the largest extant family of vertebrates and the most diverse freshwater fish family in Asia. The genus Culter consists of nine species and subspecies that are confined to Asia, particularly China, Vietnam, Korea, Mongolia, and Russia (Chen, 1998). Several studies have investigated phylogenetic and biogeography of these Culter species (Wang et al., 2008; Qi et al., 2013, 2015; Yang et al., 2016), however, the sample sizes were small and geographic distances were limited. Therefore, a much larger and more intensive survey is needed to determine the genetic structure of Culter populations. The topmouth Culter (Culter alburnus) and Mongolian redfin (C. mongolicus) are the largest carnivorous fish of the Culter genus (Chen, 1998). They are commercially important freshwater fish that widely coexist in larger rivers, lakes, and reservoirs of China. Their recent use for commercial aquaculture in some provinces of China has seen rapid development (Fan et al., 2008), but wild populations of both species have declined so quickly that conservation of their genetic diversity has become a very important issue.

Our current understanding suggests that extant Culter species originated from the Pliocene epoch in the River Plain of Eastern China (RPEC) according to fossil records (Liu and Su, 1962; Wang, 2013). The Bohai basin is a third-order reservoir constrained by the Miaodao Uplift that was formed by subsidence during the Cenozoic Era (Allen et al., 1997; Hu et al., 2001; He and Wang, 2003). It experienced a threestage evolution from the early Pliocene (3.7 million years ago, Ma) to a marine system in the late Pleistocene (0.13 Ma) (Yi et al., 2016). Fish spread from the RPEC along the Liaohe River and Nenjiang River to Heilong River Basin (LRB) during the Pleistocene. Analysis of the natural distribution of extant fish species in the Zhujiang (Pearl) River Basin (ZRB) suggested that the southward diffusion of East Asian taxa might be caused by the uplift of the Yunnan plateau in the late Pliocene. Hainan Island originally belonged to the ZRB, and was still connected to the South China mainland until the early Quaternary Holocene (10 000 years ago) (Wang et al., 2012). Today, Culter species are widely distributed in Eastern China, especially in the Changjiang River Basin (CRB), Huaihe River Basin (HRB), LRB, and ZRB. We hypothesize that Culter species were subject to drainage basin isolation, so applied the isolation by distance (IBD) model in our present study. Under the population genetics background, this is the process of generating genetic structures through geographically restricted gene flow because random genetic drift is occurring locally (Hardy and Vekemans, 1999).

Here, we employed the mitochondrial DNA control region (CR) to test hypotheses on the drivers of biogeographic process in two sympatric species, C. alburnus and C. mongolicus. Large-scale sampling allowed us to test hypotheses on the evolutionary history and population demographics of these species to determine if Pliocene and Quaternary climatic fluctuations resulted in their genetic diversity in China. We also evaluated the genetic diversity, population genetic structure, and population demographic history of C. alburnus and C. mongolicus and determined whether IBD played a role in shaping the structure of C. alburnus and C. mongolicus. Finally, we assessed the level of genetic diversity across all current populations and proposed suitable conservation plans.

2 MATERIAL AND METHOD 2.1 Samples collection

A total of 333 wild C. alburnus specimens and 225 wild C. mongolicus specimens were collected using three layer gill net from 12 lakes and reservoirs in China between November 2015 and March 2017. Sample information is provided in Table 1 and sampling sites are mapped in Fig. 1. Caudal fin tissue samples used for genomic DNA extraction were preserved at 4℃ and in 100% ethanol and deposited in the Institute of Hydrobiology, Chinese Academy of Sciences.

The Holocene is the latest phase of the geological era when the biological composition was almost identical to that of the modern era (Yan, 2006). In this study, the Songtao Reservoir (STR) on Hainan Island was classed as the ZRB because Hainan Island once belonged to the ZRB, and separated from the South China mainland during the late Quaternary (Zeng and Zeng, 1989).

2.2 DNA extraction, PCR amplification and sequencing

Total genomic DNA was extracted from topmouth Culter and Mongolian redfin tissue using phenolchloroform extraction (Russell and Sambrook, 2001). The mitochondrial CR was amplified using the universal primers MitD1-F: 5'-CACCCYTRRCTCCCAAAGCYA-3', and MitD1-R: 5'-GGTGCGGRKACTTGCATGTRTAA-3' (Huang et al., 2007; Cheng et al., 2011). Each final reaction volume of 40 μL contained 26.4 μL of sterilized ultrapure water, 4 μL 10× Ex-taq buffer, 4 μL of Mg²⁺ (25 mmol/L), 1.6 μL of dNTP-Mix (2.5 mmol/L), 0.4 μL Extaq (5 U/μL), 0.8 μL of each primer (10 mmol/L), and 2 μL of DNA template (approximately 100 ng/μL). PCR amplifications were performed using a T100 Thermal Cycler with a 5 min at 94℃, followed by 5 cycles at 94℃ for 30 s, 58℃ for 30 s, and at 72℃ for 1 min, then 35 cycles at 94℃ for 30 s, 60℃ for 30 s, 72℃ for 1 min, with a final elongation step at 72℃ for 7 min. Amplified PCR products were purified by 1.0% low-melting agarose gel electrophoresis and sequenced in both directions using primers F-SEQ-0510: 5'-TGCATTTGGCTTCTATCTCAGG-3' and R-SEQ-0510: 5'-CCTACTTAAACTCAAGCAAGG-3' and the ABI PRISM® 3730 DNA Analyzer (Applied Biosystems, Life Technologies, Shanghai, China) by the Tsing Ke Commercial Sequencing Service (Wuhan, China).

2.3 Phylogenetic analysis

The generated sequences were aligned against previously published sequences for the same genetic region, and trimmed to the same length using the software package MEGA 7.0 (Kumar et al., 2016).

Phylogenetic relationships based on CR haplotypes were estimated using Bayesian Inference (BI) and Maximum Likelihood (ML) methods. BI analysis was performed using MrBayes 3.2.6 (Ronquist et al., 2012). The best-fit nucleotide substitution models for C. alburnus and C. mongolicus were TrN+G and HKY+I+G, respectively, which was based on the Bayesian Information Criterion in jModelTest 2.1.7 (Darriba et al., 2012). Four simultaneous Markov chain Monte Carlo (MCMC) methods were run for 20 million generations sampling once every 1 000 generations, two independent parallel analyses methods were used. The average standard deviation of split frequencies was less than 0.01 and the potential scale reduction factor value was close to 1. The first 20 000 sampled trees were discarded as burn-in and 50% majority-rule consensus trees was calculated from remaining trees with posterior probability values for each node. Support for clades is expressed as posterior probabilities.

The MLs for phylogenetic analysis were reconstructed in PhyML 3.0 (Guindon and Gascuel, 2003) with 1 000 bootstrap replicates. The most suitable model of DNA substitution, TrN+I+G for C. alburnus and HKY+I+G for C. mongolicus, was determined by ModelTest 3.7 (Posada and Crandall, 1998). Grass carp (Ctenopharyngodon idellus) (GenBank: KJ614563.1) and Hemiculte leucisculus (GenBank: AF494361.1) were set as outgroups.

We constructed haplotype networks using Medianjoining network analysis by Network 4.6 (Bandelt et al., 1999) to visualize relationships among all haplotypes.

2.4 Divergence time estimates

We employed the relaxed molecular clock method incorporated in BEAST 1.8.4 to estimate the divergence time of each lineage (Drummond and Rambaut, 2007). Calibration information was according to C. mongolicus-like fossil records detected in the Pliocene (Liu and Su, 1962; Chen, 1998). The means values (4±1.9 Ma) of the normal distributions were used in the analysis. The Yule process was used as the tree topology prior, and the TrN and HKY substitution models were employed for C. alburnus and C. mongolicus, respectively. Posterior distributions of node ages were performed by MCMC for 20 million generations with sampling every 1 000 trees, and the first 10% was treated as burn-in. The effective sample size was used to determine the Bayesian statistical significance of each parameter in Tracer 1.6 (Rambaut and Drummond, 2013). Summarized and annotated phylogenetic trees with a molecular clock constraint yielded by BEAST were displayed by FigTree 1.4.3 (Rambaut, 2009).

2.5 Genetic diversity and population structure

Genetic diversity is reflected in the measures of the number of haplotypes (n), haplotype diversity (π) and nucleotide diversity (h) (Nei, 1987), using DnaSP 5.10 (Librado and Rozas, 2009) software calculate values of each population. The mean genetic distance within and between populations were calculated according to the Kimura 2-parameter (K2P) model (Kimura, 1980) as implemented in MEGA7, standard errors were estimated by bootstrapping using 10 000 replicates. Pairwise genetic differentiation between populations was estimated using F-statistics (F_ST) (Weir and Cockerham, 1984) with 10 000 permutations, based on the distance method in Arlequin 3.5 (Excoffier and Lischer, 2010). To investigate significant genetic structure between populations, we performed a two-level hierarchical analysis of molecular variance (AMOVA) using Arlequin. Significance was assessed by 10 000 random permutations. Gene flow (N_m) was evaluated from GammaST (Nei, 1982), and F_ST (Hudson et al., 1992) was evaluated using DnaSP 5.10.

Mantel tests were used to test the relationships between pairwise values of F_ST and the geographical distance. Geographic distances among populations were estimated using Google Earth 7.3.0.

2.6 Demographic history

The population demography for topmouth culter and Mongolian redfin populations were examined using three approaches. First, Tajima's D (Tajima, 1989) and Fu's Fs (Fu, 1997) were applied to check the evidence of neutral evolution in Arlequin 3.5 (Excoffier and Lischer, 2010). In the neutrality test, we assessed significance by generating null distributions from 10 000 coalescent simulations.

Second, the demographic history was investigated by comparing mismatch distributions in each basin population sample using DnaSP 5.10 and Arlequin3.5. The time of possible population expansions (t) was estimated via equation τ=2ut (Nei and Tajima, 1981; Rogers and Harpending, 1992), where u is the mutation rate for each sequence and each generation. The value of u was calculated by u=2μk, where μ is the mutation rate per nucleotide and k is the number of nucleotides in the analyzed fragment. The approximate time of expansion was calculated by multiplying t by the generation time (C. alburnus: 3 years for CRB, HRB, and ZRB, 4 years for LRB, and 3.5 years for the overall population; C.mongolicus: 2 years for CRB, 4 years for LRB, and 3 years for the overall population) (Fish Research Laboratory, Hubei Institute of Hydrobiology, 1976; Xu et al., 2009; Yang, 2009). The substitution rate of 3% per million years for freshwater fish species was generally applied in the CR (Donaldson and Wilson, 1999). The goodness-of-fit of the actual distributions were calculated by Harpending's raggedness index (Hri) (Harpending, 1994) with the expected distributions using a model of population expansion, and the sum of squared deviations (SSD) by assessing significance with 1 000 permutations.

Third, Bayesian skyline plots (BSPs) (Drummond et al., 2005) were generated using BEAST 1.8.4 for CR with 20 million generations to explore demographic histories. This analysis was performed using the TrN substitution model for C. alburnus and the HKY substitution model for C. mongolicus, and no partition into codon positions. The coalescent tree prior was specified as the BSP. The uncorrelated lognormal relaxed-clock model, with means value (4±1.9 Ma) of the normal distributions was used. The remaining settings were as default.

3 RESULT 3.1 Sequence information

The mitochondrial CR successfully amplified 327 C. alburnus from 11 populations and 223 C.mongolicus from 9 populations. C. alburnus sequences (GenBank Accession Nos. MG584736–MG584792) and C. mongolicus sequences (GenBank Accession Nos. MG584793–MG584829) were deposited in GenBank (Table S1). For C. alburnus, the aligned sequences (842 bp) contained 58 variable sites, of which 47 were parsimony-informative. For C. mongolicus, the aligned sequences (932 bp) contained 51 variable sites, of which 49 were parsimony-informative. No insertions or deletions were detected in CR.

3.2 Culter alburnus phylogenetic analysis and divergence time estimate

Based on C. alburnus CR haplotypes, BI and ML methods yielded similar tree topologies and consistently supported two lineages (Fig. 2a, only BI tree shown). Lineage I, the sister group of Lineage II, consisted of specimens from the ZRB. Lineage II included samples from the CRB, HRB, and LRB. The median-joining network analysis revealed a "weblike" topology pattern, whereas there were significant differences between the ZRB nodes versus the other three basin nodes (Fig. 3a).

We conducted divergence time estimation among different haplotypes of C. alburnus using an uncorrelated lognormal relaxed-clock model. Based on the 95% credibility intervals (95% CI) of the posterior distribution (Fig. 4a), C. alburnus diverged from outgroups such as C. idellus in the Pliocene 5.25 Ma (95% CI: 4.18–6.79). The lineage I group diverged at the same series 2.74 Ma (95% CI: 1.96–3.48), and the other three basins appeared to share an obscure common ancestor and all diverged in the Pleistocene.

3.3 Culter alburnus genetic diversity and population structure

For C. alburnus, the number of haplotypes (n), haplotype diversity (h), and nucleotide diversity (π) for each population are presented in Table 2. The overall CR h (0.913±0.008) and π (0.008 0±0.005 1) were high. A total of 57 haplotypes were detected from all sequences. The highest haplotype and nucleotide diversities were discovered in the Kuilei Lake (KLL) and Liangzi Lake (LZL) populations, whereas the lowest were discovered in the Weishan Lake (WSL) and STR. The haplotypes displayed restricted basin geographical distributions, with 98.2% of the haplotypes being unique for a single basin. Haplotype 7 (GenBank Accession No. MG584742) was the most common (shared by 64 individuals) and widespread, being shared among all CRB and HRB populations.

Based on the K2P model, the mean genetic divergence within and between each populations are presented in Table 3. The overall average genetic divergence among individuals was 0.007 6 for C. alburnus. Among the 11 populations, the overall mean intrapopulation genetic distance of the LZL population (0.004 0) was largest and Xingkai Lake (XKL) (0.000 8) population was smallest. Genetic distances among CRB and HRB populations (0.002 1–0.005 6) were very small, but genetic distances between ZRB and other basin populations (0.020 2–0.021 7) were higher.

For population structure analysis, all samples were grouped by clades and basins. An AMOVA revealed that most of the variation (79.45%, P < 0.000 1) occurred among clades (Table 4), whereas the variance among populations within groups (4.92%, P < 0.000 1) and within populations (15.64%) was relatively small (F_ST=0.843 7, P < 0.000 1), indicating a high level of clade structure. Pairwise comparisons of genetic differentiation (F_ST) values varied from -0.015 5 to 0.937 8. All pairwise F_ST comparisons between basin populations were statistically significant. For CRB, significant genetic differences were observed between Shengjin Lake (SJL) and other CRB populations. The values of N_m calculated from GammaSt (from 0.15 to 3.53) and F_ST (from 0.07 to 0.45) indicated a limited gene flow level among the four C. alburnus basin populations, especially between the ZRB and the other three basin populations, which was in accordance with their significant F_ST values (Table S2).

There was a significant positive relationship between population genetic differentiation and geographic distance (Mantel tests, R=0.502 6, P=0.023 6 < 0.05) (Fig. 5a).

3.4 Culter alburnus demographic history

For C. alburnus, neutrality tests were negative except for the HRB, and Fu's F_S test values were statistically significant for ZRB, CRB, and LRB populations (P < 0.05) (Table 5). Mismatch analysis revealed an approximately unimodal distribution of the overall population and ZRB, CRB, and LRB populations (Fig. 6 left - a, b, c, e), but the shape of the mismatch distribution of HRB was approximately bimodal (Fig. 6 left - d). Nevertheless, the model of sudden demographic expansion for the overall population was not rejected by the SSD and Hri values (Table 5). The total population of τ was calculated as between 1.2 and 3.3 by Arlequin using the evolutionary rate calibrated in this study (3% per Ma), and the expansion time for each basin and overall was estimated to be between 0.035 and 0.098 Ma (Table 5). Additionally, BSP suggested that remarkable expansion occurred in the overall population and CRB and LRB populations (Fig. 6 right-a, c, e), and that moderate expansion occurred in the ZRB population (Fig. 6 right - b), but that the effective population size was relatively stable for the HRB population (Fig. 6 right - d).

3.5 Culter mongolicus phylogenetic analysis and divergence time estimate

Based on C. mongolicus CR haplotypes, the BI and ML methods yielded similar topologies and consistently supported two lineages (Fig. 2b, only the BI tree is shown). Lineage I, the sister group of Lineage II, consisted of specimens from the ZRB. Lineage II included samples from CRB and LRB. The median-joining network analysis revealed a "web like" topology pattern, and some C. mongolicus individuals from CRB and LRB population haplotypes shared the same node (Fig. 3b). Specifically, haplotypes 6, 7, 8, and 12 are shared by CRB and LRB populations, whereas significant differences between the ZRB nodes versus the other two basin nodes are seen.

Based on the 95% CI of the posterior distribution (Fig. 4b), C. mongolicus diverged from outgroups in the Pliocene 5.28 Ma (95% CI: 4.61–6.79 Ma), and lineage I group diverged at the same time; the lineage II basin appeared to share a common ancestor 3.60 Ma (95% CI: 2.84–4.43 Ma).

3.6 Culter mongolicus genetic diversity and population structure

For C. mongolicus, the values of n, h, π for each population are presented in Table 2. The overall CR h (0.903±0.008) and π (0.012 5±0.001 2) were high. Thirty-seven haplotypes were detected among all sites. The highest haplotype and nucleotide diversities were discovered in the KLL population. The haplotype displayed restricted basin geographical distribution, with 89.2% of the haplotypes being unique for a single basin, whereas there was only one haplotype (37, GenBank accession No. MG584829) for the ZRB (STR population).

The mean genetic divergence within and between each of the populations are shown in Table 3, and the overall average genetic divergence among individuals was 0.012 9 for C. mongolicus. Among the 9 populations, the overall mean intrapopulation genetic distance of the KLL population (0.004 8) was largest and the STR (0) population was smallest. The genetic distances among CRB and HRB populations (0.000 9–0.005 1) were very small, but those between the ZRB and other basin populations (0.040 6–0.043 8) were higher.

For population structure analysis, all samples were grouped by clades and basins. An AMOVA revealed that most of the variation (92.33%, P < 0.000 1) occurred among groups (Table 4), whereas the variance within populations (5.53%, P < 0.000 1) and among populations within groups (2.14%, P < 0.000 1) was relatively small (F_ST=0.9447, P < 0.000 1), indicating a high level of clade structure.

Pairwise comparisons of genetic differentiation (F_ST) values varied from 0.006 8 to 0.990 8. All pairwise F_ST comparisons between basin populations were statistically significant. For the CRB, significant genetic differences were observed between CRB populations, with the exception of the DTH and LZL. The values of N_m calculated from GammaSt (from 0.14 to 3.32) and F_ST (from 0.06 to 1.05) indicated a limited gene flow level among C. mongolicus basin populations, especially between the ZRB and the other two basin populations, which was in accordance with their significant F_ST values (Table S2).

There was a significant positive relationship between population genetic differentiation and geographic distance (Mantel tests, R=0.592 7, P=0.007 1 < 0.05) (Fig. 5b).

3.7 Culter mongolicus demographic history

For C. mongolicus, the ZRB population was excluded from demographic analysis because of the unique haplotype. The neutrality test was negative for the overall population and the CRB and LRB. All Fu's F_S test values were significantly negative (P < 0.01) (Table 5). The mismatch distribution shape of C. mongolicus was unimodal (Fig. 6 left -f, g, h). Nevertheless, the model of sudden demographic expansion for the overall population and each basin population was not rejected by the SSD and Hri values (Table 5). The total populations of τ were calculated between 2.5 and 4.2 using Arlequin, and the expansion times for each basin and the overall population were estimated to be between 0.083 and 0.099 Ma (Table 5). Similar results were obtained with BSP, suggesting that moderate expansion occurred in the overall population and the CRB and LRB (Fig. 6 right - f, g, h).

4 DISCUSSION 4.1 Population divergence and implication of cryptic subspecies

Phylogenetic analysis based on CR revealed two lineages (I and II), which is highly consistent with the geographical patterns for C. alburnus and C. mongolicus (Fig. 2). The divergence time estimate revealed a deep divergence analogous to phylogenetic analysis. Lineage I was the first to split from Lineage II in the Pliocene and exhibits a high divergence from it; Lineage II diverged during the Pleistocene but at a later point in time. Lineage I and populations from Guangdong Province and Hainan Province were found to show a close relationship, though they are separated by the Qiongzhou Strait. However, this may be because the Qiongzhou Strait was once part of the coastal plain of the China mainland continent during the Pleistocene glaciations (Wang et al., 2012).

There is a deep genetic divergence between ZRB populations and other populations (C. alburnus, 0.020 2–0.021 7; C. mongolicus, 0.040 6–0.043 8; Table 3) based on genetic distance calculations. This level of divergence has previously been attributed to subspecies-level differences (Shaklee et al., 1982). Therefore, we suggest that C. alburnus and C. mongolicus of ZRB populations represent cryptic subspecies. BI tree analysis and web-like statistical haplotype network findings also support this. The significant pairwise genetic differentiation detected between the ZRB populations and other populations is in accordance with that of many Cyprinidae species in China such as Acheilognathus macropterus (Zhu and Liu, 2007), C. idellus (Fu et al., 2015), Hypophthalmichthys molitrix (Ji et al., 2009), and H. leucisculus (Chen et al., 2017). It is conceivable that geographical isolation is a key factor in population differentiation. For example, the ZRB was isolated from the CRB by the Nanling Mountains. Moreover, compare to the geographic distance between the two ZRB populations (XFJ and STR) is 710 km, the distance between the DJHR (CRB) and XFJ (ZRB) is only 260 km, but the two basin populations show completely different haplotypes. Gene flow has also been detected among the basin populations. Additionally, there is evidence of limited gene flow between the ZRB and other basin populations. Second, historical factors have played an important role in determining the patterns of genetic variability. Extant Culter originated from the Pliocene, and the ZRB was also diverged during the late Pliocene period (Wang et al., 2012). Therefore, the ZRB populations were isolated from their point of origin at an early stage of Culter divergence. Divergence time estimations support this hypothesis. Third, we consider that ZRB belong to peripheral area and play a role in genetic differentiation. Moreover, Culter species are a temperate fish (Froese and Pauly, 2017) but the ZRB is a tropical river (Dudgeon, 2000), and peripheral populations tend to occur in more unsuitable environments than central populations (Lesica and Allendorf, 1995).

4.2 Species-level of genetic diversity and structure

Based on haplotype and nucleotide diversities, C. alburnus and C. mongolicus both showed high overall genetic diversities as seen in previous findings (Wang et al., 2007; Qi et al., 2013, 2015; Yang et al., 2016). This could reflect the complex and variable heterogenetic habitat of Culter populations, which are widely distributed in China. However, the genetic diversity within populations shows that some populations have high levels of diversity, while others have low levels. All CRB populations typically show higher genetic variability than the other three basin populations for C. alburnus, whereas a similar level of variability was found among the C. mongolicus populations. Genetic patterns and species diversity may have been affected by speciation events and subsequent radiations (Xiao et al., 2009; Yamaguchi et al., 2010). Additionally, differences in the time of species origins have contributed to the general positive correlation with genetic diversity (Messmer et al., 2012). Populations near the center of a species' range are usually contiguous and at a high density, while populations at the margins tend to be isolated and sparse. Thus, peripheral populations are small and isolated, with low genetic variability (Brussard, 1984; Lesica and Allendorf, 1995). We therefore deduced that CRB populations were at or near the core region of their distribution. This result was also supported by the fact that the RPEC was the place of origin of Culter species (Liu and Su, 1962; Chen, 1998; Wang, 2013).

Other notable discoveries include the low genetic diversity of SJL and Dangjiangkou Reservoir (DJKR) population and the significant genetic differentiation between SJL and other CRB populations. This result is mainly caused by geographic isolation. The SJL sluice and DJKR reservoir were constructed between 1950–1970, and were isolated from the main stream of the Changjiang River. Major threats to SJL include artificial reclamation, serious sedimentation accumulation, and habitat degradation (Jiang et al., 2007). The DJKR dam prevents fish from moving through, thereby inducing geographic isolation between the DJKR and the Changjiang River. These factors can affect genetic diversity by decreasing gene flow among the CRB.

The main aim of population genetic studies is to evaluate levels of genetic connectivity among populations. Major factors that shape genetic structure among populations include population size, dispersal potential, environmental conditions, and life history strategies (Gaggiotti et al., 2009). In the present study, AMOVA and pairwise F_ST values indicated that high levels of population structure and genetic differentiation occurred in different drainage basin populations for C. alburnus and C. mongolicus. The fact that 98.2% and 89.2% of haplotypes were only detected in single C. alburnus and C. mongolicus basin localities, respectively, indicates that limited gene flow occurred among drainage basin populations.

Allopatric speciation has a central role in generating new species by reproductive isolation (April et al., 2013). Freshwater fish populations are usually confined to confined areas during a relatively short to moderate evolutionary period (Dewoody and Avise, 2000; Yamaguchi et al., 2010). During the Pleistocene glacial period, temperate freshwater fish experienced landscape changes s in intermittent connected and isolated waterways and glacial refugia, leading to complex phylogeography and many opportunities to study recent or incomplete speciation (Behnke, 1972; Bernatchez and Wilson, 1998; Taylor, 1999; Ruzzante, 2011; Ruskey and Taylor, 2016). The clear pattern of IBD in C. alburnus and C. mongolicus in the current study supports the hypothesis and was predicted by the drainage basin isolation. Geographic distance could be a reasonable explanation for the observed geographic distribution of genetic variation between populations. The two largest carnivorous Culter species mainly coexist in larger rivers, lakes, and reservoirs in China. Perhaps because of the demand for food resources they rarely inhabit small tributaries, and gene flow among basins may be rare too. Moreover, because C. alburnus and C. mongolicus spawn in flowing water (Fish Research Laboratory, Hubei Institute of Hydrobiology, 1976) and are potamodromous, gene flow among populations within the same drainage basin may be more common.

4.3 Demographic process of Culter

Climatic oscillation played an important role in the contemporary diversity of many species and communities. In particular, climate changes between ice age maxima and warm interglacials caused large changes in the distribution of species, which ultimately affect the distribution of genetic variation within species (Hewitt, 2003). During the Pleistocene, five glaciations occurred in Eastern China (Yi et al., 2007). The last interglacial period occurred at approximately 0.1 Ma, while glacial retreat occurred at 0.075 Ma. The last glaciation (LG) in Eastern China started at approximately 0.073 Ma and continued until 0.016 Ma, and the last glacial maximum (LGM) development in Eastern China occurred during the 0.024–0.016 Ma.

On the basis of estimation expansion times from mismatch distribution, most Culter populations (C. alburnus, 0.035–0.098 Ma; C. mongolicus: 0.083– 0.099 Ma) underwent a period of population expansion that began before the LGM in Eastern China. This is related to suitable environments within their areas during interglacial periods (0.1–0.075 Ma). Typical features of this type of expansion pattern are supported by high haplotype diversity but low nucleotide diversity, significantly negative F_S values, unimodal distribution of pairwise differences, and expansion BSP analysis (Nei, 1987; Tajima, 1989; Rogers and Harpending, 1992; Fu, 1997; Drummond et al., 2005).

However, BSP and Tajima's analysis suggest that effective population size was relatively low and stable for the HRB C. alburnus population. Two reasons account for this abnormal phenomenon. First, according to divergence time estimations, the HRB was the latest diverged group (1.48 Ma) for C. alburnus. Furthermore, its genetic diversity was lowest in C. alburnus. Since the Pleistocene ice age, the population bottlenecks appear to have impacted the extent of haplotype diversity in many freshwater fish populations (Billington and Hebert, 1991) and were influenced by the small founder effect later (Hoffman and Blouin, 2004). Second, habitats of the HRB (WSL, 34°34'45"N), which were distributed in the southern boundary of the permafrost area in the hilly plain region of Eastern China (near 33°20'N–33°40'N) (Pu, 1991), were unsuitable for C. alburnus during the LG period because they are temperate fish.

4.4 Conservation implication

Culter alburnus and C. mongolicus are economically valuable freshwater fish that are widely distributed in larger rivers, lakes, and reservoirs of China. However, because of their carnivorous nature, they were once considered harmful and were removed from water bodies. Today, some anthropogenic activities such as illegal fishing, random stock enhancement, and habitat degradation (e.g., dam construction and water pollution) are major threats to these species and result in the miniaturization of natural fish resources. Indeed, during the last few decades, the wild populations have declined so rapidly (Liu et al., 2007) that conservation of their genetic diversity has become a very important issue.

Maintaining genetic diversity is a fundamental concern in conservation and evolutionary biology, and is an important part of biodiversity. Populations lacking genetic variability have low adaptability to stress, while populations with high genetic variability can rapidly respond to environmental changes (Ayala, 1965; Nei et al., 1975; Frankham et al., 2002). We propose several recommendations based on the current study. Each basin (ZRB, CRB, HRB, and LRB) should be considered an individual Management Unit based on the significant genetic divergence occurring in different drainage basin populations. Additionally, the CRB population demands priority protection because it has the highest genetic diversity and provides the largest yield of Culter species, thus we should protect the integrity and continuity of their habitats and reduce human disturbance. Special attention should also be paid to SJL and DJKR populations because dam construction and artificial reclamation will lead to geographic separation and impede gene flow among basins. The most important course of action should be to reconnect the lakes with the Changjiang River, this is especially important for the restoration of genetic diversity in smaller lakes, thus we recommend installation fish passage facilities for DJKR and restoration of the connectivity with Changjiang River for SJL. The C. mongolicus population of the ZRB (STR) also requires a great deal of attention because of its genetic uniqueness, which includes exclusive mtDNA CR haplotypes. Long-term monitoring programs should also be established. Finally, to understand genetic patterns and establish effective conservation strategies for the two Culter species, more sensitive nuclear markers and larger samples should be employed in future research.

5 CONCLUSION

This study assembled a relatively larger scale dataset to infer Culter species phylogeographic patterns in China. The results support the hypothesis that geographic isolation drove divergence between lineages I and II. Based on the molecular clock, Culter species appeared to originate during the early Pliocene. Genetic diversity differences among basin populations are significant for Culter species, and demographic analysis indicated that most populations experienced population expansion during the warm interglacial period. Based on the current study, we propose several suggestions for the protection of genetic diversity: each basin should be considered an individual Management Unit and CRB populations demand priority protection; the C. mongolicus STR population possesses unique genotypes, which calls for a special conservation concern for this species at the population level.

6 DATA AVAILABILITY STATEMENT

Sequence data that support the findings of this study have been deposited in GenBank under accession numbers from MG584736 to MG584792 for C. alburnus and from MG584793 to MG584829 for C. mongolicus. Details regarding individual samples are available in Table S1.

7 ACKNOWLEDGMENT

We would like to acknowledge colleagues at the Fisheries Ecology Lab for their assistance in field sampling.