Daphnia diversity on the Tibetan Plateau measured by DNA taxonomy

Abstract Daphnia on the Tibetan Plateau has been little studied, and information on species diversity and biogeography is lacking. Here, we conducted a 4‐year survey using the barcoding fragment of the mitochondrial COI gene to determine the distribution and diversity of Daphnia species found across the Plateau. Our results show that species richness is higher than previously thought, with total described and provisional species number doubling from 5 to 10. Six of the taxonomic units recovered by DNA taxonomy agreed well with morphology, but DNA barcoding distinguished three clades each for the D. longispina (D. galeata, D. dentifera, and D. longispina) and D. pulex (D. pulex, D. cf. tenebrosa, and D. pulicaria) complexes. The sequence divergence between congeneric species varied within a large range, from 9.25% to 30.71%. The endemic D. tibetana was the most common and widespread species, occurring in 12 hyposaline to mesosaline lakes. The lineage of D. longispina is the first confirmed occurrence in west Tibet.


| INTRODUC TI ON
In the past decade, DNA sequencing has generated abundant molecular information, standard dataset platforms, and universal technical rules for modern taxonomic and biogeographical research . DNA barcoding uses a short DNA sequence in an organism's DNA to compare against that of another organism to determine the degree of relatedness between two closely related organisms. The barcoding fragment of the mitochondrial gene cytochrome c oxidase subunit I (COI) is a popular marker used to identify and differentiate closely related species that are very similar in morphology. It has assisted in species-level identity in many animal groups such as birds (Hebert, Stoeckle, Zemlak, & Francis, 2004), fishes (Ward, Zemlak, Innes, Last, & Hebert, 2005), spiders (Barrett & Hebert, 2005), butterflies (Hebert, Cywinska, Ball, & deWaard, 2003;Janzen et al., 2005), ants (Smith, Fisher, & Hebert, 2005), and crustaceans (Costa et al., 2007;Elías-Gutiérrez, Jerónimo, Ivanova, Valdez-Moreno, & Hebert, 2008), including marine decapods and euphausiids (Bucklin et al., 2007;Costa et al., 2007).
Cladocera is a monophyletic, primarily freshwater crustacean order, one of the three main components of the microcrustacean zooplankton (Dumont & Negrea, 2002). The genus Daphnia (Anomopoda: Daphniidae) has been studied in much detail (Lampert, TA B L E 1 Geographic and environmental data for lakes and ponds where Daphnia populations were sampled  2011), and the full genome of two species (D. magna and D. pulex) has been sequenced (Colbourne, Singan, & Gilbert, 2005;Colbourne et al., 2011). Daphnia is most diverse and abundant in the temperate regions, but is present in all climate zones on all continents, and is often the dominant group in freshwater zooplankton (Benzie, 2005 (Wu, 1980).
More than 25% of the total species identified are endemic (Wu, 1987). This reflects the age of the plateau, the central part of which started rising some 40 million years ago. Previous fragmentary taxonomic studies of Cladocera on the Tibetan Plateau including the genus Daphnia were based solely on morphology (Chiang, 1963;Chiang & Du, 1979;Shen & Sung, 1964). These were updated after several scientific expeditions to the area during the 1970s (Chiang & Chen, 1974;Chiang, Shen, & Gong, 1983). More recently, Möst et al. (2013) and Ma et al. (2015)  In this study, we employed DNA barcoding and DNA taxonomy through analysis of the mitochondrial marker COI to determine species diversity of the Daphnia genus in lakes and ponds on the Tibetan Plateau. We also estimated the number of endemic species in the region and generated a phylogenetic tree based on our mtCOI data and those from GenBank. Our study will greatly improve our understanding of distribution and species diversity in Cladocera and may have important implications for the conservation of the Tibetan Plateau freshwater fauna.

| Sample collection
Sample collection covered a large geographical range (>2,200,000 km 2 ) in different habitats that ranged from 2,700 m to about 5,000 m a.s.l. Zooplankton samples were collected between 2012 and 2015 from 26 permanent lakes and from several riparian temporary ponds (Table 1 and Figure 1). Samples were obtained by vertical hauls with a plankton net that has a mesh size of 100 μm.
The collected samples were fixed in 70% ethanol. Specimens were examined under a dissecting microscope in the laboratory. Sorted individuals were transferred to a fresh tube and preserved in 95% ethanol at 4°C for genetic analysis. We followed Benzie (2005) for species identification and nomenclature of Daphniidae.

| DNA extraction
Total genomic DNA was extracted using a genomic DNA isolation kit (Wizard ® Genomic DNA, Purification Kit type A1225; Promega, USA). We modified the standard protocol as follows: for DNA extraction, specimens were picked out from 95% ethanol, rinsed with double-distilled water, transferred individually to a reaction tube and stored on ice. Next, we added 200 μl warm Cell Lysis solution and 3 μl Proteinase K (20 mg/ml). We vortexed and subsequently incubated the mixture for 2 hr at 65°C and then for 2-3 days at 55°C with daily addition of 2 μl of fresh Proteinase K. Next, we added 100 μl of Precipitation Solution, vortexed vigorously at middle speed for 20 s and put on ice for 2 min, then centrifuged at RCF 15,321 g for 10 min at room temperature. We carefully removed the supernatant and transferred it to a clean 500 μl microcentrifuge tube containing F I G U R E 1 Geographic location of sample collection sites. Color dots: lakes or ponds inhabited by Daphnia species. The background map was generated using SRTM 90 m elevation data 200 μl of isopropanol at room temperature. We centrifuged at RCF 15,321 g for 1 min at room temperature and carefully decanted the supernatant. Finally, the pellet was washed with 70% ethanol, dissolved in 40 μl DNA hydration solution and stored at −20°C.

| Amplification and sequencing of the mitochondrial gene
The barcoding fragment of the mitochondrial cytochrome oxidase I (COI) gene was amplified from total genomic DNA using polymerase chain reaction (PCR) . Primers used for PCR were CO1490F and CO2198R (Folmer, Black, Hoeh, Lutz, & Vrijenhoek, 1994

| DNA taxonomy
The authenticity of all mitochondrial COI sequences was verified by a BLAST search in GenBank. The sequences were assembled and edited in BioEdit (Hall, 1999) and aligned using the CLUSTALW multiple algorithm. The first 20 and last 10 bp were not included because they were missing in some sequences. We added the COI sequences available in public databases to our analysis to ensure that our nomenclature is reliable for each Daphnia species (see Table S1). We used two different approaches to identify taxonomic units from DNA taxonomy, namely the Automatic Barcode Gap Discovery (ABGD) (Puillandre, Lambert, Brouillet, & Achaz, 2012) and Generalized Mixed Yule Coalescent (GMYC) model (Fujisawa & Barraclough, 2013;Pons et al., 2006) to infer putative species boundaries on COI dataset. The ABGD approach tests for a gap in the distribution of the pairwise genetic distances and then identifies groups of individuals united by genetic distances that are shorter than the gap. The method was performed on the COI alignment through an online tool (http:// wwwabi.snv.jussieu.fr/public/abgd/abgdweb.html) with default settings with P (prior limit to intraspecific diversity) ranged between 0.001 and 0.1 and X (gap widths) = 1 using the available models JC86 (Jukes-Cantor) and K80 (Kimura

| Genetic divergence and phylogenetic analysis
Distances between COI sequences were calculated using the Kimura   (Table 1). Sequences of COI were obtained from 93 animals (GenBank accession numbers MG544001-MG544093). Comparing each of our COI sequences with sequences in GenBank, we identified all animals as Daphnia (sequence divergence < 5%). The ABGD model detected a barcode gap in the alignment and suggested that the 93 individuals included 10 taxonomic units (named S1-S10). The

| Patterns of genetic divergence and phylogenetic analysis
The total length of the sequenced segment after alignment was 677 bp. The average base composition was A = 21.70%, C = 20.41%, G = 22.36%, T = 35.53%, and transition/transversion (ti/tv) ratio = 1.751. The uncorrected K2P pairwise distances among species in this study varied between 9.25% and 30.71% and the average pairwise distance was 25.23%. The highest distance was between D.
cf. tenebrosa and D. magna, a value which is slightly higher than the maximum congeneric distance of 30.65% recorded earlier in Daphnia by Costa et al. (2007). Two species of the D. longispina complex, viz.
The best-fitting model selected by MrModeltest 2.3 was GTR+I+G with a relative AIC weight of 0.982 and gamma distribution shape pa-

| Species diversity and genetic divergence
Cladocera have been traditionally regarded as cosmopolitan, but there is mounting evidence for the existence of numerous sibling and cryptic species. In the past 10 years, DNA barcoding has accu-  Garfias-Espejo, 2006), and in other parts of North America (Penton et al., 2004). DNA barcoding in the present study also revealed that Daphnia species diversity on the Tibetan Plateau is much higher than previously thought (Chiang, 1963;Chiang & Chen, 1974;Chiang & Du, 1979;Chiang et al., 1983;Shen & Sung, 1964), Alternatively, species in the D. longispina clade show strong morphological plasticity (Laforsch & Tollrian, 2004;Petrusek, Tollrian, Schwenk, Haas, & Laforsch, 2009) compounded by the possibility of hybridization and introgression (Ishida et al., 2011;Keller, Wolinska, Tellenbach, & Spaak, 2007;Schwenk & Spaak, 1995). Morphologybased taxonomy is insufficient for distinguishing the underlying genetic units. A lack of investigation has long delayed an appreciation of the diversity of Daphnia in Tibet. Only recently have studies begun to reveal the region's hidden species  and the impact of environmental change on cladoceran species richness and composition (Lin et al., 2017).
The genomic region of the COI gene sequence is used not only in DNA barcoding (Costa et al., 2007;Hebert, Cywinska, et al., 2003), but also in detecting speciation. The level of sequence divergence between congeneric species of crustaceans averaged 17.16%, the highest value so far in animals. As a comparison, congeneric species of lepidopterans show just 6.1% variation (Hebert, Cywinska, et al., 2003), birds 7.93% , and fishes 9.93% (Ward et al., 2005). Congeneric divergences in Daphnia are reported by Costa et al. (2007) to be extremely high at 13.18%-30.65%, which is supported by our data (9.25%-30.71%; the highest divergence being between D. cf. tenebrosa and D. magna). The average interspecific divergence between Daphnia species on the Tibetan Plateau of 25.23% is similar to that reported in Argentina (25.28%, Adamowicz, Hebert, & Marinone, 2004), but significantly higher than values reported from Churchill, Canada (14.1%, Jeffery et al., 2011). Difference in average interspecific divergence may be related to species richness: 10 congeneric species occurred in Tibet and 11 South American endemics in Argentina, against only five species in the Churchill region.
The level of intraspecific variation in crustaceans averaged 0.69%, a value that is slightly higher than those reported in other groups (most range from 0.25% to 0.30%). High intraspecific variations in our study were found in D. tibetana (1.60%) and D. pulicaria (1.72%). We collected one population of D. tibetana and D. pulicaria from Zhaling Lake, which is more than 1,600 km away from the other investigated lakes, indicating the elevated divergence values came from Zhaling Lake samples. Possibly the high values reflect limited gene flow between the two species caused by physical barriers such as mountains that separate the Tibetan lakes from Zhaling Lake, followed by adaptation to local environmental pressure in the lake.
Previous investigations showed that D. tibetana is a halobiont, living at more than 4,000 m in hyposaline to mesosaline lakes in Tibet, Mongolia, and India. Our sample area covered a large geographical range (>2,200,000 km 2 ) containing different habitats located at latitudes ranging from 2,700 m to about 5,000 m asl. Water temperatures at which the samples were collected ranged from 2°C to 20°C, salinity varied from 9 to 35 g/L, and pH ranged from 9.0 to 10.4 (Zhao, Wang, Zheng, Zhao, & Wang, 2002). However, a recent study (Lin et al., 2017) reported an even wider salinity range (6.4-46.2 g/L) for D. tibetana. The northernmost population in our TA B L E 2 Genetic diversity, assessed by Kimura two-parameter distance (median, in %) within/between the ten taxonomic units of Daphnia with uniform rates; standard error estimates obtained by neighbor-joining bootstrap procedure with 10,000 replicates  (Wei et al., 2015). Moreover, D. longispina has been documented to typically occur in alpine oligotrophic lakes (Hamrová, Krajicek, Karanovic, Černý, & Petrusek, 2012;Ventura et al., 2014). was absent from east China . Daphnia cf. himalaya is especially intriguing, as it was found in two permanent lakes and one temporary pond along the Nyenchenthanglha Mountain in the center of the Tibetan Plateau. The morphology of D. cf. himalaya in our collection is similar to the dark-pigmented Daphnia-like species described by Manca, Martin, Peñalva-Arana, and Benzie (2006) and named Daphnia himalaya from the Khumbu Region in Nepal.
However, the absence of males in our samples suggests further investigation is needed.

| CON CLUS ION
Our study is the first to use DNA barcoding as a tool to delineate species and their distribution pattern in Tibet. The technique revealed 10 described and provisional species of Daphnia on the Tibetan Plateau. This diversity is double of that shown by previous checklists. The sequence divergence among Daphnia was high and varied between 9.25% and 30.71%. Two species, D. tibetana and D. cf. himalaya, are endemic to the plateau and the Himalayas. The hygrophile D. tibetana, presumed to be the result of local speciation, was the most common species that was found in 12 hyposaline to mesosaline lakes. Our study is the first time to confirm the presence of the D. longispina lineage in western Tibet.