Phylogeography of Orinus (Poaceae), a dominant grass genus on the Qinghai-Tibet Plateau

Phylogeography of Orinus (Poaceae), a dominant grass genus on the Qinghai-Tibet Plateau Abstract To better understand the responses of arid-adapted, alpine plants to Quaternary climatic oscillations, we investigated the genetic variation and phylogeographic history of Orinus, an endemic genus of Poaceae comprising three species from the dry grasslands of the Qinghai-Tibet Plateau (QTP) in China. We measured the genetic variation of 476 individuals from 88 populations using three maternally inherited plastid DNA markers (matK, rbcL and psbA-trnH), the biparentally inherited nuclear ribosomal internal transcribed spacer (nrITS) and amplified fragment length polymorphisms (AFLPs). We found that the plastid DNA, nrITS and AFLPs show considerable, recent differentiation among the species. We detected 14 plastid haplotypes (H1–H14), of which only three were shared among all species, and 30 nrITS ribotypes (S1–S30), of which one (S10) was shared between two species, O. kokonoricus and O. intermedius, but absent in O. thoroldii. The nrITS types formed clades that were inconsistent with species boundaries. Based on these data, we propose and illustrate a complex hypothesis for the evolutionary history of Orinus involving lineage sorting and introgression, the latter of which may explain the shared S10 nrITS type. The AFLP results showed clades corresponding to current species delineation and suggest that lineage sorting in the genus is probably complete. We estimated the crown age of Orinus to be 2.85 (95% highest posterior density: 0.58–12.45) Mya (late Pliocene), and subsequent divergence occurred in the Quaternary. Early divergences were allopatric. More recently, Orinus probably underwent regional expansions corresponding to Quaternary climatic changes, especially glaciation, which is consistent with our divergence time estimates. These climatic changes could have facilitated the S10 event and other hybridization events. Our data also suggest that species of this small genus of grasses survived the Quaternary glacial period in the extremely adverse habitats of the QTP. INTRODUCTION The Quaternary period began c. 2.58 Mya and has been characterized by climatic oscillations, especially glacial and interglacial cycles in the Northern Hemisphere (Shackleton & Opdyke, 1973). The glacial cycles co-varied with and probably profoundly affected other aspects of the climate, including the intensity of the Asian monsoon, even in unglaciated regions (An et al., 2001; Jiang et al., 2011). The Quaternary climatic changes caused changes in geographical distributions and population demography of plant species and consequently left lasting, detectable genetic imprints within and among species (Abbott et al., 2000; Avise, 2000; Hewitt, 2004, 2011; Qiu, Fu & Comes, 2011; Qiu et al., 2013; Wen et al., 2014, 2016). In Europe and North America, the fossil records of plant species and phylogeographic analyses indicate that a common pattern of geographical range shift was to retreat southward and to lower elevations during glacial periods and then to rapidly recolonize the northern areas and higher elevations during the interglacial and postglacial periods (Nason, Hamrick & Fleming, 2002; Stewart et al., 2010; Sakaguchi et al., 2011; Li et al., 2012; Segovia, Pérez & Hinojosa, 2012; Voss, Eckstein & Durka, 2012; Tzedakis, Emerson & Hewitt, 2013; de Lafontaine et al., 2014). In comparison with Europe and North America, the genetic structure and variation in Quaternary phylogeographic histories of plant species are much more pronounced in topographically complex regions of China, such as the Qinghai-Tibet Plateau (QTP) (e.g. Zhang et al., 2005; Meng et al., 2007; Wang et al., 2009a, 2014; Opgenoorth et al., 2010; Xu et al., 2010; Qiu et al., 2011; Zou et al., 2012; Wen et al., 2014; Liu et al., 2015; Zhang et al., 2015; Wan et al., 2016). The QTP and adjacent mountain ranges, such as the Hengduan and the Hindu-Kush Mountains, (hereafter QTP) were formed by the collision of the Indian subcontinent with Eurasia and contain the world’s tallest mountain (Mount Everest) and have an average elevation of > 4000 m (Zheng, 1996). The QTP is regarded as a hotspot for global biodiversity and supports many endemic genera and species (Wu, 1988; Mittermeier et al., 2005). The exceptional biodiversity of the QTP has been attributed to Quaternary climatic oscillations and Miocene-Pliocene or Miocene-Quaternary phases of uplifts of the QTP (Liu et al., 2002, 2006; Liu, 2004; Wang et al., 2009c; Mao et al., 2010; Xu et al., 2010; Jia et al., 2012; Wen et al., 2014; Liu et al., 2015). Climatic oscillations and uplifts may lead to allopatric speciation, detectable as genetic imprints that date to the Miocene, Pliocene or Quaternary. Recent studies have shown evidence for the importance of Quaternary climatic oscillations in explaining genetic variability, distributions and demographic patterns observed in modern populations of QTP species (Tang & Shen, 1996; Bawa et al., 2010; Jia et al., 2011; Liu et al., 2012, 2014b, 2015; Wan et al., 2016). In particular, many alpine species retreated to eastern or south-eastern refugia during glacial periods and recolonized the central plateau platform during the interglacial or postglacial stages (e.g. Zhang et al., 2005; Meng et al., 2007; Chen et al., 2008; Yang et al., 2008; Wu et al., 2010; Tian et al., 2011; Jia et al., 2012; Liu et al., 2015; Wan et al., 2016). In contrast, some other cold-tolerant alpine species were able to survive in situ through glacial periods in multiple ice-free refugia or microrefugia (Wang et al., 2009b; Opgenoorth et al., 2010; Jia et al., 2011; Ma et al., 2014). These two major types of evolutionary hypotheses about phylogeographic relationships between plant species and Quaternary climatic oscillations have been developed and tested by previous studies on the QTP (Wu et al., 2010; Jia et al., 2011; Tian et al., 2011; Ma et al., 2014; Liu et al., 2015). Results show that there was a close relationship between geographical distribution patterns of plant species from the QTP and Quaternary climatic oscillations, and independent genetic lineages appear in different plateau regions due to isolations among developed ice sheets. Given the probable impacts of ongoing climatic change on biodiversity, phylogeographic investigations of plant species in diverse regions, such as the QTP, are increasingly important. Although many phylogeographic studies of QTP plants have been performed, these have mostly focused on woody plants or herbaceous dicot lineages that occur in forested areas (e.g. Zhang et al., 2005; Meng et al., 2007; Chen et al., 2008; Tian et al., 2011; Jia et al., 2012; Wan et al., 2016). The grasslands comprise a significant portion of the QTP and are dominated by grasses (Poaceae). Therefore, disentangling the phylogeographic history of grasses is crucial to understanding the evolution of the QTP flora. In addition, these kind of grass-dominated landscapes are interesting since many of the vegetation types thus far studied are dominated by trees or shrubs. Grasses are cold-, warm- and arid-tolerant, making their study necessary to describe evolution in the QTP. Therefore, phylogeographic studies of grasses are needed to understand the past history of biogeographic and demographic change linked to past climatic oscillations and to make informed inferences about the future of the QTP grasslands (e.g. Qiu et al., 2011; Liu et al., 2012, 2014b, 2015; Wen et al., 2014). Here, we studied the phylogeographic history of the perennial grass genus, Orinus Hitchc., that is endemic to alpine, dry grasslands of the QTP and occurs in unusually sandy habitats (Chen & Phillips, 2006; Su & Cai, 2009; Su et al., 2015). Orinus comprises three morphologically and ecologically distinct species, O. thoroldii (Stapf ex Hemsl.) Bor, O. kokonoricus (K.S.Hao) Keng and O. intermedius X. Su & J. Quan Liu (Su et al., 2015, 2017), that are distributed at elevations between 2200 and 4800 m. Recently, Orinus was placed in subtribe Orininae P.M.Peterson, Romasch. & Y.Herrera with Cleistogenes Keng (Peterson, Romaschenko & Arrieta, 2016; Soreng et al., 2017). Furthermore, Orinus reproduces mainly through clonal reproduction via root suckers. Orinus spp. are dominant in places where they occur and are an exceptional study system for understanding the effects of the Quaternary climatic oscillations on the demography of QTP grassland species. Specifically, we investigated the biogeographic and phylogeographic history in Orinus using robust sampling of the three species throughout their geographical ranges using molecular markers, including plastid DNA, nuclear ribosomal internal transcribed spacer (nrITS) and amplified fragment length polymorphism (AFLP) data sets. We address the following questions: (1) What is the phylogeographic pattern in Orinus? (2) Are the phylogeographic patterns in Orinus consistent with southern glacial refugia and with other demographic patterns observed among diverse alpine plants on the QTP? MATERIAL AND METHODS Population sampling We sampled 476 individuals belonging to 88 populations from the distributional area of the genus including the type localities of all three species (Fig. 1A; see also Supporting Information, Table S1). From each population, we collected fresh leaf blades from three to 12 individuals and preserved them in silica gel. In each population, we selected individuals that were at least 20 m away from each other to avoid clones, and we recorded the latitude, longitude and elevation of each population using an Etrex GIS (Garmin, Taiwan, China). We also prepared several voucher specimens from each population and deposited them at the Herbarium of the Northwest Plateau Institute of Biology (HNWP), the Chinese Academy of Science, Xining, Qinghai Province, China. Figure 1. View largeDownload slide Sampling sites and geographical distribution of plastid DNA haplotypes observed in three Orinus spp. A, sampling region (dashed-line box). B, documented distribution range (dotted line), sampling sites and plastid DNA haplotype frequencies per population. Each pie graph represents a plastid DNA haplotype, and size of pie is proportional to the number of individuals bearing that haplotype. C, network for 14 plastid DNA haplotypes detected from three Orinus spp. Dinebra viscida was used as the outgroup. Sizes of network circles are proportional to haplotype frequencies. Numbers on the link lines between haplotypes represent the mutational sites. Colours for haplotypes are described in (B). Figure 1. View largeDownload slide Sampling sites and geographical distribution of plastid DNA haplotypes observed in three Orinus spp. A, sampling region (dashed-line box). B, documented distribution range (dotted line), sampling sites and plastid DNA haplotype frequencies per population. Each pie graph represents a plastid DNA haplotype, and size of pie is proportional to the number of individuals bearing that haplotype. C, network for 14 plastid DNA haplotypes detected from three Orinus spp. Dinebra viscida was used as the outgroup. Sizes of network circles are proportional to haplotype frequencies. Numbers on the link lines between haplotypes represent the mutational sites. Colours for haplotypes are described in (B). DNA extraction, amplification, sequencing and AFLP fingerprinting We extracted total DNA from silica gel-dried leaves using a modified 2× CTAB procedure (Doyle & Doyle, 1987). We used universal primers to amplify nuclear ITS and three plastid DNA regions: matK, rbcL and psbA-trnH (White et al., 1990; Sang, Crawford & Stuessy, 1997; Feng et al., 1998; Sun, McLewin & Fay, 2001). We performed the amplification using PCR in 25 µL reaction volumes, containing 1 µL plant template DNA, 0.2 mM MgCl2, 0.25 mM dNTPs, 2 µM each primer, 2.50 µL 10× PCR buffer and 0.25 µL Taq DNA polymerase (5 U/µL; TakaRa, Dalian, China). Our PCR thermocycling protocol comprised an initial 5-min denaturation period at 94 °C, followed by 36 cycles of denaturation at 94 °C for 50 s, annealing at 49–58 °C for 50 s and extension at 72 °C for 1 min, with a final extension at 72 °C for 10 min. We purified the PCR products with a TIAN-quick Midi Purification Kit (Tiangen, Beijing, China) in accordance with the manufacturer’s instructions, and we sequenced the purified products using an ABI 3130XL DNA Analyzer (Applied Biosystems, Foster City, CA, USA). We aligned the sequences of the plastid DNA and nrITS separately in MUSCLE (Edgar, 2004) and refined the alignments manually in MEGA 5 (Tamura et al., 2011). We combined the three plastid DNA regions because we did not detect phylogenetic incongruence among them using an incongruence length difference test under the parsimony criterion with 1000 replicates, Tree Bisection and Reconnection branch swapping and the number of trees retained for each replicate limited to 1000 (Farris et al., 1995). All newly obtained sequences of Orinus have been submitted to GenBank under the accession numbers KP302391–KP303274. We used a slightly modified protocol (Vos et al., 1995) to examine the AFLP of all the sampled populations (Supporting Information, Table S1) and used fluorescent primers for selective amplification (Zuo et al., 2015). DNA was digested with restriction enzymes PstI and MseI (40 U/µL, Beijing Dingguo Biotechnology Co., Ltd, China) for 5 h at 37 °C in a 20 µL volume, followed by ligation (40 units T4 DNA ligase, 4 h at 8 °C) to double-stranded PstI and MseI adapters. Two rounds of PCR amplification followed. Of 15 tested PstI/MseI primer combinations, eight were selected to produce the most repeatable and unambiguously scorable profiles. Selective fluorescence-labelled amplification products were separated and analysed on an ABI PRISM 377 DNA Sequencer (Applied Biosystems, Foster City, CA, USA) using GeneScan ROX-500, with an internal size standard. We used GeneScan 3.1 (Applied Biosystems, Foster City, CA, USA) to score the bands of AFLP amplification products. We imported the raw data into Binthere (Garnhart, 2001) and MG (Zhou & Qian, 2003) to generate a 0/1 matrix for data analyses. Population descriptive statistics We calculated the haplotype diversity (HE) of the plastid DNA for each population using DnaSP 5.0 (Librado & Rozas, 2009), and we used Permut (Pons & Petit, 1996) to find the average genetic diversity within populations (HS), total genetic diversity (HT) and the coefficients of differentiation (GST, NST) for each species, separately and combined. For the plastid DNA, we used permutation of GST and NST with 1000 replicates to infer significant phylogeographic structure (Pons & Chaouche, 1995). To examine population genetic structure, we used an analysis of molecular variance (AMOVA) with 1000 permutations (Excoffier, Smouse & Quattro, 1992) in Arlequin 3.11 (Excoffier, Laval & Schneider, 2005). We calculated plastid gene flow among populations (Nem) based on the equation Nem = ([1/FST] − 1)/2, where Ne was the effective population size of female individuals and m was the migration rate (Freeland, Kirk & Petersen, 2011). In nrITS, we assumed that a site was heterozygous if it had a double peak at the same position in both strands, and the weaker signal between base calls was at least 25% of the strength of the stronger (Fuertes-Aguilar, Rosselló & Nieto-Feliner, 1999; Fuertes-Aguilar & Nieto-Feliner, 2003). We also determined nrITS ribotypes with PHASE 1.0 for heterozygous individuals (Stephens, Smith & Donnelly, 2011). We used the Bayesian implementation of the k-mean method in STRUCTURE 2.2 (Pritchard, Stephens & Donnelly, 2000) to estimate the most likely number of populations (K) in Orinus according to the plastid DNA and nrITS data sets. Specifically, we performed analyses for K values ranging from 2 to 7, with ten replicates for each value of K and a burn-in of 2 × 106. We applied the ‘no admixture model’ and independent allele frequencies for these analyses. The most likely number of clusters was estimated according to the K resulting in the best improvement in the global likelihood (Evanno, Regnaut & Goudet, 2005). Network and phylogenetic analyses We performed network and phylogenetic analyses to reconstruct relationships among populations and species. We first distinguished plastid DNA haplotypes and nrITS ribotypes with DnaSP 5.0 software package (Librado & Rozas, 2009). We constructed a median-joining network of haplotypes of the plastid DNA and ribotypes of nrITS with NETWORK 4.5.0.0 (Weir, 1996). In NETWORK, we set both site mutations and indels to evolve with equal likelihood, and we assumed that each indel originated independently of all other indels. For reconstructing phylogenetic relationships, we applied maximum parsimony (MP) and maximum likelihood (ML) to the plastid DNA haplotypes and nrITS ribotypes, independently in PAUP*4.0b10 (Swofford, 2002). We also reconstructed phylogenetic relationships with Bayesian inference (BI) in MrBayes 3.1 (Ronquist & Huelsenbeck, 2003). For the ML analyses and BI, we used a general-time-reversible model (GTR) with a gamma distribution of rate variation among sites (+ G) based on results in PAUP*4.0b10 and jModelTest according to the Akaike information criterion (Posada, 1998; Swofford, 2002). We approximated the gamma distribution using four rate categories. For other parameters, we used default settings. For the MP analyses, we performed heuristic searches with 1000 random addition sequence replicates with ten trees retained at each step during the stepwise additions and with tree-bisection-reconnection for branch swapping. We set PAUP*4.0b10 to treat all characters as unordered and equally weighted, gaps as missing and multistate data as uncertain in the MP analyses. For MP and ML, we calculated bootstrap values from 1000 replicates (Felsenstein, 1985). The BI consisted of two parallel runs with four incrementally heated chains and three million generations sampled every 1000 generations. The output was assessed for convergence using Tracer v.1.3 (Rambaut & Drummond, 2007), and summary statistics and trees were generated using the last two million generations. For network and phylogenetic analyses of plastid DNA and nrITS, we used Dinebra viscida (Scribn.) P.M.Peterson & N.Snow and four samples of Cleistogenes as outgroups, respectively (Peterson et al., 2012, 2016). For the AFLP data, we analysed a rectangular binary 0/1 data matrix using the NTSYS-pc 2.l (Rohlf, 2000) statistical package, and we generated a pairwise similarity matrix with simple matching coefficient according to the SIMQUAL procedure in NTSYS-pc. We performed cluster analysis of the binary matrix using the unweighted pair-group method to build a UPGMA dendrogram by means of SAHN package in NTSYS-pc. We performed bootstrap analysis of the UPGMA tree with 2000 replicates using the winboot computer program (Yap & Nelson, 1996). We also calculated a genetic similarity matrix using AFLP data according to the method of Nei & Li (1979). Dating the divergence between lineages We obtained the nrITS data set of Chloridoideae from Peterson et al. (2016) and used it to estimate the divergence times of Orinus. We verified equal base frequencies among accessions in PAUP*4.0b10 and tested a strict molecular clock hypothesis in MEGA 5 (Tamura et al., 2011). The strict clock was rejected (2lnLR = 434.458, d.f. = 48, P < 0.01). We performed Bayesian divergence time estimation in BEAST 2.4.5 (Bouckaert et al., 2014), which co-estimates topology and node ages. We prepared the nrITS data set for analysis in BEAST using Beauti and by directly editing the content of.xml files. We applied an HKY model of nucleotide substitutions because the more parameter-rich GTR model did not yield substantial gains according to the preliminary analyses. We approximated the gamma distribution of site-specific rates with ten rate categories and applied the ‘birth death’ tree prior process model with a standard, uniform prior. We set a ‘relaxed lognormal clock’ with a substitution rate prior of 6.85 × 10−9 based on Wolfe, Li & Sharp (1987). We calibrated the crown node age of Chloridoideae as 32.74 Mya, which was inferred in prior studies (Christin et al., 2008; Vicentini et al., 2008) and is reasonable because no isotopic surveys to-date provided evidence of an older date for C4 grass (Osborne & Beerling, 2006). Our calibration comprised a normal distribution with an SD of 1.0, yielding a 95% highest posterior density (HPD) of 28.09–37.19 Mya. We constrained the monophyly of the Chloridoideae after preliminary analyses revealed that long branches (i.e. high evolutionary rates) resulted in placement of Sporobolus R.Br. with the outgroups. Sporobolus clearly belongs in Chloridoideae based on strong evidence from prior molecular and morphological studies (Peterson et al., 2016). We ran two independent analyses in BEAST with 5 × 108 Markov chain Monte Carlo generations. We compared the posterior and likelihood distributions from each run in Tracer 1.5.0 to ensure convergence, and we selected one of two runs for downstream analyses (Rambaut & Drummond, 2009). We performed a 10% burn-in and summarized trees by selecting the maximum clade credibility tree in TreeAnnotator 2.4.5, which is part of the BEAST 2.4.5 package, and using median heights for branch lengths. We visualized the maximum clade credibility tree and other summary statistics in FigTree 1.3.1 (Rambaut, 2009). To calculate divergence times of species and populations in Orinus, we used nrITS and two methods: (1) a direct calculation from mutation rates using a strict molecular clock; and (2) a dating analysis in BEAST (Drummond & Rambaut, 2007). For the direct calculation, a strict molecular clock was appropriate, because it was not rejected by an analysis in MEGA 5 (Tamura et al., 2011) (2lnLR = 54.524, d.f. = 29, P = 0.61), and performed a coalescent Bayesian skyline tree process prior suitable for intraspecific data. The strict molecular clock comprised mutation rates of 5.8–8.1 × 10−9 substitutions/site/year, which are well established as rates of evolution of nrITS in other perennial grass genera (Wolfe et al., 1987). For the analyses in BEAST, we analysed 476 sequences of Orinus, and we used four samples of Cleistogenes as outgroups. We applied an age calibration to the Orinus crown node based on our results from the tree for Chloridoideae (above). Specifically, the calibration comprised a uniform prior with using the median date from the tree for Chloridoideae ± 1 Mya. All other parameters of this analysis were the same as for the analysis of Chloridoideae (Peterson et al., 2016), except that we used a coalescent skyline model of diversification because it is appropriate for intraspecific taxa. Glacial refugia and demographic analyses We inferred putative glacial refugia of Orinus based on occurrences of comparatively high genetic diversity among intraspecific populations and for the genus with plastid DNA and nrITS data sets. Genetic diversity is expected to be higher in refugial areas compared with surrounding areas (Hewitt, 1996, 1999, 2000) due especially to founder effects (Mayr, 1942). Although associating high diversity areas with glacial refugia may be overly simplistic (Petit et al., 2003), our intention is to generate a working hypothesis for future testing. To visualize putative refugia based on diversity across geographical areas, we generated surface plots of the number of unique haplotypes, ribotypes and combined haplotypes and ribotypes within populations using a triangulated irregular network (TIN) in ArcMap v.10 (ESRI, 2010). The from-end of link points are used as the TIN uses triangles nodes points, georeferenced populations in this case, and interpolate values, or number of genotypes in this case, between them. We trimmed the TINs according to the convex hulls of O. thoroldii and O. intermedius and for the geographically distinct northern and central metapopulations of O. kokonoricus (Fig. 6A). The resulting surfaces facilitate inference and visualization of geographical areas of high diversity, which we take as possible glacial refugia. We deposited all GIS data online with ESRI, and it is available at https://www.arcgis.com/home/item.html?id=d2f9770771ec4de58a8373d85e3695ce. Figure 6. View largeDownload slide The sampled populations and interpolated distribution of genetic diversity of plastid DNA and nrITS in three Orinus spp. A, the sampled populations and convex hulls drawn around them. B, the interpolated distribution of genetic diversity of plastid DNA. C, the interpolated distribution of genetic diversity of nrITS. D, the highest genetic diversity regions of plastid DNA and nrITS. E, the interpolated distribution of genetic diversity of plastid DNA and nrITS. F, the legend of the study area in China to provide the geographical context, and a–e correspond to the number of genotypes 1–5, respectively. Figure 6. View largeDownload slide The sampled populations and interpolated distribution of genetic diversity of plastid DNA and nrITS in three Orinus spp. A, the sampled populations and convex hulls drawn around them. B, the interpolated distribution of genetic diversity of plastid DNA. C, the interpolated distribution of genetic diversity of nrITS. D, the highest genetic diversity regions of plastid DNA and nrITS. E, the interpolated distribution of genetic diversity of plastid DNA and nrITS. F, the legend of the study area in China to provide the geographical context, and a–e correspond to the number of genotypes 1–5, respectively. We used several methods to detect demographic expansions based on plastid DNA data sets, such as from refugia, in the history of Orinus. We sought to detect recent demographic expansion in Orinus and in each species using a mismatch distribution analysis in Arlequin 3.5. The mismatch distribution analysis applies a model of computes the pairwise p-distances (absolute distances) between sequences and organizes the results in a histogram. A distribution of pairwise frequencies that is multimodal is the equilibrium state, showing expected differentiation among species or populations. In contrast, a unimodal distribution represents geneflow and a trend towards panmixia (Slatkin & Hudson, 1991; Rogers & Harpending, 1992). We analysed the results of the mismatch distribution by comparing the sum of squared deviations observed distributions and those expected following a recent demographic expansion. We also quantified modality of the distribution using Harpending’s raggedness index, where more ragged distributions are multimodal. We assessed recent demographic expansions using Tajima’s D (Tajima, 1989) or Fu’s FS (Fu, 1997) in Arlequin 3.5, with 10000 parametric bootstrap replicates. Positive values for these statistics typically indicate demographic contractions or a recent bottleneck event (Tajima, 1989; Fu, 1997; Hamilton, 2009). In addition, we estimated expansion times (t) using the relationship τ = 2ut (Rogers & Harpending, 1992), where u is the mutation rate per generation for all sequences analysed and t is the expansion time in number of generations. Here, u = μkg, where μ is the substitution rate, k is the average sequence length of the analysed DNA region and g is the generation time in years. We assumed that g = 5 years based on data from the closely related genus Triodia R.Br. (Armstrong, 2011). RESULTS Identification of plastid DNA haplotypes and nrITS ribotypes, and dating divergence times Among the 476 individuals, the total alignment length of three plastid DNA sequences (matK, rbcL and psbA-trnH) was 2898 bp; for nrITS, it was 624 bp. Nucleotide substitutions occurred at 19 polymorphic sites in plastid DNA (0.66%), and 13 of these were potentially parsimony informative; there were no indels (see Supporting Information, Table S2). The aligned nrITS sequences had 46 variable sites, 45 of which were potentially parsimony informative (see Supporting Information, Table S3). We used sequence polymorphisms to identify haplotypes in plastid DNA and ribotypes in nrITS. In plastid DNA, we identified 14 different haplotypes, H1–H14 (Fig. 1B; see Supporting Information, Table S1). Forty-nine populations (55.68%) were fixed for a single haplotype; the remaining 39 (44.32%) were polymorphic (Fig. 1B; see Supporting Information, Table S1). Two haplotypes, H4 and H6, were shared between O. kokonoricus and O. intermedius, and only one, H2, was fixed among all three Orinus spp. Of these shared haplotypes, H2 and H4 occurred in 67.05 and 17.05% of all populations, respectively, whereas H6 was limited to six populations (6.82%). H2 was widely distributed across the QTP, and H4 was widely distributed throughout the eastern and south-eastern QTP. Six haplotypes were found only in O. thoroldii, four were specific to O. kokonoricus and one was specific to O. intermedius (Fig. 1C; see Supporting Information, Fig. S1). In nrITS, we identified 30 different ribotypes, S1–S30. Of these, ten were each limited to a single population, and all but one were restricted to one species (Fig. 2; see Supporting Information, Fig. S2, Table S1). The S10 ribotype was shared between O. kokonoricus and O. intermedius. Thirty-nine populations (44.32%) only possessed one haplotype; the remaining 49 (55.68%) were polymorphic (Fig. 2; see Supporting Information, Fig. S2, Table S1). Ribotypes S2, S19 and S26 were the most common nrITS sequences and were fixed or present in nearly all sampled populations from the three Orinus spp. Sequence S2 was discovered in 84.44% of all samples for O. kokonoricus, S19 was found in all populations from O. thoroldii and S26 also occurred in 86.67% of all sampled populations of O. intermedius. Ten ribotypes (S7–S9, S15–S18, S22, S23 and S25) were private to a single population (Fig. 2; see Supporting Information, Table S1). Figure 2. View largeDownload slide Sampling sites and geographical distribution of nrITS ribotypes observed in three Orinus species. A, sampling region (dashed-line box). B, documented distribution range (dotted line), sampling sites and nrITS ribotype frequencies per population. Each pie graph represents nrITS ribotypes; size of pie is proportional to number of individuals bearing that haplotype. Figure 2. View largeDownload slide Sampling sites and geographical distribution of nrITS ribotypes observed in three Orinus species. A, sampling region (dashed-line box). B, documented distribution range (dotted line), sampling sites and nrITS ribotype frequencies per population. Each pie graph represents nrITS ribotypes; size of pie is proportional to number of individuals bearing that haplotype. The divergence time estimation showed Orinus diverged from its closest relatives c. 12.06 Mya (95% HPD: 6.87–17.50) (Miocene) and diversified c. 2.85 Mya (node 1; 95% HPD: 0.78–5.05) in the Pliocene (Fig. 3; see Supporting Information, Table S4). In Orinus, modern ribotypes diversified in the Pleistocene, 0.31–1.71 Mya (nodes 2–9) (see Supporting Information, Fig. S2, Table S4), according to the analysis using BEAST. When the divergence time of ribotypes was calculated using mutation rates, we found the similar results as the above, also in the Pleistocene (see Supporting Information, Table S4). Figure 3. View largeDownload slide Bayesian divergence time estimate of Orinus (shaded) based on the nrITS data set of Chloridoideae from Peterson et al. (2016). Clade constraint is marked by the black circle. Numbers in parentheses indicate 95% highest posterior density intervals. Figure 3. View largeDownload slide Bayesian divergence time estimate of Orinus (shaded) based on the nrITS data set of Chloridoideae from Peterson et al. (2016). Clade constraint is marked by the black circle. Numbers in parentheses indicate 95% highest posterior density intervals. Population variability, glacial refugia and demography Haplotype diversity (HE) from each population was calculated based on plastid DNA haplotype frequencies (see Supporting Information, Table S1). The average genetic diversity within populations (HS) was relatively low for O. thoroldii, O. kokonoricus and O. intermedius (0.16 ± 0.05, 0.23 ± 0.07 and 0.31 ± 0.10), and the total genetic diversity HT (0.74 ± 0.03) across the whole genus was much higher than HS (0.21 ± 0.04) based on the observed plastid DNA variations (Table 1). A permutation test suggested that population differentiation across the sampling range was high (NST = 0.83 ± 0.04) and was significantly higher than GST (0.72 ± 0.05) (0.01 < P < 0.05), indicating a strong phylogeographic structure across the entire sampling range (Table 1). Table 1. Estimates of genetic diversity for plastid DNA in three Orinus spp. on the Qinghai-Tibet Plateau Groups  HS  HT  GST  NST  O. kokonoricus (P1–45)  0.23 (0.07)  0.74 (0.06)  0.69 (0.10)  0.60 (0.12)  O. thoroldii (P46–73)  0.16 (0.05)  0.56 (0.04)  0.72 (0.08)  0.81 (0.07)*  O. intermedius (P74–88)  0.31 (0.10)  0.59 (0.05)  0.48 (0.18)  0.33 (0.15)  All populations  0.21 (0.04)  0.74 (0.03)  0.72 (0.05)  0.83 (0.04)*  Groups  HS  HT  GST  NST  O. kokonoricus (P1–45)  0.23 (0.07)  0.74 (0.06)  0.69 (0.10)  0.60 (0.12)  O. thoroldii (P46–73)  0.16 (0.05)  0.56 (0.04)  0.72 (0.08)  0.81 (0.07)*  O. intermedius (P74–88)  0.31 (0.10)  0.59 (0.05)  0.48 (0.18)  0.33 (0.15)  All populations  0.21 (0.04)  0.74 (0.03)  0.72 (0.05)  0.83 (0.04)*  HS, estimate of average genetic diversity within populations; HT, total genetic diversity; GST, coefficient of genetic variation over all populations; NST, coefficient of genetic variation influenced by both haplotype frequencies and genetic distance between haplotypes. *P < 0.05. View Large Hierarchical AMOVA based on the plastid DNA data set revealed that there was higher interspecific genetic differentiation across the whole distribution of this genus (FST = 0.82), and the interpopulation genetic differentiation in O. thoroldii and O. kokonoricus was much more pronounced (FST = 0.72, FST = 0.73), but there was a little interpopulation differentiation in O. intermedius (FST = 0.34) (Table 2). These results suggested that the grouping pattern of hierarchical AMOVA was basically consistent with relationships among plastid haplotypes and their geographical distribution (Fig. 1C, Table 2). Table 2. Results of analysis of molecular variance for plastid DNA and nrITS sequence data from populations of Orinus on the Qinghai-Tibet Plateau Groups  Source of variation  d.f.  SS  VC  Variation (%)  Fixation index    Plastid DNA  All samples  Among groups  2  121.90  0.74 Va  44.78  FST = 0.82*  Among populations within groups  53  152.13  0.62 Vb  37.72  FSC = 0.68*  Within populations  420  52.46  0.29 Vc  17.50  FCT = 0.45*  Total  475  326.49  1.65      O. kokonoricus  Among populations  20  50.09  0.54 Va  72.69  FST = 0.73*  Within populations  161  14.17  0.20 Vb  27.31    Total  181  64.26  0.74      O. thoroldii  Among populations  26  92.52  0.92 Va  71.90  FST = 0.72*  Within populations  161  24.18  0.36 Vb  28.10    Total  187  116.70  1.28      O. intermedius  Among populations  7  9.52  0.16 Va  34.27  FST = 0.34*  Within populations  98  14.10  0.32 Vb  65.73    Total  105  23.62  0.48        nrITS  All samples  Among groups  2  1977.08  6.32 Va  81.90  FST = 0.96*  Among populations within groups  53  491.63  1.08 Vb  13.96  FSC = 0.77*  Within populations  420  134.40  0.32 Vc  4.15  FCT = 0.82*  Total  475  2603.11  7.72      O. kokonoricus  Among populations  20  204.28  1.13 Va  67.89  FST = 0.68*  Within populations  161  86.25  0.54 Vb  32.11    Total  181  290.53  1.67      O. thoroldii  Among populations  26  21.58  0.09 Va  30.17  FST = 0.30*  Within populations  161  33.48  0.21 Vb  69.83    Total  187  55.06  0.30      O. intermedius  Among populations  7  265.77  2.95 Va  95.17  FST = 0.95*  Within populations  98  14.67  0.15 Vb  4.83    Total  105  280.43  3.10      Groups  Source of variation  d.f.  SS  VC  Variation (%)  Fixation index    Plastid DNA  All samples  Among groups  2  121.90  0.74 Va  44.78  FST = 0.82*  Among populations within groups  53  152.13  0.62 Vb  37.72  FSC = 0.68*  Within populations  420  52.46  0.29 Vc  17.50  FCT = 0.45*  Total  475  326.49  1.65      O. kokonoricus  Among populations  20  50.09  0.54 Va  72.69  FST = 0.73*  Within populations  161  14.17  0.20 Vb  27.31    Total  181  64.26  0.74      O. thoroldii  Among populations  26  92.52  0.92 Va  71.90  FST = 0.72*  Within populations  161  24.18  0.36 Vb  28.10    Total  187  116.70  1.28      O. intermedius  Among populations  7  9.52  0.16 Va  34.27  FST = 0.34*  Within populations  98  14.10  0.32 Vb  65.73    Total  105  23.62  0.48        nrITS  All samples  Among groups  2  1977.08  6.32 Va  81.90  FST = 0.96*  Among populations within groups  53  491.63  1.08 Vb  13.96  FSC = 0.77*  Within populations  420  134.40  0.32 Vc  4.15  FCT = 0.82*  Total  475  2603.11  7.72      O. kokonoricus  Among populations  20  204.28  1.13 Va  67.89  FST = 0.68*  Within populations  161  86.25  0.54 Vb  32.11    Total  181  290.53  1.67      O. thoroldii  Among populations  26  21.58  0.09 Va  30.17  FST = 0.30*  Within populations  161  33.48  0.21 Vb  69.83    Total  187  55.06  0.30      O. intermedius  Among populations  7  265.77  2.95 Va  95.17  FST = 0.95*  Within populations  98  14.67  0.15 Vb  4.83    Total  105  280.43  3.10      d.f., degrees of freedom; SS, sum of squares; VC, variance components; FCT, correlation within populations relative to groups; FSC, correlation of haplotypes within groups relative to the total; FST, correlation within populations relative to the total. *P < 0.001 (10000 permutations). View Large Our analyses revealed extensive gene flow among populations of Orinus across its geographical range. In the STRUCTURE analyses (Fig. 4A), the highest likelihood of the plastid DNA data was obtained when all samples were weakly clustered into three groups according to the recognized species (K = 3). However, this result did not correspond to separate geographical regions because there was strong gene flow across the range of Orinus. In addition, a Nem value of 0.61 also indicated extensive gene flow existed between groups based on plastid DNA data. Figure 4. View largeDownload slide Estimated genetic structure for K = 2–3 obtained with the program STRUCTURE (Pritchard et al., 2000) for 56 populations of Orinus based on plastid DNA (A) and nrITS (B) sequence variations. Each vertical bar represents an individual and its assignment proportion (not shown). Red, Orinus thoroldii; Green, O. kokonoricus; Blue-purple, O. intermedius. Figure 4. View largeDownload slide Estimated genetic structure for K = 2–3 obtained with the program STRUCTURE (Pritchard et al., 2000) for 56 populations of Orinus based on plastid DNA (A) and nrITS (B) sequence variations. Each vertical bar represents an individual and its assignment proportion (not shown). Red, Orinus thoroldii; Green, O. kokonoricus; Blue-purple, O. intermedius. For nrITS, we found that 81.90% of the total genetic variation occurred between groups, 13.96% occurred among populations within groups and only 4.15% occurred within populations. The pairwise FST value between three groups was 0.96 (Table 2). It suggested that there was high interspecific differentiation among three species, and interpopulation differentiation of O. kokonoricus and O. intermedius was also extremely high (FST = 0.68, FST = 0.95). In contrast, O. thoroldii exhibited lower interpopulation variation (FST = 0.30). The grouping pattern of hierarchical AMOVA based on nrITS data was in accordance with relationships of ribotypes and their geographical distribution (Fig. 5, Table 2). Similarly, the STRUCTURE analyses indicated that almost all samples of O. intermedius and O. kokonoricus were clustered into a lineage, and there was only obscure gene flow between them and O. thoroldii (K = 2) (Fig. 4B). However, all samples of this genus were clustered into three independent lineages when the highest likelihood of the data was obtained (K = 3), and there was prominent introgression between O. kokonoricus and O. intermedius (Fig. 4B). In addition, the gene flow between the three groups was low according to the Nem value was 0.11. Figure 5. View largeDownload slide Hypothesis of incomplete lineage sorting and reticulation in nrITS of Orinus shown within a species tree of the genus. Lineage A is fixed in O. thoroldii, whereas descendent types from lineage B are incompletely sorted in O. kokonoricus and O. intermedius. nrITS type S10 evolved in lineage F of O. kokonoricus and was later shared with O. intermedius during a secondary contact event. Black nodes indicate divergence events among types that we can infer from our data, and extant types sampled in this study are numbered 1–30 at terminal nodes. Two particularly large clades of nrITS types in lineage 4 are collapsed into single terminal nodes, and all types represented by those nodes are listed numerically. Figure 5. View largeDownload slide Hypothesis of incomplete lineage sorting and reticulation in nrITS of Orinus shown within a species tree of the genus. Lineage A is fixed in O. thoroldii, whereas descendent types from lineage B are incompletely sorted in O. kokonoricus and O. intermedius. nrITS type S10 evolved in lineage F of O. kokonoricus and was later shared with O. intermedius during a secondary contact event. Black nodes indicate divergence events among types that we can infer from our data, and extant types sampled in this study are numbered 1–30 at terminal nodes. Two particularly large clades of nrITS types in lineage 4 are collapsed into single terminal nodes, and all types represented by those nodes are listed numerically. Our surface analysis in ArcMap revealed four centres of diversity for haplotypes, ribotypes and the combined diversity of haplotypes and ribotypes (Fig. 6). For O. thoroldii, the haplotype and ribotype data are congruent in showing a centre of diversity in the eastern part of the range (Fig. 6B, C). The haplotype and ribotype data also agree that the northern metapopulation of O. kokonoricus has a centre of diversity at the southern portion of its range (Fig. 6B, C). In contrast, the haplotypes show centres of diversity at the northern ends of the ranges of O. intermedius and the central metapopulation of O. kokonoricus, whereas the ribotypes show centres of diversity at the southern ends of these ranges (Fig. 6B, C). Overall, greater diversity among ribotypes than haplotypes (Figs 1, 2) yielded a geographical pattern in the combined analysis of genotypes (Fig. 6E) consistent with ribotypes alone. The mismatch distributions of pairwise differences for all samples from the identified haplotypes of Orinus and O. thoroldii were close to bimodal, whereas those of O. kokonoricus and O. intermedius were close to be unimodal (see Supporting Information, Fig. S3), indicating that the former should experience recent bottlenecks and the latter should experience recent and rapid population expansion. However, further analyses using the non-significant variance and raggedness index tests indicated that none of the observed distributions differed significantly from those expected under a sudden expansion model (Table 3). Moreover, the significant negative calculated values for Tajima’s D (P < 0.01) and non-significant positive calculated values for Fu’s FS (P > 0.05) were also supported this earlier, rapid demographic or geographical expansion and subsequent population growth (Table 3). Besides, the mean expansion times were estimated to be c. 25, 29, 40 and 53 kya for O. intermedius, O. kokonoricus, O. thoroldii and Orinus, respectively (Table 3). These results support a recent and rapid demographic population growth in the genus and each of the three species, and all expansions occurred before the Last Glacial Maximum (LGM; c. 20 kya). Possibly, Orinus and O. thoroldii have experienced recent bottlenecks. Table 3. Results of mismatch distribution analyses and neutrality tests (Tajima’s D and Fu’s FS) based on cpDNA data   Mismatch distribution  Neutrality tests  Groups  τ (95% CI)  SSD (P value)  RAG (P value)  tmin–tmax (kya)  tave (kya)  Tajima’s D (P value)  Fu’s Fs (P value)  O. thoroldii  0.772 (0.046–6.021)  0.057 (0.079)  0.160 (0.208)  26.639–88.797  40.983  −0.239 (0.001)  0.393 (0.000)  O. kokonoricus  0.554 (0.050–1.261)  0.050 (0.099)  0.170 (0.217)  19.117–63.722  29.410  −0.197 (0.004)  0.342 (0.006)  O. intermedius  0.480 (0.047–1.331)  0.025 (0.357)  0.196 (0.513)  16.563–55.210  25.482  −0.247 (0.006)  0.465 (0.014)  Total  0.994 (0.050–3.577)  0.050 (0.128)  0.169 (0.254)  34.300–114.332  52.768  −0.224 (0.001)  0.384 (0.002)    Mismatch distribution  Neutrality tests  Groups  τ (95% CI)  SSD (P value)  RAG (P value)  tmin–tmax (kya)  tave (kya)  Tajima’s D (P value)  Fu’s Fs (P value)  O. thoroldii  0.772 (0.046–6.021)  0.057 (0.079)  0.160 (0.208)  26.639–88.797  40.983  −0.239 (0.001)  0.393 (0.000)  O. kokonoricus  0.554 (0.050–1.261)  0.050 (0.099)  0.170 (0.217)  19.117–63.722  29.410  −0.197 (0.004)  0.342 (0.006)  O. intermedius  0.480 (0.047–1.331)  0.025 (0.357)  0.196 (0.513)  16.563–55.210  25.482  −0.247 (0.006)  0.465 (0.014)  Total  0.994 (0.050–3.577)  0.050 (0.128)  0.169 (0.254)  34.300–114.332  52.768  −0.224 (0.001)  0.384 (0.002)  τ, time in number of generations elapsed since the sudden expansion episode; SSD, sum of squared deviations; RAG, Harpending’s raggedness index; tmin, absolute expansion time estimated with a high substitution rate, 1.0 × 10−9 substitutions/site/year (s/s/year); tmax, absolute expansion time estimated with a low substitution rate, 0.3 × 10−9 s/s/year; tave, mean expansion time; kya, thousand years ago; CI, confidence interval. View Large DISCUSSION Phylogeographic patterns in Orinus In recent years, many studies have focused on the early stages of speciation (Abbott et al., 2013; Liu et al., 2014b). Some of these have shown that divergence among species can occur even in the presence of gene flow and that species boundaries are maintained by divergent selection (Hewitt, 1996; Petit & Excoffier, 2009; Schluter, 2009; Nosil, 2012; Sousa & Hey, 2013). Our plastid DNA, nrITS and AFLP data suggest that Orinus may be at an advanced stage of speciation in which gene flow is limited, by geographical isolation or divergent selection, and fixation is underway so that most sequences form clades comprising only one species, but these species are paraphyletic (Fig. 5) (Funk & Omland, 2003) Population bottlenecks have also played a role in promoting divergence among closely related taxa but has rarely been studied (Butlin et al., 2012; Bi et al., 2016). Our results suggest a subtle but complex geographical pattern in the genetic diversity of Orinus. High levels of interspecific genetic variation occurred among samples of plastid DNA (44.78%) and nrITS (81.90%) (Table 2). In addition, the differentiation among species was high overall evidenced by significant FST of 0.86 and 0.96 for plastid DNA and nrITS (Table 2), indicating nearly complete genetic isolation (Rieseberg, Church & Morjan, 2004; Hauquier et al., 2017). Moreover, only one plastid DNA haplotype and no nrITS types were shared among all populations (Figs 1, 2). However, the pattern is more complex than isolation by species; the eastern and south-eastern species, O. kokonoricus and O. intermedius, clearly have more gene flow between them than either does with the western species, O. thoroldii. Orinus kokonoricus and O. intermedius share two additional plastid DNA haplotypes in common (Fig. 1), exclusive of O. thoroldii, and one nrITS type (Fig. 2), strongly suggestive of recent or ongoing gene flow (Parchman, Benkman & Britch, 2006). Moreover, the phylogenetic reconstruction of nrITS sequences in O. kokonoricus and O. intermedius shows that the haplotypes unique to each species do not form clades; they are intermixed in their phylogenetic positions, and this is highly consistent with genetic connectivity and subsequent isolation and fixation. In contrast, O. thoroldii does not share any additional haplotypes or nrITS ribotypes in common with O. kokonoricus or O. intermedius. Overall, the demographic pattern suggests greater and/or more frequent connectivity between the eastern and south-eastern species than with the western one. It is noteworthy that Orinus spp. are largely differentiated and possess few shared ribotypes or haplotypes, but the ribotypes or haplotypes do not form the same clades as the species (Fig. 5; see Supporting Information, Figs S1, S2). At least two, non-mutually exclusive hypotheses may explain this phylogeographic pattern: (1) somewhat recent species divergence with incomplete lineage sorting (Choleva et al., 2014; Zhou et al., 2017) or (2) recent and/or repeated secondary contact and introgression (Avise, 2000; Sardell & Uy, 2016). In Orinus, both mechanisms are needed to explain our results (Figs 1, 2, 4). Based on our nrITS data, we propose nrITS split into two lineages, A and B before or after the divergence of O. thoroldii from the rest of Orinus (Fig. 5). Lineage A became fixed in O. thoroldii, in which lineage sorting is complete (Fig. 5). Lineage B diversified into lineages C, D, E and F in an ancestor of O. kokonoricus and O. intermedius (Fig. 5). Lineages C and E are limited to O. intermedius, whereas D and F occur only in O. kokonoricus (Fig. 5). Thus, sorting of the descendent genotypes of lineage B is incomplete between O. kokonoricus and O. intermedius. The S10 nrITS type is shared by both O. kokonoricus and O. intermedius, but is clearly derived in lineage F (Fig. 5). Therefore, it seems probable that the S10 type was shared by O. kokonoricus with O. intermedius during a secondary contact and hybridization event (Sardell & Uy, 2016; Zhou et al., 2017). The plastid DNA data and the AFLPs are not in conflict with this hypothesis. Plastid DNA loci sort and become fixed more slowly than nuclear markers (Wolfe et al., 1987; de Villiers et al., 2013; Shaw et al., 2017), such as nrITS, and the shared haplotypes among all three species and between O. kokonoricus and O. intermedius may represent incomplete lineage sorting (Choleva et al., 2014; Yasuda et al., 2015; Kuritzin et al., 2016). Alternatively, the haplotype shared among all three species may represent an older secondary contact (Miraldo et al., 2011; Sardell & Uy, 2016), predating the S10 event and too old to be preserved in nrITS data. AFLPs revealed that species and populations within species comprise clades. Therefore, most but not all lineage sorting among species is complete (Choleva et al., 2014; Yasuda et al., 2015; Kuritzin et al., 2016). The demographic pattern that we observed in Orinus seems highly consistent with the history of geomorphism in the QTP and adjacent areas. The evidence for QTP uplifts comes primarily from detecting its effects on regional aridity and the monsoon climatic using dust sedimentation and indicator fossils of marine foraminifera (e.g. An et al., 2001; Sun & An, 2001; Guo et al., 2002; Zheng et al., 2004). Although the exact timing of uplift events remains somewhat in question, the QTP probably underwent early orogenesis c. 40 Mya along its western margin as the Indian subcontinent collided with Eurasia and more extensive, rapid uplift 20–7 Mya (Harrison et al., 1992; Ruddiman, 1998). In addition, many geological and historical biogeographic studies infer that that additional uplifts occurred east of the QTP around c. 4.0–2.5 Mya in adjacent mountain ranges, especially the Hengduan Mountains (Chen, 1992, 1996; Zheng et al., 2004; Gao et al., 2013; Favre et al., 2015). Thus, the Hengduan Mountains and other adjacent ranges have experienced more recent, Pliocene geological changes and associated climatic perturbations than the QTP proper (Xing & Ree, 2017). Contemporaneously, glacial cycling also impacted the area, and its effects can be difficult to distinguish from the genetic signature of geomorphism (Shackleton & Opdyke, 1973; Liu et al., 2014b). Nevertheless, it is possible that regional geomorphism to the east of the QTP facilitated the dynamic history of recent allopatric speciation and secondary contact suggested by our nrITS data for O. kokonoricus and O. intermedius. Similarly, the relatively recent geomorphism to the east of the QTP by comparison to the QTP may also be reflected in the evolutionary and demographic histories of other diverse plant groups. For example, Yang et al. (2012) found that species of Meconopsis Vig. (Papaveraceae) occurring in the QTP diverged earlier than species occurring in adjacent eastern areas. Allopatric divergence, glacial refugia and historical demography The role of climatic oscillations in shaping the regional biota is important and should not be ruled out as a driver of evolutionary events in Orinus. Other species in the region have histories affected by climatic change. For example, Anisodus tanguticus Pascher (Solanaceae) has western and eastern clades that diverged during the Pleistocene, probably a result of genetic isolation in glacial refugia, rather than a geological event (Wan et al., 2016). Similar patterns of intraspecific divergence due to Quaternary climatic change have also been found in other plant species in the QTP and adjacent regions, such as Aconitum gymnandrum Maxim. (Ranunculaceae) (Wang et al., 2009a), Sophora davidii (Franch.) Skeels (Fabaceae) (Fan et al., 2013) and Oxyria sinensis Hemsl. (Polygonaceae) (Meng et al., 2015). Based on the divergence time of 32.74 Mya for Chloridoideae, the crown age was dated to be 2.85 (95% HPD: 0.58–12.45) Mya, with interspecific divergences occurred in the range of 0.31–1.70 Mya (see Supporting Information, Fig. S2, Table S4). This relatively recent interspecific diversification atop an unbranched stem lineage is sometimes regarded as reflecting ancient extinctions (Harvey, May & Nee, 1994). Deep divergence in Orinus apparently began in the mid-Pleistocene (see Supporting Information, Fig. S2, Table S4). Although these estimated timescales for deep divergence must be interpreted with caution, they (the major nodes) correspond relatively well with recent QTP uplifts (0.01–1.80 Mya) (Li, Shi & Li, 1995; Shi, Li & Li, 1998; Zheng, Xu & Shen, 2002). Development of aridification in the QTP was a response to climatic change at the time, such as long-term cooling and drying trends. By the early Pleistocene, the QTP dramatically uplifted to present heights (Wu et al., 2008), which made the climate of this area become extremely cool and dry. Especially, strong episodic cooling resulted in sharp increases in aridity of the QTP regions (Fang et al., 2002). Aridification played a significant role in the increase of genetic diversity, and even species divergence and speciation of alpine plants. Habitat fragmentation, vicariance and population isolation on the QTP provided opportunities for allopatric speciation through the action of selection and/or genetic drift (Wang et al., 2013; Meng et al., 2014; Wen et al., 2014). Thus, deep intraspecific divergence in Orinus may be a result of adaptation to arid habitats, with the mechanisms probably being similar to those in other arid/desert species (Wang et al., 2013; Yu et al., 2013; Liu et al., 2014b; Meng et al., 2014). Therefore, it is likely that the deep divergence in Orinus was caused by allopatric differentiation and reduction in gene flow following plateau uplifts and climatic oscillation. Plant populations in refugia generally have high genetic divergence (Petit et al., 2003). Our surface analysis in ArcMap revealed four centres of diversity for haplotypes, ribotypes and the combined diversity of haplotypes and ribotypes, which were consistent with the result of glacial refugia based on both plastid DNA and nrITS data sets (Figs 1B, 2B, 6B–E). Among them, there were three areas with relatively high plastid DNA and nrITS genetic diversity and one region with unique plastid DNA haplotypes and high nrITS genotypes (Figs 1B, 2B, 6B–E). One region with plastid DNA haplotype uniqueness and high nrITS genetic diversity (O. intermedius) is located on the southern edge of the QTP, an area well known as an important refugium for other alpine plant species (Figs 1B, 2B, 6E) (Zhang et al., 2005; Meng et al., 2007; Chen et al., 2008; Yang et al., 2008; Wang et al., 2009a; Gao et al., 2012; Liu et al., 2012, 2014b; Ma et al., 2014), which was further supported by the result of surface analysis in ArcMap (Fig. 6C–E). Furthermore, three (S28–S30) of four (S10, S28–S30) nrITS ribotypes from this region were endemic, indicating a divergence centre for Orinus. Plastid and nrITS genetic diversities are higher in two regions (O. kokonoricus) in the north-eastern and eastern plateau (Figs 1B, 2B, 6E), which may represent two important refugia throughout the Quaternary. Our surface analysis in ArcMap also indicated the northern metapopulation of O. kokonoricus has a centre of diversity at the southern portion of its range (Fig. 6B, C), and the southern metapopulation of O. kokonoricus has another centre of diversity at the southern portion of its range (Fig. 6E). We have two alternative biogeographic hypotheses regarding high levels of genetic diversity in O. kokonoricus. A glacial refugium allowed O. kokonoricus to survive climatic oscillations and to accumulate genetic diversity (Tzedakis et al., 2002) or this region was composed of a mixture of individuals with different origins, resulting in higher genetic diversity than the original source (Petit et al., 2003). If the latter hypothesis is true, then haplotypes detected in this region should form a subset of haplotypes of the multiple original sources, although four (H1, H3, H5 and H7) out of the seven plastid DNA haplotypes in these populations of O. kokonoricus are endemic. Furthermore, the majority of nrITS ribotypes (S1–S9 and S11–S18) were restricted to this area, indicating other divergence centres of O. kokonoricus. In addition, widespread haplotypes (H2 and H4) and ribotypes (S2 and S3) coexisting among the great number of plastid DNA haplotypes and nrITS ribotypes also occurred in these two regions. Thus, the most-parsimonious explanation is consequently the presence of glacial refugia for O. kokonoricus. A similar result was found in previous studies of endemics species from the QTP (e.g. Gao et al., 2012; Ma et al., 2014). In addition, the western region of the QTP was only occupied by O. thoroldii probably represent the fourth separate refugia during the large glaciation stage (Figs 1, 2, 6B–E). The populations of O. thoroldii in eastern part of the range had relatively high plastid DNA and nrITS genetic diversity (Figs 1, 2), which were congruent in showing a centre of diversity according to the surface analysis in ArcMap (Fig. 6B–E). Therefore, we thought greater diversity among ribotypes than haplotypes yielded a geographical pattern in the combined analysis of genotypes consistent with ribotypes alone. Based on the plastid DNA data set, the three Orinus spp. might have experienced recent and repeated population expansion after allopatric divergence (e.g. Avise, 2004; Wang et al., 2009b; Gao et al., 2012; Ma et al., 2014; Liu et al., 2015), and O. thoroldii might also have experienced recent bottlenecks due to its bimodal mismatch distribution (Fig. S3B). A star-like phylogeny pattern for O. thoroldii and O. kokonoricus is shown in the haplotype network (Fig. 1C), suggesting the haplotypes might have originated from a sudden expansion event (Hudson, 1990). This is supported by significantly negative Tajima’s D values (Table 3), as well as unimodal mismatch distributions for O. kokonoricus and O. intermedius (see Supporting Information, Fig. S3C, D). Overall, these findings support expansions of Orinus in western, eastern and north-eastern regions, perhaps between 25.48 and 52.77 kya, prior to the LGM (18–24 kya) (Shi et al., 1998) (Table 3), but later than other previously reported QTP endemic genera (Ma et al., 2014) and species (e.g. Yang et al., 2008; Jia et al., 2011; Li et al., 2012). Early expansion of Orinus would have yielded all haplotypes in the QTP region. However, some haplotypes are fixed in specific populations. For example, H3, H4 and H6 correspond to populations 3–5, 42–43 and 17–20, respectively. Thus, the existence of multiple high-elevation refugia that would have facilitated survival during the LGM (Wang et al., 2009c; Opgenoorth et al., 2010; Jia et al., 2011; Gao et al., 2012; Liu et al., 2012; Ma et al., 2014). Perhaps, only small and isolated populations, such as P87 and P88, survived in restricted ice-free regions, and range expansion in the refugia appears to be restricted during glacial and interglacial periods. However, another distinct genetic signature via large-scale range expansions, that is wide distribution of a single haplotype (Avise, 1987; Comes & Kadereit, 1998; Hewitt, 2000). In the present study, two widely distributed haplotypes (H2, H9) were found in Orinus. For instance, H2 was found as a single haplotype in 59 of the 88 populations of the QTP and fixed in 25 of them, whereas H9 was fixed in 11 of the 28 populations in O. thoroldii from the western QTP (e.g. Jia et al., 2011; Li et al., 2012; Wang et al., 2014; Liu et al., 2015). Therefore, Orinus might have experienced a second recent range expansion, probably after the LGM, which led to the wide distribution of the two haplotypes in the QTP region. Orinus intermedius is a newly described species (Su et al., 2017) and is sister to O. kokonoricus (Figs 1C, 3, 5). The morphological characters of O. intermedius are intermediate between its two congeners (Su et al., 2017). We initially hypothesized that O. intermedius may have formed via hybrid speciation on the QTP and its adjacent areas (Su et al., 2015, 2017) in a similar manner to Ostryopsis intermedia B.Tian & J.Q.Liu (Betulaceae) (Liu et al., 2014a), Picea purpurea Masters (Pinaceae) (Sun et al., 2014), Roscoea humeana Balf.f. & W.W.Sm. and R. cautleoides Gagnep. (Zingiberaceae) (Zhao et al., 2016). However, our plastid and nrITS data support recent divergences between O. intermedius and O. kokonoricus, which may explain their morphological similarities. In addition, the AFLP data suggest that the three Orinus spp. each form a separate clade indicating phylogenetic distinctiveness (Supporting Information, Fig. S4). As the next steps, we plan to use more rapidly evolving markers, such as single- or low-copy nuclear genes, or even next- generation sequencing approaches in Orinus (Zimmer & Wen, 2012, 2015) to unravel the genetic structure and phylogeographic patterns and to explore speciation mechanisms (Abbott, 2017; Crawford & Archibald, 2017). SUPPORTING INFORMATION Additional Supporting Information may be found in the online version of this article at the publisher’s web-site: Figure S1. Strict MP consensus tree of 14 plastid DNA haplotypes of three Orinus spp. Each tree is topologically congruent with ML tree and Bayesian consensus tree based on same data set. Bootstrap support values from Bayesian and MP analyses, and corresponding posterior probabilities from MP analyses, are given above branches (values > 50%). Dinebra viscida was used as an outgroup. Figure S2. Strict MP consensus tree of 30 nrITS ribotypes of three Orinus species. Each tree is topologically congruent with ML tree and Bayesian consensus tree based on same data set. Bootstrap support values from Bayesian and ML analyses, and corresponding posterior probabilities from MP analyses, are given above branches (values > 50%). Four Cleistogenes spp. (C. bulgarica, C. festucacea, C. mucronata and C. songorica) were used as outgroups. Figure S3. Mismatch distributions for plastid DNA data in genus Orinus (A), O. thoroldii (B), O. kokonoricus (C) and O. intermedius (D). The thin line represents the expected (Exp) values, and the dashed line represents the observed (Obs) values. Figure S4. UPGMA dendrograms for three Orinus spp. generated by cluster analysis based on AFLP data set. Table S1. Sampling data, estimates of haplotype diversity (HE), and haplotype composition for plastid DNA and nrITS data sets from 88 populations of Orinus. Table S2. Variable sites for 14 plastid DNA haplotypes of three Orinus spp. on the Qinghai-Tibet Plateau (QTP). Table S3. Variable sites for 30 nrITS ribotypes of three Orinus spp. on the Qinghai-Tibet Plateau (QTP). Table S4. Estimates of divergence times in Orinus for the major nodes based on nrITS data set. ACKNOWLEDGEMENTS We are grateful to Prof. Jianquan Liu for providing space and expendables in the laboratory and two anonymous reviewers for their comments on the manuscript. This work was supported by grants from the National Natural Science Foundation of China (31260052 and 41761009), the Natural Science Foundation of Qinghai Province (2017-ZJ-904, 2016-ZJ-937Q and 2014-ZJ-947Q) and the State Scholarship Fund of China Scholarship Council (201508635003). REFERENCES Abbott RJ. 2017. Plant speciation across environmental gradients and the occurrence and nature of hybrid zones. Journal of Systematics and Evolution  55: 238– 258. Google Scholar CrossRef Search ADS   Abbott R, Albach D, Ansell S, Arntzen JW, Baird SJ, Bierne N, Boughman J, Brelsford A, Buerkle CA, Buggs R, Butlin RK, Dieckmann U, Eroukhmanoff F, Grill A, Cahan SH, Hermansen JS, Hewitt G, Hudson AG, Jiggins C, Jones J, Keller B, Marczewski T, Mallet J, Martinez-Rodriguez P, Möst M, Mullen S, Nichols R, Nolte AW, Parisod C, Pfennig K, Rice AM, Ritchie MG, Seifert B, Smadja CM, Stelkens R, Szymura JM, Väinölä R, Wolf JB, Zinner D. 2013. Hybridization and speciation. Journal of Evolutionary Biology  26: 229– 246. Google Scholar CrossRef Search ADS PubMed  Abbott RJ, Smith LC, Milne RI, Crawford RM, Wolff K, Balfour J. 2000. Molecular analysis of plant migration and refugia in the Arctic. Science  289: 1343– 1346. Google Scholar CrossRef Search ADS PubMed  An ZS, Kutzbach JE, Prell WL, Porter SC. 2001. Evolution of Asian monsoons and phased uplift of the Himalaya-Tibetan plateau since late Miocene times. Nature  411: 62– 66. Google Scholar CrossRef Search ADS PubMed  Armstrong G. 2011. Evidence for the equal resilience of Triodia spp. (Poaceae), from different functional groups, to frequent fire dating back to the late Pleistocene. Heredity  107: 558– 564. Google Scholar CrossRef Search ADS PubMed  Avise JC. 1987. Identification and interpretation of mitochondrial DNA stocks in marine species. In: Kumpf H, Nakamura EL, eds. Proceedings of the Stock Identification Workshop . Panama City: National Oceanographic and Atmospheric Administration, 105– 136. Avise JC. 2000. Phylogeography: the history and formation of species . Cambridge: Harvard University Press. Avise JC. 2004. Molecular markers, natural history and evolution, 2nd edn . Sunderland: Sinauer Associates. Bawa KS, Koh LP, Lee TM, Liu J, Ramakrishnan PS, Yu DW, Zhang YP, Raven PH. 2010. Ecology. China, India, and the environment. Science  327: 1457– 1459. Google Scholar CrossRef Search ADS PubMed  Bi H, Yue W, Wang X, Zou J, Li L, Liu J, Sun Y. 2016. Late Pleistocene climate change promoted divergence between Picea asperata and P. crassifolia on the Qinghai-Tibet Plateau through recent bottlenecks. Ecology and Evolution  6: 4435– 4444. Google Scholar CrossRef Search ADS PubMed  Bouckaert R, Heled J, Kühnert D, Vaughan T, Wu CH, Xie D, Suchard MA, Rambaut A, Drummond AJ. 2014. BEAST 2: a software platform for Bayesian evolutionary analysis. PLoS Computational Biology  10: e1003537. Google Scholar CrossRef Search ADS PubMed  Butlin R, Debelle A, Kerth C, Snook RR, Beukeboom LW, Castillo CRF, Diao W, Maan ME, Paolucci S, Weissing FJ, van de Zande L, Hoikkala A, Geuverink E, Jennings J, Kankare M, Knott KE, Tyukmaeva VI, Zoumadakis C, Ritchie MG, Barker D, Immonen E, Kirkpatrick M, Noor M, Macias Garcia C, Schmitt T, Schilthuizen M. 2012. What do we need to know about speciation? Trends in Ecology & Evolution  27: 27– 39. Google Scholar CrossRef Search ADS PubMed  Chen FB. 1992. H-D event: an important tectonic event of the late Cenozoic in Eastern Asia. Mountain Research  10: 195– 202. Chen FB. 1996. Second discussion on the H-D movement. Mountain Research  17: 14– 22. Chen K, Abbott RJ, Milne RI, Tian XM, Liu J. 2008. Phylogeography of Pinus tabulaeformis Carr. (Pinaceae), a dominant species of coniferous forest in northern China. Molecular Ecology  17: 4276– 4288. Google Scholar CrossRef Search ADS PubMed  Chen SL, Phillips SM. 2006. Orinus (Poaceae). In: Wu ZY, Raven PH, ed. Flora of China, Vol. 22 . Beijing/St. Louis: Science Press/Missouri Botanical Garden, 464– 465. Choleva L, Musilova Z, Kohoutova-Sediva A, Paces J, Rab P, Janko K. 2014. Distinguishing between incomplete lineage sorting and genomic introgressions: complete fixation of allospecific mitochondrial DNA in a sexually reproducing fish (Cobitis; Teleostei), despite clonal reproduction of hybrids. PLoS One  9: e80641. Google Scholar CrossRef Search ADS PubMed  Christin PA, Besnard G, Samaritani E, Duvall MR, Hodkinson TR, Savolainen V, Salamin N. 2008. Oligocene CO2 decline promoted C4 photosynthesis in grasses. Current Biology: CB  18: 37– 43. Google Scholar CrossRef Search ADS PubMed  Comes HP, Kadereit JW. 1998. The effect of Quaternary climatic changes on plant distribution and evolution. Trends in Plant Science  3: 432– 438. Google Scholar CrossRef Search ADS   Crawford DJ, Archibald JK. 2017. Island floras as model systems for studies of plant speciation: prospects and challenges. Journal of Systematics and Evolution  55: 1– 15. Google Scholar CrossRef Search ADS   Doyle JJ, Doyle JL. 1987. A rapid DNA isolation procedure for small quantities of fresh leaf material. Phytochemical Bulletin, Botanical Society of America  19: 11– 15. Drummond AJ, Rambaut A. 2007. BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evolutionary Biology  7: 214. Google Scholar CrossRef Search ADS PubMed  Edgar RC. 2004. MUSCLE: multiple sequence alignment with high accuracy and high throughput. Nucleic Acids Research  32: 1792– 1797. Google Scholar CrossRef Search ADS PubMed  ESRI. 2010. ESRI info: Company history . Available at: http://www.esri.com/about-esri/about/history.html (accessed 6 May 2010). Evanno G, Regnaut S, Goudet J. 2005. Detecting the number of clusters of individuals using the software structure: a simulation study. Molecular Ecology  14: 2611– 2620. Google Scholar CrossRef Search ADS PubMed  Excoffier L, Laval G, Schneider S. 2005. ARLEQUIN (version 3.0): an integrated software package for population genetics data analysis. Evolutionary Bioinformatics  1: 47– 50. Google Scholar CrossRef Search ADS   Excoffier L, Smouse PE, Quattro JM. 1992. Analysis of molecular variance inferred from metric distances among DNA haplotypes: application to human mitochondrial DNA restriction data. Genetics  131: 479– 491. Google Scholar PubMed  Fan DM, Yue JP, Nie ZL, Li ZM, Comes HP, Sun H. 2013. Phylogeography of Sophora davidii (Leguminosae) across the ‘Tanaka-Kaiyong Line’, an important phytogeographic boundary in Southwest China. Molecular Ecology  22: 4270– 4288. Google Scholar CrossRef Search ADS PubMed  Fang XM, Lü LQ, Yang SL, Li JJ, An ZS, Jiang PA, Chen XL. 2002. Loess in Kunlun Mountains and its implications on desert development and Tibetan Plateau uplift in west China. Science in China Series D  45: 289– 299. Google Scholar CrossRef Search ADS   Farris JS, Källersjö M, Kluge AG, Bult C. 1995. Testing the significance of incongruence. Cladistics  10: 315– 319. Google Scholar CrossRef Search ADS   Favre A, Päckert M, Pauls SU, Jähnig SC, Uhl D, Michalak I, Muellner-Riehl AN. 2015. The role of the uplift of the Qinghai-Tibetan Plateau for the evolution of Tibetan biotas. Biological Reviews of the Cambridge Philosophical Society  90: 236– 253. Google Scholar CrossRef Search ADS PubMed  Felsenstein J. 1985. Confidence limits on phylogenies: an approach using the bootstrap. Evolution  39: 783– 791. Google Scholar CrossRef Search ADS PubMed  Feng YX, Wang XQ, Pan KY, Hong DY. 1998. A revaluation of the systematic positions of the Cercidiphyllaceae and Daphniphyllaceae based on rbcL gene sequence analysis, with reference to the relationship in the “lower” Hamamelidae. Acta Phytotaxonomica Sinica  36: 411– 422. Freeland JR, Kirk H, Petersen SD. 2011. Molecular ecology, 2nd edn . Chichester: John Wiley & Sons. Google Scholar CrossRef Search ADS   Fuertes-Aguilar J, Nieto-Feliner G. 2003. Additive polymorphisms and reticulation in an ITS phylogeny of thrifts (Armeria, Plumbaginaceae). Molecular Phylogenetics and Evolution  28: 430– 447. Google Scholar CrossRef Search ADS PubMed  Fuertes-Aguilar J, Rosselló JA, Nieto-Feliner G. 1999. Nuclear ribosomal DNA (nrDNA) concerted evolution in natural and artificial hybrids of Armeria (Plumbaginaceae). Molecular Ecology  8: 1341– 1346. Google Scholar CrossRef Search ADS PubMed  Fu YX. 1997. Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics  147: 915– 925. Google Scholar PubMed  Funk DJ, Omland KE. 2003. Species level paraphyly and polyphyly: frequency, causes, and consequences, with insights from animal mitochondrial DNA. Annual Review of Ecology, Evolution, and Systematics  34: 397– 423. Google Scholar CrossRef Search ADS   Gao YD, Harris AJ, Zhou SD, He XJ. 2013. Evolutionary events in Lilium (including Nomocharis, Liliaceae) are temporally correlated with orogenies of the Q-T plateau and the Hengduan Mountains. Molecular Phylogenetics and Evolution  68: 443– 460. Google Scholar CrossRef Search ADS PubMed  Gao QB, Zhang DJ, Duan YZ, Zhang FQ, Li YH, Fu PC, Chen SL. 2012. Intraspecific divergences of Rhodiola alsia (Crassulaceae) based on plastid DNA and internal transcribed spacer fragments. Botanical Journal of the Linnean Society  168: 204– 215. Google Scholar CrossRef Search ADS   Garnhart NJ. 2001. Binthere V1.0, a program to bin AFLP data . Durham: University of New Hampshire. Guo ZT, Ruddiman WF, Hao QZ, Wu HB, Qiao YS, Zhu RX, Peng SZ, Wei JJ, Yuan BY, Liu TS. 2002. Onset of Asian desertification by 22 Myr ago inferred from loess deposits in China. Nature  416: 159– 163. Google Scholar CrossRef Search ADS PubMed  Hamilton MB. 2009. Population genetics . Chichester: John Wiley & Sons Ltd. Harrison TM, Copeland P, Kidd WS, Yin A. 1992. Raising Tibet. Science  255: 1663– 1670. Google Scholar CrossRef Search ADS PubMed  Harvey PH, May RM, Nee S. 1994. Phylogenies without fossils. Evolution  48: 523– 529. Google Scholar CrossRef Search ADS PubMed  Hauquier F, Leliaert F, Rigaux A, Derycke S, Vanreusel A. 2017. Distinct genetic differentiation and species diversification within two marine nematodes with different habitat preference in Antarctic sediments. BMC Evolutionary Biology  17: 120. Google Scholar CrossRef Search ADS PubMed  Hewitt GM. 1996. Some genetic consequences of ice ages, and their role in divergence and speciation. Biological Journal of the Linnean Society  58: 247– 276. Google Scholar CrossRef Search ADS   Hewitt GM. 1999. Post-glacial re-colonization of European Biota. Biological Journal of the Linnean Society  68: 87– 112. Google Scholar CrossRef Search ADS   Hewitt GM. 2000. The genetic legacy of the Quaternary ice ages. Nature  405: 907– 913. Google Scholar CrossRef Search ADS PubMed  Hewitt GM. 2004. Genetic consequences of climatic oscillations in the Quaternary. Philosophical Transactions of the Royal Society of London B: Biological Sciences  359: 183– 195. Google Scholar CrossRef Search ADS   Hewitt GM. 2011. Quaternary phylogeography: the roots of hybrid zones. Genetica  139: 617– 638. Google Scholar CrossRef Search ADS PubMed  Hudson RR. 1990. Gene genealogies and the coalescent process. In: Futuyma D, Antonovics J, eds. Oxford surveys in evolutionary biology . Oxford: Oxford University Press, 1– 44. Jia DR, Abbott RJ, Liu TL, Mao KS, Bartish IV, Liu JQ. 2012. Out of the Qinghai-Tibet Plateau: evidence for the origin and dispersal of Eurasian temperate plants from a phylogeographic study of Hippophaë rhamnoides (Elaeagnaceae). The New Phytologist  194: 1123– 1133. Google Scholar CrossRef Search ADS PubMed  Jia DR, Liu TL, Wang LY, Zhou DW, Liu JQ. 2011. Evolutionary history of an alpine shrub Hippophae tibetana (Elaeagnaceae): allopatric divergence and regional expansion. Botanical Journal of the Linnean Society  102: 37– 50. Google Scholar CrossRef Search ADS   Jiang DB, Lang XM, Tian ZP, Guo DL. 2011. Last glacial maximum climate over China from PMIP simulations. Palaeogeography, Palaeoclimatology, Palaeoecology  309: 347– 357. Google Scholar CrossRef Search ADS   Kuritzin A, Kischka T, Schmitz J, Churakov G. 2016. Incomplete lineage sorting and hybridization statistics for large-scale retroposon insertion data. PLoS Computational Biology  12: e1004812. Google Scholar CrossRef Search ADS PubMed  de Lafontaine G, Amasifuen Guerra CA, Ducousso A, Sánchez-Goñi MF, Petit RJ. 2014. Beyond skepticism: uncovering cryptic refugia using multiple lines of evidence. The New Phytologist  204: 450– 454. Google Scholar CrossRef Search ADS PubMed  Li JJ, Shi YF, Li BY. 1995. Uplift of the Qinghai-Xizang (Tibet) Plateau and global change . Lanzhou: Lanzhou University Press. Li ZH, Zou JB, Mao KS, Lin K, Li HP, Liu JQ, Källman T. 2012. Population genetic evidence for complex evolutionary histories of four high altitude juniper species in the Qinghai-Tibetan Plateau. Evolution  66: 831– 845. Google Scholar CrossRef Search ADS PubMed  Librado P, Rozas J. 2009. DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics  25: 1451– 1452. Google Scholar CrossRef Search ADS PubMed  Liu JQ. 2004. Uniformity of karyotypes in Ligularia (Asteraceae: Senecioneae), a highly-diversified genus of the eastern Qinghai-Tibet Plateau highlands and adjacent areas. Botanical Journal of the Linnean Society  144: 329– 342. Google Scholar CrossRef Search ADS   Liu B, Abbott RJ, Lu Z, Tian B, Liu J. 2014a. Diploid hybrid origin of Ostryopsis intermedia (Betulaceae) in the Qinghai-Tibet Plateau triggered by Quaternary climate change. Molecular Ecology  23: 3013– 3027. Google Scholar CrossRef Search ADS   Liu JQ, Duan YW, Hao G, Ge XJ, Sun H. 2014b. Evolutionary history and underlying adaptation of alpine plants on the Qinghai-Tibet Plateau. Journal of Systematics and Evolution  52: 241– 249. Google Scholar CrossRef Search ADS   Liu JQ, Gao TG, Chen ZD, Lu AM. 2002. Molecular phylogeny and biogeography of the Qinghai-Tibet Plateau endemic Nannoglottis (Asteraceae). Molecular Phylogenetics and Evolution  23: 307– 325. Google Scholar CrossRef Search ADS PubMed  Liu YP, Su X, He YH, Han LM, Huang YY, Wang ZZ. 2015. Evolutionary history of Orinus thoroldii (Poaceae), endemic to the western Qinghai-Tibetan Plateau in China. Biochemical Systematics and Ecology  59: 159– 167. Google Scholar CrossRef Search ADS   Liu JQ, Sun YS, Ge XJ, Gao LM, Qiu YX. 2012. Phylogeographic studies of plants in China: advances in the past and directions in the future. Journal of Systematics and Evolution  50: 267– 275. Google Scholar CrossRef Search ADS   Liu JQ, Wang YJ, Wang AL, Hideaki O, Abbott RJ. 2006. Radiation and diversification within the Ligularia-Cremanthodium-Parasenecio complex (Asteraceae) triggered by uplift of the Qinghai-Tibetan Plateau. Molecular Phylogenetics and Evolution  38: 31– 49. Google Scholar CrossRef Search ADS PubMed  Ma YZ, Li ZH, Wang X, Shang BL, Wu GL, Wang YJ. 2014. Phylogeography of the genus Dasiphora (Rosaceae) in the Qinghai-Tibetan Plateau: divergence blurred by expansion. Botanical Journal of the Linnean Society  111: 777– 788. Google Scholar CrossRef Search ADS   Mao K, Hao G, Liu J, Adams RP, Milne RI. 2010. Diversification and biogeography of Juniperus (Cupressaceae): variable diversification rates and multiple intercontinental dispersals. The New Phytologist  188: 254– 272. Google Scholar CrossRef Search ADS PubMed  Mayr E. 1942. Systematics and the origin of species, from the viewpoint of a zoologist . New York: Columbia University Press. Meng L, Chen G, Li Z, Yang Y, Wang Z, Wang L. 2015. Refugial isolation and range expansions drive the genetic structure of Oxyria sinensis (Polygonaceae) in the Himalaya-Hengduan Mountains. Scientific Reports  5: 10396. Google Scholar CrossRef Search ADS PubMed  Meng HH, Gao XY, Huang JF, Zhang ML. 2014. Plant phylogeography in arid northwest China: retrospectives and perspectives. Journal of Systematics and Evolution  52: 1– 16. Google Scholar CrossRef Search ADS   Meng L, Yang R, Abbott RJ, Miehe G, Hu T, Liu J. 2007. Mitochondrial and chloroplast phylogeography of Picea crassifolia Kom. (Pinaceae) in the Qinghai-Tibetan Plateau and adjacent highlands. Molecular Ecology  16: 4128– 4137. Google Scholar CrossRef Search ADS PubMed  Miraldo A, Hewitt GM, Paulo OS, Emerson BC. 2011. Phylogeography and demographic history of Lacerta lepida in the Iberian Peninsula: multiple refugia, range expansions and secondary contact zones. BMC Evolutionary Biology  11: 170. Google Scholar CrossRef Search ADS PubMed  Mittermeier RA, Gil PR, Hoffman M, Pilgrim J, Brooks T, Mittermeier CG, Lamoreux J, da Fonseca GAB. 2005. Hotspots revisited: earth’s biologically richest and most endangered terrestrial ecoregions . Chicago: University of Chicago Press. Nason JD, Hamrick JL, Fleming TH. 2002. Historical vicariance and postglacial colonization effects on the evolution of genetic structure in Lophocereus, a Sonoran Desert columnar cactus. Evolution  56: 2214– 2226. Google Scholar CrossRef Search ADS PubMed  Nei M, Li WH. 1979. Mathematical model for studying genetical variation in terms of restriction endonucleases. Proceedings of the National Academy of Sciences of the United States of America  76: 5269– 5273. Google Scholar CrossRef Search ADS PubMed  Nosil P. 2012. Ecological speciation . New York: Oxford University Press. Google Scholar CrossRef Search ADS   Opgenoorth L, Vendramin GG, Mao K, Miehe G, Miehe S, Liepelt S, Liu J, Ziegenhagen B. 2010. Tree endurance on the Tibetan Plateau marks the world’s highest known tree line of the Last Glacial Maximum. The New Phytologist  185: 332– 342. Google Scholar CrossRef Search ADS PubMed  Osborne CP, Beerling DJ. 2006. Nature’s green revolution: the remarkable evolutionary rise of C4 plants. Philosophical Transactions of the Royal Society of London B: Biological Sciences  361: 173– 194. Google Scholar CrossRef Search ADS   Parchman TL, Benkman CW, Britch SC. 2006. Patterns of genetic variation in the adaptive radiation of New World crossbills (Aves: Loxia). Molecular Ecology  15: 1873– 1887. Google Scholar CrossRef Search ADS PubMed  Peterson PM, Romaschenko K, Arrieta YH. 2016. A molecular phylogeny and classification of the Cynodonteae (Poaceae: Chloridoideae) with four new genera: Orthacanthus, Triplasiella, Tripogonella, and Zaqiqah; three new subtribes: Dactylocteniinae, Orininae, and Zaqiqahinae; and a subgeneric classification of Distichlis. Taxon  65: 1263– 1287. Google Scholar CrossRef Search ADS   Peterson PM, Romaschenko K, Snow N, Johnson G. 2012. A molecular phylogeny and classification of Leptochloa (Poaceae: Chloridoideae: Chlorideae) sensu lato and related genera. Annals of Botany  109: 1317– 1330. Google Scholar CrossRef Search ADS PubMed  Petit R, Aguinagalde I, de Beaulieu JL, Bittkau C, Brewer S, Cheddadi R, Ennos R, Fineschi S, Grivet D, Lascoux M, Mohanty A, Müller-Starck G, Demesure-Musch B, Palmé A, Martín JP, Rendell S, Vendramin GG. 2003. Glacial refugia: hotspots but not melting pots of genetic diversity. Science  300: 1563– 1565. Google Scholar CrossRef Search ADS PubMed  Petit RJ, Excoffier L. 2009. Gene flow and species delimitation. Trends in Ecology & Evolution  24: 386– 393. Google Scholar CrossRef Search ADS PubMed  Pons O, Chaouche K. 1995. Estimation, variance and optimal sampling of gene diversity II. Diploid locus. Theoretical and Applied Genetics  91: 122– 130. Google Scholar CrossRef Search ADS PubMed  Pons O, Petit RJ. 1996. Measuring and testing genetic differentiation with ordered versus unordered alleles. Genetics  144: 1237– 1245. Google Scholar PubMed  Posada D. 1998. jModelTest: phylogenetic model averaging. Molecular Biology and Evolution  25: 1235– 1256. Pritchard JK, Stephens M, Donnelly P. 2000. Inference of population structure using multilocus genotype data. Genetics  155: 945– 959. Google Scholar PubMed  Qiu YX, Fu CX, Comes HP. 2011. Plant molecular phylogeography in China and adjacent regions: tracing the genetic imprints of Quaternary climate and environmental change in the world’s most diverse temperate flora. Molecular Phylogenetics and Evolution  59: 225– 244. Google Scholar CrossRef Search ADS PubMed  Qiu YX, Liu YF, Kang M, Yi GM, Huang HW. 2013. Spatial and temporal population genetic variation and structure of Nothotsuga longibracteata (Pinaceae), a relic conifer species endemic to subtropical China. Genetics and Molecular Biology  36: 598– 607. Google Scholar CrossRef Search ADS PubMed  Rambaut A. 2009. FIGTREE 1.3.1. Available at: http://tree.bio.ed.ac.uk/software/figtree (accessed 7 March 2014). Rambaut A, Drummond AJ. 2007. Tracer. Available at: http://beast.bio.ed.ac.uk/tracer/ (accessed 11 December 2013). Rambaut A, Drummond AJ. 2009. TRACER v1.5. Available at: http://tree.bio.ed.ac.uk/software/tracer/ (accessed 7 March 2014). Rieseberg LH, Church SA, Morjan CL. 2004. Integration of populations and differentiation of species. The New Phytologist  161: 59– 69. Google Scholar CrossRef Search ADS PubMed  Rogers AR, Harpending H. 1992. Population growth makes waves in the distribution of pairwise genetic differences. Molecular Biology and Evolution  9: 552– 569. Google Scholar PubMed  Rohlf FJ. 2000. NTSYS-pc version 2.1: numerical taxonomy and multivariance analysis system . Setauket: Applied Biostatistics Inc., Exeter Software. Ronquist F, Huelsenbeck JP. 2003. MrBayes 3: Bayesian phylogenetic inference under mixed models. Bioinformatics  19: 1572– 1574. Google Scholar CrossRef Search ADS PubMed  Ruddiman W. 1998. Geology: early uplift in Tibet? Nature  394: 723– 725. Google Scholar CrossRef Search ADS   Sakaguchi S, Takeuchi Y, Yamasaki M, Sakurai S, Isagi Y. 2011. Lineage admixture during postglacial range expansion is responsible for the increased gene diversity of Kalopanax septemlobus in a recently colonised territory. Heredity  107: 338– 348. Google Scholar CrossRef Search ADS PubMed  Sang T, Crawford D, Stuessy T. 1997. Chloroplast DNA phylogeny, reticulate evolution, and biogeography of Paeonia (Paeoniaceae). American Journal of Botany  84: 1120. Google Scholar CrossRef Search ADS PubMed  Sardell JM, Uy JAC. 2016. Hybridization following recent secondary contact results in asymmetric genotypic and phenotypic introgression between island species of Myzomela honeyeaters. Evolution  70: 257– 269. Google Scholar CrossRef Search ADS PubMed  Schluter D. 2009. Evidence for ecological speciation and its alternative. Science  323: 737– 741. Google Scholar CrossRef Search ADS PubMed  Segovia RA, Pérez MF, Hinojosa LF. 2012. Genetic evidence for glacial refugia of the temperate tree Eucryphia cordifolia (Cunoniaceae) in southern South America. American Journal of Botany  99: 121– 129. Google Scholar CrossRef Search ADS PubMed  Shackleton NJ, Opdyke ND. 1973. Oxygen isotope and palaeomagnetic stratigraphy of Equatorial Pacific core V28-238: oxygen isotope temperatures and ice volumes on a 105 year and 106 year scale. Quaternary Research  3: 39– 55. Google Scholar CrossRef Search ADS   Shaw J, Shafer HL, Leonard OR, Margaret JK, Schorr M, Morris AB. 2017. Chloroplast DNA sequence utility for the lowest phylogenetic and phylogeographic inferences in angiosperms: the tortoise and the hare IV. American Journal of Botany  101: 1987– 2004. Google Scholar CrossRef Search ADS   Shi YF, Li JJ, Li BY. 1998. Uplift and environmental changes of Qinghai-Tibetan Plateau in the late Cenozoic . Guangzhou: Guangdong Science and Technology Press. Slatkin M, Hudson RR. 1991. Pairwise comparisons of mitochondrial DNA sequences in stable and exponentially growing populations. Genetics  129: 555– 562. Google Scholar PubMed  Soreng RJ, Peterson PM, Romaschenko K, Davidse G, Teisher JK, Clark LG, Barbéra P, Gillespie LJ, Zuloaga FO. 2017. A worldwide phylogenetic classification of the Poaceae (Gramineae) II: an update and comparison of two 2015 classifications. Journal of Systematics and Evolution  55: 259– 290. Google Scholar CrossRef Search ADS   Sousa V, Hey J. 2013. Understanding the origin of species with genome-scale data: modelling gene flow. Nature Reviews Genetics  14: 404– 414. Google Scholar CrossRef Search ADS PubMed  Stephens M, Smith NJ, Donnelly P. 2011. A new statistical method for haplotype reconstruction from population data. American Journal of Human Genetics  68: 978– 989. Google Scholar CrossRef Search ADS   Stewart JR, Lister AM, Barnes I, Dalen L. 2010. Refugia revisited: individualistic responses of species in space and time. Proceedings of the Royal Society of London B: Biological Sciences  277: 661– 671. Google Scholar CrossRef Search ADS   Su X, Cai LB. 2009. Orinus longiglumis (Poaceae: Chloridoideae), a new species from Xizang (Tibet), China. Annales Botanici Fennici  46: 143– 147. Google Scholar CrossRef Search ADS   Su X, Liu YP, Wu GL, Luo WC, Liu JQ. 2017. A taxonomic revision of Orinus (Poaceae) with a new species, O. intermedius, from the Qinghai-Tibet Plateau. Novon  25: 206– 213. Google Scholar CrossRef Search ADS   Su X, Wu G, Li L, Liu J. 2015. Species delimitation in plants using the Qinghai-Tibet Plateau endemic Orinus (Poaceae: Tridentinae) as an example. Annals of Botany  116: 35– 48. Google Scholar CrossRef Search ADS PubMed  Sun Y, Abbott RJ, Li L, Li L, Zou J, Liu J. 2014. Evolutionary history of Purple cone spruce (Picea purpurea) in the Qinghai-Tibet Plateau: homoploid hybrid origin and Pleistocene expansion. Molecular Ecology  23: 343– 359. Google Scholar CrossRef Search ADS PubMed  Sun YB, An ZS. 2001. History and variability of Asian interior aridity recorded by eolian flux in the Chinese Loess Plateau during the past 7 Ma. Science in China Series D: Earth Sciences  45: 420– 429. Google Scholar CrossRef Search ADS   Sun H, McLewin W, Fay MF. 2001. Molecular phylogeny of Helleborus (Ranunculaceae), with an emphasis on the East Asia-Mediterranean disjunction. Taxon  50: 1001– 1018. Google Scholar CrossRef Search ADS   Swofford DL. 2002. PAUP*: Phylogenetic analysis using parsimony (*and other methods), version 4 . Sunderland: Sinauer Associates. Tajima F. 1989. Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics  123: 585– 595. Google Scholar PubMed  Tamura K, Peterson D, Peterson N, Stecher G, Nei M, Kumar S. 2011. MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Molecular Biology and Evolution  28: 2731– 2739. Google Scholar CrossRef Search ADS PubMed  Tang LY, Shen CM. 1996. Late Cenozoic vegetational history and climatic characteristics of Qinghai-Xizang Plateau. Acta Micropalaeontologica Sinica  13: 321– 337. Tian X, Luo J, Wang A, Mao K, Liu J. 2011. On the origin of the woody buckwheat Fagopyrum tibeticum (=Parapteropyrum tibeticum) in the Qinghai-Tibetan Plateau. Molecular Phylogenetics and Evolution  61: 515– 520. Google Scholar CrossRef Search ADS PubMed  Tzedakis PC, Emerson BC, Hewitt GM. 2013. Cryptic or mystic? Glacial tree refugia in northern Europe. Trends in Ecology & Evolution  28: 696– 704. Google Scholar CrossRef Search ADS PubMed  Tzedakis PC, Lawson IT, Forgley MR, Hewitt GM. 2002. Buffered tree population changes in Quaternary refugium: evolutionary implications. Science  292: 267– 269. Vicentini A, Barber JC, Aliscioni SS, Giussani LM, Kellogg EA. 2008. The age of the grasses and clusters of origins of C4 photosynthesis. Global Change Biology  14: 2963– 2977. Google Scholar CrossRef Search ADS   de Villiers MJ, Pirie MD, Hughes M, Möller M, Edwards TJ, Bellstedt DU. 2013. An approach to identify putative hybrids in the ‘coalescent stochasticity zone’, as exemplified in the African plant genus Streptocarpus (Gesneriaceae). The New Phytologist  198: 284– 300. Google Scholar CrossRef Search ADS PubMed  Vos P, Hogers R, Bleeker M, Reijans M, van de Lee T, Hornes M, Frijters A, Pot J, Peleman J, Kuiper M. 1995. AFLP: a new technique for DNA fingerprinting. Nucleic Acids Research  23: 4407– 4414. Google Scholar CrossRef Search ADS PubMed  Voss N, Eckstein RL, Durka W. 2012. Range expansion of a selfing polyploid plant despite widespread genetic uniformity. Annals of Botany  110: 585– 593. Google Scholar CrossRef Search ADS PubMed  Wan DS, Feng JJ, Jiang DC, Mao KS, Duan YW, Miehe G, Opgenoorth L. 2016. The Quaternary evolutionary history, potential distribution dynamics, and conservation implications for a Qinghai-Tibet Plateau endemic herbaceous perennial, Anisodus tanguticus (Solanaceae). Ecology and Evolution  6: 1977– 1995. Google Scholar CrossRef Search ADS PubMed  Wang Q, Abbott RJ, Yu QS, Lin K, Liu JQ. 2013. Pleistocene climate change and the origin of two desert plant species, Pugionium cornutum and Pugionium dolabratum (Brassicaceae), in northwest China. The New Phytologist  199: 277– 287. Google Scholar CrossRef Search ADS PubMed  Wang L, Abbott RJ, Zheng W, Chen P, Wang Y, Liu J. 2009a. History and evolution of alpine plants endemic to the Qinghai-Tibetan Plateau: Aconitum gymnandrum (Ranunculaceae). Molecular Ecology  18: 709– 721. Google Scholar CrossRef Search ADS   Wang GN, He XY, Miehe G, Mao KS. 2014. Phylogeography of the Qinghai-Tibet Plateau endemic alpine herb Pomatosace filicula (Primulaceae). Journal of Systematics and Evolution  52: 1– 14. Google Scholar CrossRef Search ADS   Wang LY, Ikeda H, Liu TL, Wang YJ, Liu JQ. 2009b. Repeated range expansion and glacial endurance of Potentilla glabra (Rosaceae) in the Qinghai-Tibetan Plateau. Journal of Integrative Plant Biology  51: 698– 706. Google Scholar CrossRef Search ADS   Wang YJ, Susanna A, Von Raab-Straube E, Milne R, Liu JQ. 2009c. Island-like radiation of Saussurea (Asteraceae: Cardueae) triggered by uplifts of the Qinghai-Tibetan Plateau. Botanical Journal of the Linnean Society  97: 893– 903. Google Scholar CrossRef Search ADS   Weir BS. 1996. Genetic data analysis II . Sunderland: Sinauer Associates. Wen J, Nie ZL, Ickert-Bond SM. 2016. Intercontinental disjunctions between eastern Asia and western North America in vascular plants highlight the biogeographic importance of the Bering land bridge from late Cretaceous to Neogene. Journal of Systematics and Evolution  54: 469– 490. Google Scholar CrossRef Search ADS   Wen J, Zhang JQ, Nie ZL, Zhong Y, Sun H. 2014. Evolutionary diversifications of plants on the Qinghai-Tibetan Plateau. Frontiers in Genetics  5: 4. Google Scholar PubMed  White TJ, Bruns T, Lee S, Taylor J. 1990. Amplification and direct sequencing of fungal ribosomal RNA genes for phylogenetics. In: Innis M, Gelfand D, Sninsky J, White T, eds. PCR protocols . San Diego: Academic Press, 315– 322. Google Scholar CrossRef Search ADS   Wolfe KH, Li WH, Sharp PM. 1987. Rates of nucleotide substitution vary greatly among plant mitochondrial, chloroplast, and nuclear DNAs. Proceedings of the National Academy of Sciences of the United States of America  84: 9054– 9058. Google Scholar CrossRef Search ADS PubMed  Wu CY. 1988. Hengduan mountains flora and their significance. Journal of Japanese Botany  63: 297– 311. Wu ZH, Barosh PJ, Wu ZH, Hu DG, Zhao X, Ye PS. 2008. Vast early Miocene lakes of the central Tibetan Plateau. Geological Society of America Bulletin  120: 1326– 1337. Google Scholar CrossRef Search ADS   Wu LL, Cui XK, Milne RI, Sun YS, Liu JQ. 2010. Multiple autopolyploidizations and range expansion of Allium przewalskianum Regel. (Alliaceae) in the Qinghai-Tibetan Plateau. Molecular Ecology  19: 1691– 1704. Google Scholar CrossRef Search ADS PubMed  Xing YW, Ree RH. 2017. Uplift-driven diversification in the Hengduan Mountains, a temperate biodiversity hotspot. Proceedings of the National Academy of Sciences of the United States of America  114: E3444– E3451. Google Scholar CrossRef Search ADS PubMed  Xu T, Abbott RJ, Milne RI, Mao K, Du FK, Wu G, Ciren Z, Miehe G, Liu J. 2010. Phylogeography and allopatric divergence of cypress species (Cupressus L.) in the Qinghai-Tibetan Plateau and adjacent regions. BMC Evolutionary Biology  10: 194. Google Scholar CrossRef Search ADS PubMed  Yang FS, Li YF, Ding X, Wang XQ. 2008. Extensive population expansion of Pedicularis longiflora (Orobanchaceae) on the Qinghai-Tibetan Plateau and its correlation with the Quaternary climate change. Molecular Ecology  17: 5135– 5145. Google Scholar CrossRef Search ADS PubMed  Yang FS, Qin AL, Li YF, Wang XQ. 2012. Great genetic differentiation among populations of Meconopsis integrifolia and its implication for plant speciation in the Qinghai-Tibetan Plateau. PLoS One  7: e37196. Google Scholar CrossRef Search ADS PubMed  Yap IV, Nelson RJ. 1996. Winboot: a program for performing bootstrap analysis of binary data to determine the confidence limits of UPGMA-based dendrograms . IRRI Discussion Paper Series No. 14. Manila: International Rice Research Institute. Yasuda N, Taquet C, Nagai S, Fortes M, Fan TY, Harii S, Yoshida T, Sito Y, Nadaoka K. 2015. Genetic diversity, paraphyly and incomplete lineage sorting of mtDNA, ITS2 and microsatellite flanking region in closely related Heliopora species (Octocorallia). Molecular Phylogenetics and Evolution  93: 161– 171. Google Scholar CrossRef Search ADS PubMed  Yu QS, Wang Q, Wu GL, Ma YZ, He XY, Wang X, Xie PH, Hu LH, Liu JQ. 2013. Genetic differentiation and delimitation of Pugionium dolabratum and Pugionium cornutum (Brassicaceae). Plant Systematics and Evolution  299: 1355– 1365. Google Scholar CrossRef Search ADS   Zhang Q, Chiang TY, George M, Liu JQ, Abbott RJ. 2005. Phylogeography of the Qinghai-Tibetan Plateau endemic Juniperus przewalskii (Cupressaceae) inferred from chloroplast DNA sequence variation. Molecular Ecology  14: 3513– 3524. Google Scholar CrossRef Search ADS PubMed  Zhang JQ, Meng SY, Wen J, Rao GY. 2015. DNA barcoding of Rhodiola (Crassulaceae): a case study on a group of recently diversified medicinal plants from the Qinghai-Tibetan Plateau. PLoS One  10: 1– 15. Zhao JL, Gugger PF, Xia YM, Li QJ. 2016. Ecological divergence of two closely related Roscoea species associated with late Quaternary climate change. Journal of Biogeography  43: 1990– 2001. Google Scholar CrossRef Search ADS   Zheng D. 1996. The system of physico-geographical regions of the Qinghai-Xizang (Tibet) Plateau. Science China: Earth Sciences  39: 410– 417. Zheng HB, Powell CM, Rea DK, Wang JL, Wang PX. 2004. Late Miocene and mid-Pliocene enhancement of the East Asian monsoon as viewed from the land and sea. Global and Planetary Change  41: 147– 155. Google Scholar CrossRef Search ADS   Zheng BX, Xu QQ, Shen YP. 2002. The relationship between climate change and Quaternary glacial cycles on the Qinghai-Tibetan Plateau: review and speculation. Quaternary International  97–98: 93– 101. Google Scholar CrossRef Search ADS   Zhou Y, Duvaux L, Ren G, Zhang L, Savolainen O, Liu J. 2017. Importance of incomplete lineage sorting and introgression in the origin of shared genetic variation between two closely related pines with overlapping distributions. Heredity  118: 211– 220. Google Scholar CrossRef Search ADS PubMed  Zhou SL, Qian P. 2003. Matrix Generator (MG): a program for creating 0/1 matrix from DNA fragments. Acta Botanica Sinica  45: 766. Zimmer EA, Wen J. 2012. Using nuclear gene data for plant phylogenetics: progress and prospects. Molecular Phylogenetics and Evolution  65: 774– 785. Google Scholar CrossRef Search ADS PubMed  Zimmer EA, Wen J. 2015. Using nuclear gene data for plant phylogenetics: progress and prospects II. Next-gen approaches. Journal of Systematics and Evolution  53: 371– 379. Google Scholar CrossRef Search ADS   Zou JB, Peng XL, Li L, Liu JQ, Miehe G, Opgenoorth L. 2012. Molecular phylogeography and evolutionary history of Picea likiangensis in the Qinghai-Tibetan Plateau inferred from mitochondrial and chloroplast DNA sequence variation. Journal of Systematics and Evolution  50: 341– 350. Google Scholar CrossRef Search ADS   Zuo YJ, Wen J, Ma JS, Zhou SL. 2015. Evolutionary radiation of the Panax bipinnatifidus species complex (Araliaceae) in the Sino-Himalayan region of eastern Asia as inferred from AFLP analysis. Journal of Systematics and Evolution  53: 210– 220. Google Scholar CrossRef Search ADS   © 2018 The Linnean Society of London, Botanical Journal of the Linnean Society http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Botanical Journal of the Linnean Society Oxford University Press

Phylogeography of Orinus (Poaceae), a dominant grass genus on the Qinghai-Tibet Plateau

Loading next page...
 
/lp/ou_press/phylogeography-of-orinus-poaceae-a-dominant-grass-genus-on-the-qinghai-5iruUR8TYZ
Publisher
The Linnean Society of London
Copyright
© 2018 The Linnean Society of London, Botanical Journal of the Linnean Society
ISSN
0024-4074
eISSN
1095-8339
D.O.I.
10.1093/botlinnean/box091
Publisher site
See Article on Publisher Site

Abstract

Abstract To better understand the responses of arid-adapted, alpine plants to Quaternary climatic oscillations, we investigated the genetic variation and phylogeographic history of Orinus, an endemic genus of Poaceae comprising three species from the dry grasslands of the Qinghai-Tibet Plateau (QTP) in China. We measured the genetic variation of 476 individuals from 88 populations using three maternally inherited plastid DNA markers (matK, rbcL and psbA-trnH), the biparentally inherited nuclear ribosomal internal transcribed spacer (nrITS) and amplified fragment length polymorphisms (AFLPs). We found that the plastid DNA, nrITS and AFLPs show considerable, recent differentiation among the species. We detected 14 plastid haplotypes (H1–H14), of which only three were shared among all species, and 30 nrITS ribotypes (S1–S30), of which one (S10) was shared between two species, O. kokonoricus and O. intermedius, but absent in O. thoroldii. The nrITS types formed clades that were inconsistent with species boundaries. Based on these data, we propose and illustrate a complex hypothesis for the evolutionary history of Orinus involving lineage sorting and introgression, the latter of which may explain the shared S10 nrITS type. The AFLP results showed clades corresponding to current species delineation and suggest that lineage sorting in the genus is probably complete. We estimated the crown age of Orinus to be 2.85 (95% highest posterior density: 0.58–12.45) Mya (late Pliocene), and subsequent divergence occurred in the Quaternary. Early divergences were allopatric. More recently, Orinus probably underwent regional expansions corresponding to Quaternary climatic changes, especially glaciation, which is consistent with our divergence time estimates. These climatic changes could have facilitated the S10 event and other hybridization events. Our data also suggest that species of this small genus of grasses survived the Quaternary glacial period in the extremely adverse habitats of the QTP. INTRODUCTION The Quaternary period began c. 2.58 Mya and has been characterized by climatic oscillations, especially glacial and interglacial cycles in the Northern Hemisphere (Shackleton & Opdyke, 1973). The glacial cycles co-varied with and probably profoundly affected other aspects of the climate, including the intensity of the Asian monsoon, even in unglaciated regions (An et al., 2001; Jiang et al., 2011). The Quaternary climatic changes caused changes in geographical distributions and population demography of plant species and consequently left lasting, detectable genetic imprints within and among species (Abbott et al., 2000; Avise, 2000; Hewitt, 2004, 2011; Qiu, Fu & Comes, 2011; Qiu et al., 2013; Wen et al., 2014, 2016). In Europe and North America, the fossil records of plant species and phylogeographic analyses indicate that a common pattern of geographical range shift was to retreat southward and to lower elevations during glacial periods and then to rapidly recolonize the northern areas and higher elevations during the interglacial and postglacial periods (Nason, Hamrick & Fleming, 2002; Stewart et al., 2010; Sakaguchi et al., 2011; Li et al., 2012; Segovia, Pérez & Hinojosa, 2012; Voss, Eckstein & Durka, 2012; Tzedakis, Emerson & Hewitt, 2013; de Lafontaine et al., 2014). In comparison with Europe and North America, the genetic structure and variation in Quaternary phylogeographic histories of plant species are much more pronounced in topographically complex regions of China, such as the Qinghai-Tibet Plateau (QTP) (e.g. Zhang et al., 2005; Meng et al., 2007; Wang et al., 2009a, 2014; Opgenoorth et al., 2010; Xu et al., 2010; Qiu et al., 2011; Zou et al., 2012; Wen et al., 2014; Liu et al., 2015; Zhang et al., 2015; Wan et al., 2016). The QTP and adjacent mountain ranges, such as the Hengduan and the Hindu-Kush Mountains, (hereafter QTP) were formed by the collision of the Indian subcontinent with Eurasia and contain the world’s tallest mountain (Mount Everest) and have an average elevation of > 4000 m (Zheng, 1996). The QTP is regarded as a hotspot for global biodiversity and supports many endemic genera and species (Wu, 1988; Mittermeier et al., 2005). The exceptional biodiversity of the QTP has been attributed to Quaternary climatic oscillations and Miocene-Pliocene or Miocene-Quaternary phases of uplifts of the QTP (Liu et al., 2002, 2006; Liu, 2004; Wang et al., 2009c; Mao et al., 2010; Xu et al., 2010; Jia et al., 2012; Wen et al., 2014; Liu et al., 2015). Climatic oscillations and uplifts may lead to allopatric speciation, detectable as genetic imprints that date to the Miocene, Pliocene or Quaternary. Recent studies have shown evidence for the importance of Quaternary climatic oscillations in explaining genetic variability, distributions and demographic patterns observed in modern populations of QTP species (Tang & Shen, 1996; Bawa et al., 2010; Jia et al., 2011; Liu et al., 2012, 2014b, 2015; Wan et al., 2016). In particular, many alpine species retreated to eastern or south-eastern refugia during glacial periods and recolonized the central plateau platform during the interglacial or postglacial stages (e.g. Zhang et al., 2005; Meng et al., 2007; Chen et al., 2008; Yang et al., 2008; Wu et al., 2010; Tian et al., 2011; Jia et al., 2012; Liu et al., 2015; Wan et al., 2016). In contrast, some other cold-tolerant alpine species were able to survive in situ through glacial periods in multiple ice-free refugia or microrefugia (Wang et al., 2009b; Opgenoorth et al., 2010; Jia et al., 2011; Ma et al., 2014). These two major types of evolutionary hypotheses about phylogeographic relationships between plant species and Quaternary climatic oscillations have been developed and tested by previous studies on the QTP (Wu et al., 2010; Jia et al., 2011; Tian et al., 2011; Ma et al., 2014; Liu et al., 2015). Results show that there was a close relationship between geographical distribution patterns of plant species from the QTP and Quaternary climatic oscillations, and independent genetic lineages appear in different plateau regions due to isolations among developed ice sheets. Given the probable impacts of ongoing climatic change on biodiversity, phylogeographic investigations of plant species in diverse regions, such as the QTP, are increasingly important. Although many phylogeographic studies of QTP plants have been performed, these have mostly focused on woody plants or herbaceous dicot lineages that occur in forested areas (e.g. Zhang et al., 2005; Meng et al., 2007; Chen et al., 2008; Tian et al., 2011; Jia et al., 2012; Wan et al., 2016). The grasslands comprise a significant portion of the QTP and are dominated by grasses (Poaceae). Therefore, disentangling the phylogeographic history of grasses is crucial to understanding the evolution of the QTP flora. In addition, these kind of grass-dominated landscapes are interesting since many of the vegetation types thus far studied are dominated by trees or shrubs. Grasses are cold-, warm- and arid-tolerant, making their study necessary to describe evolution in the QTP. Therefore, phylogeographic studies of grasses are needed to understand the past history of biogeographic and demographic change linked to past climatic oscillations and to make informed inferences about the future of the QTP grasslands (e.g. Qiu et al., 2011; Liu et al., 2012, 2014b, 2015; Wen et al., 2014). Here, we studied the phylogeographic history of the perennial grass genus, Orinus Hitchc., that is endemic to alpine, dry grasslands of the QTP and occurs in unusually sandy habitats (Chen & Phillips, 2006; Su & Cai, 2009; Su et al., 2015). Orinus comprises three morphologically and ecologically distinct species, O. thoroldii (Stapf ex Hemsl.) Bor, O. kokonoricus (K.S.Hao) Keng and O. intermedius X. Su & J. Quan Liu (Su et al., 2015, 2017), that are distributed at elevations between 2200 and 4800 m. Recently, Orinus was placed in subtribe Orininae P.M.Peterson, Romasch. & Y.Herrera with Cleistogenes Keng (Peterson, Romaschenko & Arrieta, 2016; Soreng et al., 2017). Furthermore, Orinus reproduces mainly through clonal reproduction via root suckers. Orinus spp. are dominant in places where they occur and are an exceptional study system for understanding the effects of the Quaternary climatic oscillations on the demography of QTP grassland species. Specifically, we investigated the biogeographic and phylogeographic history in Orinus using robust sampling of the three species throughout their geographical ranges using molecular markers, including plastid DNA, nuclear ribosomal internal transcribed spacer (nrITS) and amplified fragment length polymorphism (AFLP) data sets. We address the following questions: (1) What is the phylogeographic pattern in Orinus? (2) Are the phylogeographic patterns in Orinus consistent with southern glacial refugia and with other demographic patterns observed among diverse alpine plants on the QTP? MATERIAL AND METHODS Population sampling We sampled 476 individuals belonging to 88 populations from the distributional area of the genus including the type localities of all three species (Fig. 1A; see also Supporting Information, Table S1). From each population, we collected fresh leaf blades from three to 12 individuals and preserved them in silica gel. In each population, we selected individuals that were at least 20 m away from each other to avoid clones, and we recorded the latitude, longitude and elevation of each population using an Etrex GIS (Garmin, Taiwan, China). We also prepared several voucher specimens from each population and deposited them at the Herbarium of the Northwest Plateau Institute of Biology (HNWP), the Chinese Academy of Science, Xining, Qinghai Province, China. Figure 1. View largeDownload slide Sampling sites and geographical distribution of plastid DNA haplotypes observed in three Orinus spp. A, sampling region (dashed-line box). B, documented distribution range (dotted line), sampling sites and plastid DNA haplotype frequencies per population. Each pie graph represents a plastid DNA haplotype, and size of pie is proportional to the number of individuals bearing that haplotype. C, network for 14 plastid DNA haplotypes detected from three Orinus spp. Dinebra viscida was used as the outgroup. Sizes of network circles are proportional to haplotype frequencies. Numbers on the link lines between haplotypes represent the mutational sites. Colours for haplotypes are described in (B). Figure 1. View largeDownload slide Sampling sites and geographical distribution of plastid DNA haplotypes observed in three Orinus spp. A, sampling region (dashed-line box). B, documented distribution range (dotted line), sampling sites and plastid DNA haplotype frequencies per population. Each pie graph represents a plastid DNA haplotype, and size of pie is proportional to the number of individuals bearing that haplotype. C, network for 14 plastid DNA haplotypes detected from three Orinus spp. Dinebra viscida was used as the outgroup. Sizes of network circles are proportional to haplotype frequencies. Numbers on the link lines between haplotypes represent the mutational sites. Colours for haplotypes are described in (B). DNA extraction, amplification, sequencing and AFLP fingerprinting We extracted total DNA from silica gel-dried leaves using a modified 2× CTAB procedure (Doyle & Doyle, 1987). We used universal primers to amplify nuclear ITS and three plastid DNA regions: matK, rbcL and psbA-trnH (White et al., 1990; Sang, Crawford & Stuessy, 1997; Feng et al., 1998; Sun, McLewin & Fay, 2001). We performed the amplification using PCR in 25 µL reaction volumes, containing 1 µL plant template DNA, 0.2 mM MgCl2, 0.25 mM dNTPs, 2 µM each primer, 2.50 µL 10× PCR buffer and 0.25 µL Taq DNA polymerase (5 U/µL; TakaRa, Dalian, China). Our PCR thermocycling protocol comprised an initial 5-min denaturation period at 94 °C, followed by 36 cycles of denaturation at 94 °C for 50 s, annealing at 49–58 °C for 50 s and extension at 72 °C for 1 min, with a final extension at 72 °C for 10 min. We purified the PCR products with a TIAN-quick Midi Purification Kit (Tiangen, Beijing, China) in accordance with the manufacturer’s instructions, and we sequenced the purified products using an ABI 3130XL DNA Analyzer (Applied Biosystems, Foster City, CA, USA). We aligned the sequences of the plastid DNA and nrITS separately in MUSCLE (Edgar, 2004) and refined the alignments manually in MEGA 5 (Tamura et al., 2011). We combined the three plastid DNA regions because we did not detect phylogenetic incongruence among them using an incongruence length difference test under the parsimony criterion with 1000 replicates, Tree Bisection and Reconnection branch swapping and the number of trees retained for each replicate limited to 1000 (Farris et al., 1995). All newly obtained sequences of Orinus have been submitted to GenBank under the accession numbers KP302391–KP303274. We used a slightly modified protocol (Vos et al., 1995) to examine the AFLP of all the sampled populations (Supporting Information, Table S1) and used fluorescent primers for selective amplification (Zuo et al., 2015). DNA was digested with restriction enzymes PstI and MseI (40 U/µL, Beijing Dingguo Biotechnology Co., Ltd, China) for 5 h at 37 °C in a 20 µL volume, followed by ligation (40 units T4 DNA ligase, 4 h at 8 °C) to double-stranded PstI and MseI adapters. Two rounds of PCR amplification followed. Of 15 tested PstI/MseI primer combinations, eight were selected to produce the most repeatable and unambiguously scorable profiles. Selective fluorescence-labelled amplification products were separated and analysed on an ABI PRISM 377 DNA Sequencer (Applied Biosystems, Foster City, CA, USA) using GeneScan ROX-500, with an internal size standard. We used GeneScan 3.1 (Applied Biosystems, Foster City, CA, USA) to score the bands of AFLP amplification products. We imported the raw data into Binthere (Garnhart, 2001) and MG (Zhou & Qian, 2003) to generate a 0/1 matrix for data analyses. Population descriptive statistics We calculated the haplotype diversity (HE) of the plastid DNA for each population using DnaSP 5.0 (Librado & Rozas, 2009), and we used Permut (Pons & Petit, 1996) to find the average genetic diversity within populations (HS), total genetic diversity (HT) and the coefficients of differentiation (GST, NST) for each species, separately and combined. For the plastid DNA, we used permutation of GST and NST with 1000 replicates to infer significant phylogeographic structure (Pons & Chaouche, 1995). To examine population genetic structure, we used an analysis of molecular variance (AMOVA) with 1000 permutations (Excoffier, Smouse & Quattro, 1992) in Arlequin 3.11 (Excoffier, Laval & Schneider, 2005). We calculated plastid gene flow among populations (Nem) based on the equation Nem = ([1/FST] − 1)/2, where Ne was the effective population size of female individuals and m was the migration rate (Freeland, Kirk & Petersen, 2011). In nrITS, we assumed that a site was heterozygous if it had a double peak at the same position in both strands, and the weaker signal between base calls was at least 25% of the strength of the stronger (Fuertes-Aguilar, Rosselló & Nieto-Feliner, 1999; Fuertes-Aguilar & Nieto-Feliner, 2003). We also determined nrITS ribotypes with PHASE 1.0 for heterozygous individuals (Stephens, Smith & Donnelly, 2011). We used the Bayesian implementation of the k-mean method in STRUCTURE 2.2 (Pritchard, Stephens & Donnelly, 2000) to estimate the most likely number of populations (K) in Orinus according to the plastid DNA and nrITS data sets. Specifically, we performed analyses for K values ranging from 2 to 7, with ten replicates for each value of K and a burn-in of 2 × 106. We applied the ‘no admixture model’ and independent allele frequencies for these analyses. The most likely number of clusters was estimated according to the K resulting in the best improvement in the global likelihood (Evanno, Regnaut & Goudet, 2005). Network and phylogenetic analyses We performed network and phylogenetic analyses to reconstruct relationships among populations and species. We first distinguished plastid DNA haplotypes and nrITS ribotypes with DnaSP 5.0 software package (Librado & Rozas, 2009). We constructed a median-joining network of haplotypes of the plastid DNA and ribotypes of nrITS with NETWORK 4.5.0.0 (Weir, 1996). In NETWORK, we set both site mutations and indels to evolve with equal likelihood, and we assumed that each indel originated independently of all other indels. For reconstructing phylogenetic relationships, we applied maximum parsimony (MP) and maximum likelihood (ML) to the plastid DNA haplotypes and nrITS ribotypes, independently in PAUP*4.0b10 (Swofford, 2002). We also reconstructed phylogenetic relationships with Bayesian inference (BI) in MrBayes 3.1 (Ronquist & Huelsenbeck, 2003). For the ML analyses and BI, we used a general-time-reversible model (GTR) with a gamma distribution of rate variation among sites (+ G) based on results in PAUP*4.0b10 and jModelTest according to the Akaike information criterion (Posada, 1998; Swofford, 2002). We approximated the gamma distribution using four rate categories. For other parameters, we used default settings. For the MP analyses, we performed heuristic searches with 1000 random addition sequence replicates with ten trees retained at each step during the stepwise additions and with tree-bisection-reconnection for branch swapping. We set PAUP*4.0b10 to treat all characters as unordered and equally weighted, gaps as missing and multistate data as uncertain in the MP analyses. For MP and ML, we calculated bootstrap values from 1000 replicates (Felsenstein, 1985). The BI consisted of two parallel runs with four incrementally heated chains and three million generations sampled every 1000 generations. The output was assessed for convergence using Tracer v.1.3 (Rambaut & Drummond, 2007), and summary statistics and trees were generated using the last two million generations. For network and phylogenetic analyses of plastid DNA and nrITS, we used Dinebra viscida (Scribn.) P.M.Peterson & N.Snow and four samples of Cleistogenes as outgroups, respectively (Peterson et al., 2012, 2016). For the AFLP data, we analysed a rectangular binary 0/1 data matrix using the NTSYS-pc 2.l (Rohlf, 2000) statistical package, and we generated a pairwise similarity matrix with simple matching coefficient according to the SIMQUAL procedure in NTSYS-pc. We performed cluster analysis of the binary matrix using the unweighted pair-group method to build a UPGMA dendrogram by means of SAHN package in NTSYS-pc. We performed bootstrap analysis of the UPGMA tree with 2000 replicates using the winboot computer program (Yap & Nelson, 1996). We also calculated a genetic similarity matrix using AFLP data according to the method of Nei & Li (1979). Dating the divergence between lineages We obtained the nrITS data set of Chloridoideae from Peterson et al. (2016) and used it to estimate the divergence times of Orinus. We verified equal base frequencies among accessions in PAUP*4.0b10 and tested a strict molecular clock hypothesis in MEGA 5 (Tamura et al., 2011). The strict clock was rejected (2lnLR = 434.458, d.f. = 48, P < 0.01). We performed Bayesian divergence time estimation in BEAST 2.4.5 (Bouckaert et al., 2014), which co-estimates topology and node ages. We prepared the nrITS data set for analysis in BEAST using Beauti and by directly editing the content of.xml files. We applied an HKY model of nucleotide substitutions because the more parameter-rich GTR model did not yield substantial gains according to the preliminary analyses. We approximated the gamma distribution of site-specific rates with ten rate categories and applied the ‘birth death’ tree prior process model with a standard, uniform prior. We set a ‘relaxed lognormal clock’ with a substitution rate prior of 6.85 × 10−9 based on Wolfe, Li & Sharp (1987). We calibrated the crown node age of Chloridoideae as 32.74 Mya, which was inferred in prior studies (Christin et al., 2008; Vicentini et al., 2008) and is reasonable because no isotopic surveys to-date provided evidence of an older date for C4 grass (Osborne & Beerling, 2006). Our calibration comprised a normal distribution with an SD of 1.0, yielding a 95% highest posterior density (HPD) of 28.09–37.19 Mya. We constrained the monophyly of the Chloridoideae after preliminary analyses revealed that long branches (i.e. high evolutionary rates) resulted in placement of Sporobolus R.Br. with the outgroups. Sporobolus clearly belongs in Chloridoideae based on strong evidence from prior molecular and morphological studies (Peterson et al., 2016). We ran two independent analyses in BEAST with 5 × 108 Markov chain Monte Carlo generations. We compared the posterior and likelihood distributions from each run in Tracer 1.5.0 to ensure convergence, and we selected one of two runs for downstream analyses (Rambaut & Drummond, 2009). We performed a 10% burn-in and summarized trees by selecting the maximum clade credibility tree in TreeAnnotator 2.4.5, which is part of the BEAST 2.4.5 package, and using median heights for branch lengths. We visualized the maximum clade credibility tree and other summary statistics in FigTree 1.3.1 (Rambaut, 2009). To calculate divergence times of species and populations in Orinus, we used nrITS and two methods: (1) a direct calculation from mutation rates using a strict molecular clock; and (2) a dating analysis in BEAST (Drummond & Rambaut, 2007). For the direct calculation, a strict molecular clock was appropriate, because it was not rejected by an analysis in MEGA 5 (Tamura et al., 2011) (2lnLR = 54.524, d.f. = 29, P = 0.61), and performed a coalescent Bayesian skyline tree process prior suitable for intraspecific data. The strict molecular clock comprised mutation rates of 5.8–8.1 × 10−9 substitutions/site/year, which are well established as rates of evolution of nrITS in other perennial grass genera (Wolfe et al., 1987). For the analyses in BEAST, we analysed 476 sequences of Orinus, and we used four samples of Cleistogenes as outgroups. We applied an age calibration to the Orinus crown node based on our results from the tree for Chloridoideae (above). Specifically, the calibration comprised a uniform prior with using the median date from the tree for Chloridoideae ± 1 Mya. All other parameters of this analysis were the same as for the analysis of Chloridoideae (Peterson et al., 2016), except that we used a coalescent skyline model of diversification because it is appropriate for intraspecific taxa. Glacial refugia and demographic analyses We inferred putative glacial refugia of Orinus based on occurrences of comparatively high genetic diversity among intraspecific populations and for the genus with plastid DNA and nrITS data sets. Genetic diversity is expected to be higher in refugial areas compared with surrounding areas (Hewitt, 1996, 1999, 2000) due especially to founder effects (Mayr, 1942). Although associating high diversity areas with glacial refugia may be overly simplistic (Petit et al., 2003), our intention is to generate a working hypothesis for future testing. To visualize putative refugia based on diversity across geographical areas, we generated surface plots of the number of unique haplotypes, ribotypes and combined haplotypes and ribotypes within populations using a triangulated irregular network (TIN) in ArcMap v.10 (ESRI, 2010). The from-end of link points are used as the TIN uses triangles nodes points, georeferenced populations in this case, and interpolate values, or number of genotypes in this case, between them. We trimmed the TINs according to the convex hulls of O. thoroldii and O. intermedius and for the geographically distinct northern and central metapopulations of O. kokonoricus (Fig. 6A). The resulting surfaces facilitate inference and visualization of geographical areas of high diversity, which we take as possible glacial refugia. We deposited all GIS data online with ESRI, and it is available at https://www.arcgis.com/home/item.html?id=d2f9770771ec4de58a8373d85e3695ce. Figure 6. View largeDownload slide The sampled populations and interpolated distribution of genetic diversity of plastid DNA and nrITS in three Orinus spp. A, the sampled populations and convex hulls drawn around them. B, the interpolated distribution of genetic diversity of plastid DNA. C, the interpolated distribution of genetic diversity of nrITS. D, the highest genetic diversity regions of plastid DNA and nrITS. E, the interpolated distribution of genetic diversity of plastid DNA and nrITS. F, the legend of the study area in China to provide the geographical context, and a–e correspond to the number of genotypes 1–5, respectively. Figure 6. View largeDownload slide The sampled populations and interpolated distribution of genetic diversity of plastid DNA and nrITS in three Orinus spp. A, the sampled populations and convex hulls drawn around them. B, the interpolated distribution of genetic diversity of plastid DNA. C, the interpolated distribution of genetic diversity of nrITS. D, the highest genetic diversity regions of plastid DNA and nrITS. E, the interpolated distribution of genetic diversity of plastid DNA and nrITS. F, the legend of the study area in China to provide the geographical context, and a–e correspond to the number of genotypes 1–5, respectively. We used several methods to detect demographic expansions based on plastid DNA data sets, such as from refugia, in the history of Orinus. We sought to detect recent demographic expansion in Orinus and in each species using a mismatch distribution analysis in Arlequin 3.5. The mismatch distribution analysis applies a model of computes the pairwise p-distances (absolute distances) between sequences and organizes the results in a histogram. A distribution of pairwise frequencies that is multimodal is the equilibrium state, showing expected differentiation among species or populations. In contrast, a unimodal distribution represents geneflow and a trend towards panmixia (Slatkin & Hudson, 1991; Rogers & Harpending, 1992). We analysed the results of the mismatch distribution by comparing the sum of squared deviations observed distributions and those expected following a recent demographic expansion. We also quantified modality of the distribution using Harpending’s raggedness index, where more ragged distributions are multimodal. We assessed recent demographic expansions using Tajima’s D (Tajima, 1989) or Fu’s FS (Fu, 1997) in Arlequin 3.5, with 10000 parametric bootstrap replicates. Positive values for these statistics typically indicate demographic contractions or a recent bottleneck event (Tajima, 1989; Fu, 1997; Hamilton, 2009). In addition, we estimated expansion times (t) using the relationship τ = 2ut (Rogers & Harpending, 1992), where u is the mutation rate per generation for all sequences analysed and t is the expansion time in number of generations. Here, u = μkg, where μ is the substitution rate, k is the average sequence length of the analysed DNA region and g is the generation time in years. We assumed that g = 5 years based on data from the closely related genus Triodia R.Br. (Armstrong, 2011). RESULTS Identification of plastid DNA haplotypes and nrITS ribotypes, and dating divergence times Among the 476 individuals, the total alignment length of three plastid DNA sequences (matK, rbcL and psbA-trnH) was 2898 bp; for nrITS, it was 624 bp. Nucleotide substitutions occurred at 19 polymorphic sites in plastid DNA (0.66%), and 13 of these were potentially parsimony informative; there were no indels (see Supporting Information, Table S2). The aligned nrITS sequences had 46 variable sites, 45 of which were potentially parsimony informative (see Supporting Information, Table S3). We used sequence polymorphisms to identify haplotypes in plastid DNA and ribotypes in nrITS. In plastid DNA, we identified 14 different haplotypes, H1–H14 (Fig. 1B; see Supporting Information, Table S1). Forty-nine populations (55.68%) were fixed for a single haplotype; the remaining 39 (44.32%) were polymorphic (Fig. 1B; see Supporting Information, Table S1). Two haplotypes, H4 and H6, were shared between O. kokonoricus and O. intermedius, and only one, H2, was fixed among all three Orinus spp. Of these shared haplotypes, H2 and H4 occurred in 67.05 and 17.05% of all populations, respectively, whereas H6 was limited to six populations (6.82%). H2 was widely distributed across the QTP, and H4 was widely distributed throughout the eastern and south-eastern QTP. Six haplotypes were found only in O. thoroldii, four were specific to O. kokonoricus and one was specific to O. intermedius (Fig. 1C; see Supporting Information, Fig. S1). In nrITS, we identified 30 different ribotypes, S1–S30. Of these, ten were each limited to a single population, and all but one were restricted to one species (Fig. 2; see Supporting Information, Fig. S2, Table S1). The S10 ribotype was shared between O. kokonoricus and O. intermedius. Thirty-nine populations (44.32%) only possessed one haplotype; the remaining 49 (55.68%) were polymorphic (Fig. 2; see Supporting Information, Fig. S2, Table S1). Ribotypes S2, S19 and S26 were the most common nrITS sequences and were fixed or present in nearly all sampled populations from the three Orinus spp. Sequence S2 was discovered in 84.44% of all samples for O. kokonoricus, S19 was found in all populations from O. thoroldii and S26 also occurred in 86.67% of all sampled populations of O. intermedius. Ten ribotypes (S7–S9, S15–S18, S22, S23 and S25) were private to a single population (Fig. 2; see Supporting Information, Table S1). Figure 2. View largeDownload slide Sampling sites and geographical distribution of nrITS ribotypes observed in three Orinus species. A, sampling region (dashed-line box). B, documented distribution range (dotted line), sampling sites and nrITS ribotype frequencies per population. Each pie graph represents nrITS ribotypes; size of pie is proportional to number of individuals bearing that haplotype. Figure 2. View largeDownload slide Sampling sites and geographical distribution of nrITS ribotypes observed in three Orinus species. A, sampling region (dashed-line box). B, documented distribution range (dotted line), sampling sites and nrITS ribotype frequencies per population. Each pie graph represents nrITS ribotypes; size of pie is proportional to number of individuals bearing that haplotype. The divergence time estimation showed Orinus diverged from its closest relatives c. 12.06 Mya (95% HPD: 6.87–17.50) (Miocene) and diversified c. 2.85 Mya (node 1; 95% HPD: 0.78–5.05) in the Pliocene (Fig. 3; see Supporting Information, Table S4). In Orinus, modern ribotypes diversified in the Pleistocene, 0.31–1.71 Mya (nodes 2–9) (see Supporting Information, Fig. S2, Table S4), according to the analysis using BEAST. When the divergence time of ribotypes was calculated using mutation rates, we found the similar results as the above, also in the Pleistocene (see Supporting Information, Table S4). Figure 3. View largeDownload slide Bayesian divergence time estimate of Orinus (shaded) based on the nrITS data set of Chloridoideae from Peterson et al. (2016). Clade constraint is marked by the black circle. Numbers in parentheses indicate 95% highest posterior density intervals. Figure 3. View largeDownload slide Bayesian divergence time estimate of Orinus (shaded) based on the nrITS data set of Chloridoideae from Peterson et al. (2016). Clade constraint is marked by the black circle. Numbers in parentheses indicate 95% highest posterior density intervals. Population variability, glacial refugia and demography Haplotype diversity (HE) from each population was calculated based on plastid DNA haplotype frequencies (see Supporting Information, Table S1). The average genetic diversity within populations (HS) was relatively low for O. thoroldii, O. kokonoricus and O. intermedius (0.16 ± 0.05, 0.23 ± 0.07 and 0.31 ± 0.10), and the total genetic diversity HT (0.74 ± 0.03) across the whole genus was much higher than HS (0.21 ± 0.04) based on the observed plastid DNA variations (Table 1). A permutation test suggested that population differentiation across the sampling range was high (NST = 0.83 ± 0.04) and was significantly higher than GST (0.72 ± 0.05) (0.01 < P < 0.05), indicating a strong phylogeographic structure across the entire sampling range (Table 1). Table 1. Estimates of genetic diversity for plastid DNA in three Orinus spp. on the Qinghai-Tibet Plateau Groups  HS  HT  GST  NST  O. kokonoricus (P1–45)  0.23 (0.07)  0.74 (0.06)  0.69 (0.10)  0.60 (0.12)  O. thoroldii (P46–73)  0.16 (0.05)  0.56 (0.04)  0.72 (0.08)  0.81 (0.07)*  O. intermedius (P74–88)  0.31 (0.10)  0.59 (0.05)  0.48 (0.18)  0.33 (0.15)  All populations  0.21 (0.04)  0.74 (0.03)  0.72 (0.05)  0.83 (0.04)*  Groups  HS  HT  GST  NST  O. kokonoricus (P1–45)  0.23 (0.07)  0.74 (0.06)  0.69 (0.10)  0.60 (0.12)  O. thoroldii (P46–73)  0.16 (0.05)  0.56 (0.04)  0.72 (0.08)  0.81 (0.07)*  O. intermedius (P74–88)  0.31 (0.10)  0.59 (0.05)  0.48 (0.18)  0.33 (0.15)  All populations  0.21 (0.04)  0.74 (0.03)  0.72 (0.05)  0.83 (0.04)*  HS, estimate of average genetic diversity within populations; HT, total genetic diversity; GST, coefficient of genetic variation over all populations; NST, coefficient of genetic variation influenced by both haplotype frequencies and genetic distance between haplotypes. *P < 0.05. View Large Hierarchical AMOVA based on the plastid DNA data set revealed that there was higher interspecific genetic differentiation across the whole distribution of this genus (FST = 0.82), and the interpopulation genetic differentiation in O. thoroldii and O. kokonoricus was much more pronounced (FST = 0.72, FST = 0.73), but there was a little interpopulation differentiation in O. intermedius (FST = 0.34) (Table 2). These results suggested that the grouping pattern of hierarchical AMOVA was basically consistent with relationships among plastid haplotypes and their geographical distribution (Fig. 1C, Table 2). Table 2. Results of analysis of molecular variance for plastid DNA and nrITS sequence data from populations of Orinus on the Qinghai-Tibet Plateau Groups  Source of variation  d.f.  SS  VC  Variation (%)  Fixation index    Plastid DNA  All samples  Among groups  2  121.90  0.74 Va  44.78  FST = 0.82*  Among populations within groups  53  152.13  0.62 Vb  37.72  FSC = 0.68*  Within populations  420  52.46  0.29 Vc  17.50  FCT = 0.45*  Total  475  326.49  1.65      O. kokonoricus  Among populations  20  50.09  0.54 Va  72.69  FST = 0.73*  Within populations  161  14.17  0.20 Vb  27.31    Total  181  64.26  0.74      O. thoroldii  Among populations  26  92.52  0.92 Va  71.90  FST = 0.72*  Within populations  161  24.18  0.36 Vb  28.10    Total  187  116.70  1.28      O. intermedius  Among populations  7  9.52  0.16 Va  34.27  FST = 0.34*  Within populations  98  14.10  0.32 Vb  65.73    Total  105  23.62  0.48        nrITS  All samples  Among groups  2  1977.08  6.32 Va  81.90  FST = 0.96*  Among populations within groups  53  491.63  1.08 Vb  13.96  FSC = 0.77*  Within populations  420  134.40  0.32 Vc  4.15  FCT = 0.82*  Total  475  2603.11  7.72      O. kokonoricus  Among populations  20  204.28  1.13 Va  67.89  FST = 0.68*  Within populations  161  86.25  0.54 Vb  32.11    Total  181  290.53  1.67      O. thoroldii  Among populations  26  21.58  0.09 Va  30.17  FST = 0.30*  Within populations  161  33.48  0.21 Vb  69.83    Total  187  55.06  0.30      O. intermedius  Among populations  7  265.77  2.95 Va  95.17  FST = 0.95*  Within populations  98  14.67  0.15 Vb  4.83    Total  105  280.43  3.10      Groups  Source of variation  d.f.  SS  VC  Variation (%)  Fixation index    Plastid DNA  All samples  Among groups  2  121.90  0.74 Va  44.78  FST = 0.82*  Among populations within groups  53  152.13  0.62 Vb  37.72  FSC = 0.68*  Within populations  420  52.46  0.29 Vc  17.50  FCT = 0.45*  Total  475  326.49  1.65      O. kokonoricus  Among populations  20  50.09  0.54 Va  72.69  FST = 0.73*  Within populations  161  14.17  0.20 Vb  27.31    Total  181  64.26  0.74      O. thoroldii  Among populations  26  92.52  0.92 Va  71.90  FST = 0.72*  Within populations  161  24.18  0.36 Vb  28.10    Total  187  116.70  1.28      O. intermedius  Among populations  7  9.52  0.16 Va  34.27  FST = 0.34*  Within populations  98  14.10  0.32 Vb  65.73    Total  105  23.62  0.48        nrITS  All samples  Among groups  2  1977.08  6.32 Va  81.90  FST = 0.96*  Among populations within groups  53  491.63  1.08 Vb  13.96  FSC = 0.77*  Within populations  420  134.40  0.32 Vc  4.15  FCT = 0.82*  Total  475  2603.11  7.72      O. kokonoricus  Among populations  20  204.28  1.13 Va  67.89  FST = 0.68*  Within populations  161  86.25  0.54 Vb  32.11    Total  181  290.53  1.67      O. thoroldii  Among populations  26  21.58  0.09 Va  30.17  FST = 0.30*  Within populations  161  33.48  0.21 Vb  69.83    Total  187  55.06  0.30      O. intermedius  Among populations  7  265.77  2.95 Va  95.17  FST = 0.95*  Within populations  98  14.67  0.15 Vb  4.83    Total  105  280.43  3.10      d.f., degrees of freedom; SS, sum of squares; VC, variance components; FCT, correlation within populations relative to groups; FSC, correlation of haplotypes within groups relative to the total; FST, correlation within populations relative to the total. *P < 0.001 (10000 permutations). View Large Our analyses revealed extensive gene flow among populations of Orinus across its geographical range. In the STRUCTURE analyses (Fig. 4A), the highest likelihood of the plastid DNA data was obtained when all samples were weakly clustered into three groups according to the recognized species (K = 3). However, this result did not correspond to separate geographical regions because there was strong gene flow across the range of Orinus. In addition, a Nem value of 0.61 also indicated extensive gene flow existed between groups based on plastid DNA data. Figure 4. View largeDownload slide Estimated genetic structure for K = 2–3 obtained with the program STRUCTURE (Pritchard et al., 2000) for 56 populations of Orinus based on plastid DNA (A) and nrITS (B) sequence variations. Each vertical bar represents an individual and its assignment proportion (not shown). Red, Orinus thoroldii; Green, O. kokonoricus; Blue-purple, O. intermedius. Figure 4. View largeDownload slide Estimated genetic structure for K = 2–3 obtained with the program STRUCTURE (Pritchard et al., 2000) for 56 populations of Orinus based on plastid DNA (A) and nrITS (B) sequence variations. Each vertical bar represents an individual and its assignment proportion (not shown). Red, Orinus thoroldii; Green, O. kokonoricus; Blue-purple, O. intermedius. For nrITS, we found that 81.90% of the total genetic variation occurred between groups, 13.96% occurred among populations within groups and only 4.15% occurred within populations. The pairwise FST value between three groups was 0.96 (Table 2). It suggested that there was high interspecific differentiation among three species, and interpopulation differentiation of O. kokonoricus and O. intermedius was also extremely high (FST = 0.68, FST = 0.95). In contrast, O. thoroldii exhibited lower interpopulation variation (FST = 0.30). The grouping pattern of hierarchical AMOVA based on nrITS data was in accordance with relationships of ribotypes and their geographical distribution (Fig. 5, Table 2). Similarly, the STRUCTURE analyses indicated that almost all samples of O. intermedius and O. kokonoricus were clustered into a lineage, and there was only obscure gene flow between them and O. thoroldii (K = 2) (Fig. 4B). However, all samples of this genus were clustered into three independent lineages when the highest likelihood of the data was obtained (K = 3), and there was prominent introgression between O. kokonoricus and O. intermedius (Fig. 4B). In addition, the gene flow between the three groups was low according to the Nem value was 0.11. Figure 5. View largeDownload slide Hypothesis of incomplete lineage sorting and reticulation in nrITS of Orinus shown within a species tree of the genus. Lineage A is fixed in O. thoroldii, whereas descendent types from lineage B are incompletely sorted in O. kokonoricus and O. intermedius. nrITS type S10 evolved in lineage F of O. kokonoricus and was later shared with O. intermedius during a secondary contact event. Black nodes indicate divergence events among types that we can infer from our data, and extant types sampled in this study are numbered 1–30 at terminal nodes. Two particularly large clades of nrITS types in lineage 4 are collapsed into single terminal nodes, and all types represented by those nodes are listed numerically. Figure 5. View largeDownload slide Hypothesis of incomplete lineage sorting and reticulation in nrITS of Orinus shown within a species tree of the genus. Lineage A is fixed in O. thoroldii, whereas descendent types from lineage B are incompletely sorted in O. kokonoricus and O. intermedius. nrITS type S10 evolved in lineage F of O. kokonoricus and was later shared with O. intermedius during a secondary contact event. Black nodes indicate divergence events among types that we can infer from our data, and extant types sampled in this study are numbered 1–30 at terminal nodes. Two particularly large clades of nrITS types in lineage 4 are collapsed into single terminal nodes, and all types represented by those nodes are listed numerically. Our surface analysis in ArcMap revealed four centres of diversity for haplotypes, ribotypes and the combined diversity of haplotypes and ribotypes (Fig. 6). For O. thoroldii, the haplotype and ribotype data are congruent in showing a centre of diversity in the eastern part of the range (Fig. 6B, C). The haplotype and ribotype data also agree that the northern metapopulation of O. kokonoricus has a centre of diversity at the southern portion of its range (Fig. 6B, C). In contrast, the haplotypes show centres of diversity at the northern ends of the ranges of O. intermedius and the central metapopulation of O. kokonoricus, whereas the ribotypes show centres of diversity at the southern ends of these ranges (Fig. 6B, C). Overall, greater diversity among ribotypes than haplotypes (Figs 1, 2) yielded a geographical pattern in the combined analysis of genotypes (Fig. 6E) consistent with ribotypes alone. The mismatch distributions of pairwise differences for all samples from the identified haplotypes of Orinus and O. thoroldii were close to bimodal, whereas those of O. kokonoricus and O. intermedius were close to be unimodal (see Supporting Information, Fig. S3), indicating that the former should experience recent bottlenecks and the latter should experience recent and rapid population expansion. However, further analyses using the non-significant variance and raggedness index tests indicated that none of the observed distributions differed significantly from those expected under a sudden expansion model (Table 3). Moreover, the significant negative calculated values for Tajima’s D (P < 0.01) and non-significant positive calculated values for Fu’s FS (P > 0.05) were also supported this earlier, rapid demographic or geographical expansion and subsequent population growth (Table 3). Besides, the mean expansion times were estimated to be c. 25, 29, 40 and 53 kya for O. intermedius, O. kokonoricus, O. thoroldii and Orinus, respectively (Table 3). These results support a recent and rapid demographic population growth in the genus and each of the three species, and all expansions occurred before the Last Glacial Maximum (LGM; c. 20 kya). Possibly, Orinus and O. thoroldii have experienced recent bottlenecks. Table 3. Results of mismatch distribution analyses and neutrality tests (Tajima’s D and Fu’s FS) based on cpDNA data   Mismatch distribution  Neutrality tests  Groups  τ (95% CI)  SSD (P value)  RAG (P value)  tmin–tmax (kya)  tave (kya)  Tajima’s D (P value)  Fu’s Fs (P value)  O. thoroldii  0.772 (0.046–6.021)  0.057 (0.079)  0.160 (0.208)  26.639–88.797  40.983  −0.239 (0.001)  0.393 (0.000)  O. kokonoricus  0.554 (0.050–1.261)  0.050 (0.099)  0.170 (0.217)  19.117–63.722  29.410  −0.197 (0.004)  0.342 (0.006)  O. intermedius  0.480 (0.047–1.331)  0.025 (0.357)  0.196 (0.513)  16.563–55.210  25.482  −0.247 (0.006)  0.465 (0.014)  Total  0.994 (0.050–3.577)  0.050 (0.128)  0.169 (0.254)  34.300–114.332  52.768  −0.224 (0.001)  0.384 (0.002)    Mismatch distribution  Neutrality tests  Groups  τ (95% CI)  SSD (P value)  RAG (P value)  tmin–tmax (kya)  tave (kya)  Tajima’s D (P value)  Fu’s Fs (P value)  O. thoroldii  0.772 (0.046–6.021)  0.057 (0.079)  0.160 (0.208)  26.639–88.797  40.983  −0.239 (0.001)  0.393 (0.000)  O. kokonoricus  0.554 (0.050–1.261)  0.050 (0.099)  0.170 (0.217)  19.117–63.722  29.410  −0.197 (0.004)  0.342 (0.006)  O. intermedius  0.480 (0.047–1.331)  0.025 (0.357)  0.196 (0.513)  16.563–55.210  25.482  −0.247 (0.006)  0.465 (0.014)  Total  0.994 (0.050–3.577)  0.050 (0.128)  0.169 (0.254)  34.300–114.332  52.768  −0.224 (0.001)  0.384 (0.002)  τ, time in number of generations elapsed since the sudden expansion episode; SSD, sum of squared deviations; RAG, Harpending’s raggedness index; tmin, absolute expansion time estimated with a high substitution rate, 1.0 × 10−9 substitutions/site/year (s/s/year); tmax, absolute expansion time estimated with a low substitution rate, 0.3 × 10−9 s/s/year; tave, mean expansion time; kya, thousand years ago; CI, confidence interval. View Large DISCUSSION Phylogeographic patterns in Orinus In recent years, many studies have focused on the early stages of speciation (Abbott et al., 2013; Liu et al., 2014b). Some of these have shown that divergence among species can occur even in the presence of gene flow and that species boundaries are maintained by divergent selection (Hewitt, 1996; Petit & Excoffier, 2009; Schluter, 2009; Nosil, 2012; Sousa & Hey, 2013). Our plastid DNA, nrITS and AFLP data suggest that Orinus may be at an advanced stage of speciation in which gene flow is limited, by geographical isolation or divergent selection, and fixation is underway so that most sequences form clades comprising only one species, but these species are paraphyletic (Fig. 5) (Funk & Omland, 2003) Population bottlenecks have also played a role in promoting divergence among closely related taxa but has rarely been studied (Butlin et al., 2012; Bi et al., 2016). Our results suggest a subtle but complex geographical pattern in the genetic diversity of Orinus. High levels of interspecific genetic variation occurred among samples of plastid DNA (44.78%) and nrITS (81.90%) (Table 2). In addition, the differentiation among species was high overall evidenced by significant FST of 0.86 and 0.96 for plastid DNA and nrITS (Table 2), indicating nearly complete genetic isolation (Rieseberg, Church & Morjan, 2004; Hauquier et al., 2017). Moreover, only one plastid DNA haplotype and no nrITS types were shared among all populations (Figs 1, 2). However, the pattern is more complex than isolation by species; the eastern and south-eastern species, O. kokonoricus and O. intermedius, clearly have more gene flow between them than either does with the western species, O. thoroldii. Orinus kokonoricus and O. intermedius share two additional plastid DNA haplotypes in common (Fig. 1), exclusive of O. thoroldii, and one nrITS type (Fig. 2), strongly suggestive of recent or ongoing gene flow (Parchman, Benkman & Britch, 2006). Moreover, the phylogenetic reconstruction of nrITS sequences in O. kokonoricus and O. intermedius shows that the haplotypes unique to each species do not form clades; they are intermixed in their phylogenetic positions, and this is highly consistent with genetic connectivity and subsequent isolation and fixation. In contrast, O. thoroldii does not share any additional haplotypes or nrITS ribotypes in common with O. kokonoricus or O. intermedius. Overall, the demographic pattern suggests greater and/or more frequent connectivity between the eastern and south-eastern species than with the western one. It is noteworthy that Orinus spp. are largely differentiated and possess few shared ribotypes or haplotypes, but the ribotypes or haplotypes do not form the same clades as the species (Fig. 5; see Supporting Information, Figs S1, S2). At least two, non-mutually exclusive hypotheses may explain this phylogeographic pattern: (1) somewhat recent species divergence with incomplete lineage sorting (Choleva et al., 2014; Zhou et al., 2017) or (2) recent and/or repeated secondary contact and introgression (Avise, 2000; Sardell & Uy, 2016). In Orinus, both mechanisms are needed to explain our results (Figs 1, 2, 4). Based on our nrITS data, we propose nrITS split into two lineages, A and B before or after the divergence of O. thoroldii from the rest of Orinus (Fig. 5). Lineage A became fixed in O. thoroldii, in which lineage sorting is complete (Fig. 5). Lineage B diversified into lineages C, D, E and F in an ancestor of O. kokonoricus and O. intermedius (Fig. 5). Lineages C and E are limited to O. intermedius, whereas D and F occur only in O. kokonoricus (Fig. 5). Thus, sorting of the descendent genotypes of lineage B is incomplete between O. kokonoricus and O. intermedius. The S10 nrITS type is shared by both O. kokonoricus and O. intermedius, but is clearly derived in lineage F (Fig. 5). Therefore, it seems probable that the S10 type was shared by O. kokonoricus with O. intermedius during a secondary contact and hybridization event (Sardell & Uy, 2016; Zhou et al., 2017). The plastid DNA data and the AFLPs are not in conflict with this hypothesis. Plastid DNA loci sort and become fixed more slowly than nuclear markers (Wolfe et al., 1987; de Villiers et al., 2013; Shaw et al., 2017), such as nrITS, and the shared haplotypes among all three species and between O. kokonoricus and O. intermedius may represent incomplete lineage sorting (Choleva et al., 2014; Yasuda et al., 2015; Kuritzin et al., 2016). Alternatively, the haplotype shared among all three species may represent an older secondary contact (Miraldo et al., 2011; Sardell & Uy, 2016), predating the S10 event and too old to be preserved in nrITS data. AFLPs revealed that species and populations within species comprise clades. Therefore, most but not all lineage sorting among species is complete (Choleva et al., 2014; Yasuda et al., 2015; Kuritzin et al., 2016). The demographic pattern that we observed in Orinus seems highly consistent with the history of geomorphism in the QTP and adjacent areas. The evidence for QTP uplifts comes primarily from detecting its effects on regional aridity and the monsoon climatic using dust sedimentation and indicator fossils of marine foraminifera (e.g. An et al., 2001; Sun & An, 2001; Guo et al., 2002; Zheng et al., 2004). Although the exact timing of uplift events remains somewhat in question, the QTP probably underwent early orogenesis c. 40 Mya along its western margin as the Indian subcontinent collided with Eurasia and more extensive, rapid uplift 20–7 Mya (Harrison et al., 1992; Ruddiman, 1998). In addition, many geological and historical biogeographic studies infer that that additional uplifts occurred east of the QTP around c. 4.0–2.5 Mya in adjacent mountain ranges, especially the Hengduan Mountains (Chen, 1992, 1996; Zheng et al., 2004; Gao et al., 2013; Favre et al., 2015). Thus, the Hengduan Mountains and other adjacent ranges have experienced more recent, Pliocene geological changes and associated climatic perturbations than the QTP proper (Xing & Ree, 2017). Contemporaneously, glacial cycling also impacted the area, and its effects can be difficult to distinguish from the genetic signature of geomorphism (Shackleton & Opdyke, 1973; Liu et al., 2014b). Nevertheless, it is possible that regional geomorphism to the east of the QTP facilitated the dynamic history of recent allopatric speciation and secondary contact suggested by our nrITS data for O. kokonoricus and O. intermedius. Similarly, the relatively recent geomorphism to the east of the QTP by comparison to the QTP may also be reflected in the evolutionary and demographic histories of other diverse plant groups. For example, Yang et al. (2012) found that species of Meconopsis Vig. (Papaveraceae) occurring in the QTP diverged earlier than species occurring in adjacent eastern areas. Allopatric divergence, glacial refugia and historical demography The role of climatic oscillations in shaping the regional biota is important and should not be ruled out as a driver of evolutionary events in Orinus. Other species in the region have histories affected by climatic change. For example, Anisodus tanguticus Pascher (Solanaceae) has western and eastern clades that diverged during the Pleistocene, probably a result of genetic isolation in glacial refugia, rather than a geological event (Wan et al., 2016). Similar patterns of intraspecific divergence due to Quaternary climatic change have also been found in other plant species in the QTP and adjacent regions, such as Aconitum gymnandrum Maxim. (Ranunculaceae) (Wang et al., 2009a), Sophora davidii (Franch.) Skeels (Fabaceae) (Fan et al., 2013) and Oxyria sinensis Hemsl. (Polygonaceae) (Meng et al., 2015). Based on the divergence time of 32.74 Mya for Chloridoideae, the crown age was dated to be 2.85 (95% HPD: 0.58–12.45) Mya, with interspecific divergences occurred in the range of 0.31–1.70 Mya (see Supporting Information, Fig. S2, Table S4). This relatively recent interspecific diversification atop an unbranched stem lineage is sometimes regarded as reflecting ancient extinctions (Harvey, May & Nee, 1994). Deep divergence in Orinus apparently began in the mid-Pleistocene (see Supporting Information, Fig. S2, Table S4). Although these estimated timescales for deep divergence must be interpreted with caution, they (the major nodes) correspond relatively well with recent QTP uplifts (0.01–1.80 Mya) (Li, Shi & Li, 1995; Shi, Li & Li, 1998; Zheng, Xu & Shen, 2002). Development of aridification in the QTP was a response to climatic change at the time, such as long-term cooling and drying trends. By the early Pleistocene, the QTP dramatically uplifted to present heights (Wu et al., 2008), which made the climate of this area become extremely cool and dry. Especially, strong episodic cooling resulted in sharp increases in aridity of the QTP regions (Fang et al., 2002). Aridification played a significant role in the increase of genetic diversity, and even species divergence and speciation of alpine plants. Habitat fragmentation, vicariance and population isolation on the QTP provided opportunities for allopatric speciation through the action of selection and/or genetic drift (Wang et al., 2013; Meng et al., 2014; Wen et al., 2014). Thus, deep intraspecific divergence in Orinus may be a result of adaptation to arid habitats, with the mechanisms probably being similar to those in other arid/desert species (Wang et al., 2013; Yu et al., 2013; Liu et al., 2014b; Meng et al., 2014). Therefore, it is likely that the deep divergence in Orinus was caused by allopatric differentiation and reduction in gene flow following plateau uplifts and climatic oscillation. Plant populations in refugia generally have high genetic divergence (Petit et al., 2003). Our surface analysis in ArcMap revealed four centres of diversity for haplotypes, ribotypes and the combined diversity of haplotypes and ribotypes, which were consistent with the result of glacial refugia based on both plastid DNA and nrITS data sets (Figs 1B, 2B, 6B–E). Among them, there were three areas with relatively high plastid DNA and nrITS genetic diversity and one region with unique plastid DNA haplotypes and high nrITS genotypes (Figs 1B, 2B, 6B–E). One region with plastid DNA haplotype uniqueness and high nrITS genetic diversity (O. intermedius) is located on the southern edge of the QTP, an area well known as an important refugium for other alpine plant species (Figs 1B, 2B, 6E) (Zhang et al., 2005; Meng et al., 2007; Chen et al., 2008; Yang et al., 2008; Wang et al., 2009a; Gao et al., 2012; Liu et al., 2012, 2014b; Ma et al., 2014), which was further supported by the result of surface analysis in ArcMap (Fig. 6C–E). Furthermore, three (S28–S30) of four (S10, S28–S30) nrITS ribotypes from this region were endemic, indicating a divergence centre for Orinus. Plastid and nrITS genetic diversities are higher in two regions (O. kokonoricus) in the north-eastern and eastern plateau (Figs 1B, 2B, 6E), which may represent two important refugia throughout the Quaternary. Our surface analysis in ArcMap also indicated the northern metapopulation of O. kokonoricus has a centre of diversity at the southern portion of its range (Fig. 6B, C), and the southern metapopulation of O. kokonoricus has another centre of diversity at the southern portion of its range (Fig. 6E). We have two alternative biogeographic hypotheses regarding high levels of genetic diversity in O. kokonoricus. A glacial refugium allowed O. kokonoricus to survive climatic oscillations and to accumulate genetic diversity (Tzedakis et al., 2002) or this region was composed of a mixture of individuals with different origins, resulting in higher genetic diversity than the original source (Petit et al., 2003). If the latter hypothesis is true, then haplotypes detected in this region should form a subset of haplotypes of the multiple original sources, although four (H1, H3, H5 and H7) out of the seven plastid DNA haplotypes in these populations of O. kokonoricus are endemic. Furthermore, the majority of nrITS ribotypes (S1–S9 and S11–S18) were restricted to this area, indicating other divergence centres of O. kokonoricus. In addition, widespread haplotypes (H2 and H4) and ribotypes (S2 and S3) coexisting among the great number of plastid DNA haplotypes and nrITS ribotypes also occurred in these two regions. Thus, the most-parsimonious explanation is consequently the presence of glacial refugia for O. kokonoricus. A similar result was found in previous studies of endemics species from the QTP (e.g. Gao et al., 2012; Ma et al., 2014). In addition, the western region of the QTP was only occupied by O. thoroldii probably represent the fourth separate refugia during the large glaciation stage (Figs 1, 2, 6B–E). The populations of O. thoroldii in eastern part of the range had relatively high plastid DNA and nrITS genetic diversity (Figs 1, 2), which were congruent in showing a centre of diversity according to the surface analysis in ArcMap (Fig. 6B–E). Therefore, we thought greater diversity among ribotypes than haplotypes yielded a geographical pattern in the combined analysis of genotypes consistent with ribotypes alone. Based on the plastid DNA data set, the three Orinus spp. might have experienced recent and repeated population expansion after allopatric divergence (e.g. Avise, 2004; Wang et al., 2009b; Gao et al., 2012; Ma et al., 2014; Liu et al., 2015), and O. thoroldii might also have experienced recent bottlenecks due to its bimodal mismatch distribution (Fig. S3B). A star-like phylogeny pattern for O. thoroldii and O. kokonoricus is shown in the haplotype network (Fig. 1C), suggesting the haplotypes might have originated from a sudden expansion event (Hudson, 1990). This is supported by significantly negative Tajima’s D values (Table 3), as well as unimodal mismatch distributions for O. kokonoricus and O. intermedius (see Supporting Information, Fig. S3C, D). Overall, these findings support expansions of Orinus in western, eastern and north-eastern regions, perhaps between 25.48 and 52.77 kya, prior to the LGM (18–24 kya) (Shi et al., 1998) (Table 3), but later than other previously reported QTP endemic genera (Ma et al., 2014) and species (e.g. Yang et al., 2008; Jia et al., 2011; Li et al., 2012). Early expansion of Orinus would have yielded all haplotypes in the QTP region. However, some haplotypes are fixed in specific populations. For example, H3, H4 and H6 correspond to populations 3–5, 42–43 and 17–20, respectively. Thus, the existence of multiple high-elevation refugia that would have facilitated survival during the LGM (Wang et al., 2009c; Opgenoorth et al., 2010; Jia et al., 2011; Gao et al., 2012; Liu et al., 2012; Ma et al., 2014). Perhaps, only small and isolated populations, such as P87 and P88, survived in restricted ice-free regions, and range expansion in the refugia appears to be restricted during glacial and interglacial periods. However, another distinct genetic signature via large-scale range expansions, that is wide distribution of a single haplotype (Avise, 1987; Comes & Kadereit, 1998; Hewitt, 2000). In the present study, two widely distributed haplotypes (H2, H9) were found in Orinus. For instance, H2 was found as a single haplotype in 59 of the 88 populations of the QTP and fixed in 25 of them, whereas H9 was fixed in 11 of the 28 populations in O. thoroldii from the western QTP (e.g. Jia et al., 2011; Li et al., 2012; Wang et al., 2014; Liu et al., 2015). Therefore, Orinus might have experienced a second recent range expansion, probably after the LGM, which led to the wide distribution of the two haplotypes in the QTP region. Orinus intermedius is a newly described species (Su et al., 2017) and is sister to O. kokonoricus (Figs 1C, 3, 5). The morphological characters of O. intermedius are intermediate between its two congeners (Su et al., 2017). We initially hypothesized that O. intermedius may have formed via hybrid speciation on the QTP and its adjacent areas (Su et al., 2015, 2017) in a similar manner to Ostryopsis intermedia B.Tian & J.Q.Liu (Betulaceae) (Liu et al., 2014a), Picea purpurea Masters (Pinaceae) (Sun et al., 2014), Roscoea humeana Balf.f. & W.W.Sm. and R. cautleoides Gagnep. (Zingiberaceae) (Zhao et al., 2016). However, our plastid and nrITS data support recent divergences between O. intermedius and O. kokonoricus, which may explain their morphological similarities. In addition, the AFLP data suggest that the three Orinus spp. each form a separate clade indicating phylogenetic distinctiveness (Supporting Information, Fig. S4). As the next steps, we plan to use more rapidly evolving markers, such as single- or low-copy nuclear genes, or even next- generation sequencing approaches in Orinus (Zimmer & Wen, 2012, 2015) to unravel the genetic structure and phylogeographic patterns and to explore speciation mechanisms (Abbott, 2017; Crawford & Archibald, 2017). SUPPORTING INFORMATION Additional Supporting Information may be found in the online version of this article at the publisher’s web-site: Figure S1. Strict MP consensus tree of 14 plastid DNA haplotypes of three Orinus spp. Each tree is topologically congruent with ML tree and Bayesian consensus tree based on same data set. Bootstrap support values from Bayesian and MP analyses, and corresponding posterior probabilities from MP analyses, are given above branches (values > 50%). Dinebra viscida was used as an outgroup. Figure S2. Strict MP consensus tree of 30 nrITS ribotypes of three Orinus species. Each tree is topologically congruent with ML tree and Bayesian consensus tree based on same data set. Bootstrap support values from Bayesian and ML analyses, and corresponding posterior probabilities from MP analyses, are given above branches (values > 50%). Four Cleistogenes spp. (C. bulgarica, C. festucacea, C. mucronata and C. songorica) were used as outgroups. Figure S3. Mismatch distributions for plastid DNA data in genus Orinus (A), O. thoroldii (B), O. kokonoricus (C) and O. intermedius (D). The thin line represents the expected (Exp) values, and the dashed line represents the observed (Obs) values. Figure S4. UPGMA dendrograms for three Orinus spp. generated by cluster analysis based on AFLP data set. Table S1. Sampling data, estimates of haplotype diversity (HE), and haplotype composition for plastid DNA and nrITS data sets from 88 populations of Orinus. Table S2. Variable sites for 14 plastid DNA haplotypes of three Orinus spp. on the Qinghai-Tibet Plateau (QTP). Table S3. Variable sites for 30 nrITS ribotypes of three Orinus spp. on the Qinghai-Tibet Plateau (QTP). Table S4. Estimates of divergence times in Orinus for the major nodes based on nrITS data set. ACKNOWLEDGEMENTS We are grateful to Prof. Jianquan Liu for providing space and expendables in the laboratory and two anonymous reviewers for their comments on the manuscript. This work was supported by grants from the National Natural Science Foundation of China (31260052 and 41761009), the Natural Science Foundation of Qinghai Province (2017-ZJ-904, 2016-ZJ-937Q and 2014-ZJ-947Q) and the State Scholarship Fund of China Scholarship Council (201508635003). REFERENCES Abbott RJ. 2017. Plant speciation across environmental gradients and the occurrence and nature of hybrid zones. Journal of Systematics and Evolution  55: 238– 258. Google Scholar CrossRef Search ADS   Abbott R, Albach D, Ansell S, Arntzen JW, Baird SJ, Bierne N, Boughman J, Brelsford A, Buerkle CA, Buggs R, Butlin RK, Dieckmann U, Eroukhmanoff F, Grill A, Cahan SH, Hermansen JS, Hewitt G, Hudson AG, Jiggins C, Jones J, Keller B, Marczewski T, Mallet J, Martinez-Rodriguez P, Möst M, Mullen S, Nichols R, Nolte AW, Parisod C, Pfennig K, Rice AM, Ritchie MG, Seifert B, Smadja CM, Stelkens R, Szymura JM, Väinölä R, Wolf JB, Zinner D. 2013. Hybridization and speciation. Journal of Evolutionary Biology  26: 229– 246. Google Scholar CrossRef Search ADS PubMed  Abbott RJ, Smith LC, Milne RI, Crawford RM, Wolff K, Balfour J. 2000. Molecular analysis of plant migration and refugia in the Arctic. Science  289: 1343– 1346. Google Scholar CrossRef Search ADS PubMed  An ZS, Kutzbach JE, Prell WL, Porter SC. 2001. Evolution of Asian monsoons and phased uplift of the Himalaya-Tibetan plateau since late Miocene times. Nature  411: 62– 66. Google Scholar CrossRef Search ADS PubMed  Armstrong G. 2011. Evidence for the equal resilience of Triodia spp. (Poaceae), from different functional groups, to frequent fire dating back to the late Pleistocene. Heredity  107: 558– 564. Google Scholar CrossRef Search ADS PubMed  Avise JC. 1987. Identification and interpretation of mitochondrial DNA stocks in marine species. In: Kumpf H, Nakamura EL, eds. Proceedings of the Stock Identification Workshop . Panama City: National Oceanographic and Atmospheric Administration, 105– 136. Avise JC. 2000. Phylogeography: the history and formation of species . Cambridge: Harvard University Press. Avise JC. 2004. Molecular markers, natural history and evolution, 2nd edn . Sunderland: Sinauer Associates. Bawa KS, Koh LP, Lee TM, Liu J, Ramakrishnan PS, Yu DW, Zhang YP, Raven PH. 2010. Ecology. China, India, and the environment. Science  327: 1457– 1459. Google Scholar CrossRef Search ADS PubMed  Bi H, Yue W, Wang X, Zou J, Li L, Liu J, Sun Y. 2016. Late Pleistocene climate change promoted divergence between Picea asperata and P. crassifolia on the Qinghai-Tibet Plateau through recent bottlenecks. Ecology and Evolution  6: 4435– 4444. Google Scholar CrossRef Search ADS PubMed  Bouckaert R, Heled J, Kühnert D, Vaughan T, Wu CH, Xie D, Suchard MA, Rambaut A, Drummond AJ. 2014. BEAST 2: a software platform for Bayesian evolutionary analysis. PLoS Computational Biology  10: e1003537. Google Scholar CrossRef Search ADS PubMed  Butlin R, Debelle A, Kerth C, Snook RR, Beukeboom LW, Castillo CRF, Diao W, Maan ME, Paolucci S, Weissing FJ, van de Zande L, Hoikkala A, Geuverink E, Jennings J, Kankare M, Knott KE, Tyukmaeva VI, Zoumadakis C, Ritchie MG, Barker D, Immonen E, Kirkpatrick M, Noor M, Macias Garcia C, Schmitt T, Schilthuizen M. 2012. What do we need to know about speciation? Trends in Ecology & Evolution  27: 27– 39. Google Scholar CrossRef Search ADS PubMed  Chen FB. 1992. H-D event: an important tectonic event of the late Cenozoic in Eastern Asia. Mountain Research  10: 195– 202. Chen FB. 1996. Second discussion on the H-D movement. Mountain Research  17: 14– 22. Chen K, Abbott RJ, Milne RI, Tian XM, Liu J. 2008. Phylogeography of Pinus tabulaeformis Carr. (Pinaceae), a dominant species of coniferous forest in northern China. Molecular Ecology  17: 4276– 4288. Google Scholar CrossRef Search ADS PubMed  Chen SL, Phillips SM. 2006. Orinus (Poaceae). In: Wu ZY, Raven PH, ed. Flora of China, Vol. 22 . Beijing/St. Louis: Science Press/Missouri Botanical Garden, 464– 465. Choleva L, Musilova Z, Kohoutova-Sediva A, Paces J, Rab P, Janko K. 2014. Distinguishing between incomplete lineage sorting and genomic introgressions: complete fixation of allospecific mitochondrial DNA in a sexually reproducing fish (Cobitis; Teleostei), despite clonal reproduction of hybrids. PLoS One  9: e80641. Google Scholar CrossRef Search ADS PubMed  Christin PA, Besnard G, Samaritani E, Duvall MR, Hodkinson TR, Savolainen V, Salamin N. 2008. Oligocene CO2 decline promoted C4 photosynthesis in grasses. Current Biology: CB  18: 37– 43. Google Scholar CrossRef Search ADS PubMed  Comes HP, Kadereit JW. 1998. The effect of Quaternary climatic changes on plant distribution and evolution. Trends in Plant Science  3: 432– 438. Google Scholar CrossRef Search ADS   Crawford DJ, Archibald JK. 2017. Island floras as model systems for studies of plant speciation: prospects and challenges. Journal of Systematics and Evolution  55: 1– 15. Google Scholar CrossRef Search ADS   Doyle JJ, Doyle JL. 1987. A rapid DNA isolation procedure for small quantities of fresh leaf material. Phytochemical Bulletin, Botanical Society of America  19: 11– 15. Drummond AJ, Rambaut A. 2007. BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evolutionary Biology  7: 214. Google Scholar CrossRef Search ADS PubMed  Edgar RC. 2004. MUSCLE: multiple sequence alignment with high accuracy and high throughput. Nucleic Acids Research  32: 1792– 1797. Google Scholar CrossRef Search ADS PubMed  ESRI. 2010. ESRI info: Company history . Available at: http://www.esri.com/about-esri/about/history.html (accessed 6 May 2010). Evanno G, Regnaut S, Goudet J. 2005. Detecting the number of clusters of individuals using the software structure: a simulation study. Molecular Ecology  14: 2611– 2620. Google Scholar CrossRef Search ADS PubMed  Excoffier L, Laval G, Schneider S. 2005. ARLEQUIN (version 3.0): an integrated software package for population genetics data analysis. Evolutionary Bioinformatics  1: 47– 50. Google Scholar CrossRef Search ADS   Excoffier L, Smouse PE, Quattro JM. 1992. Analysis of molecular variance inferred from metric distances among DNA haplotypes: application to human mitochondrial DNA restriction data. Genetics  131: 479– 491. Google Scholar PubMed  Fan DM, Yue JP, Nie ZL, Li ZM, Comes HP, Sun H. 2013. Phylogeography of Sophora davidii (Leguminosae) across the ‘Tanaka-Kaiyong Line’, an important phytogeographic boundary in Southwest China. Molecular Ecology  22: 4270– 4288. Google Scholar CrossRef Search ADS PubMed  Fang XM, Lü LQ, Yang SL, Li JJ, An ZS, Jiang PA, Chen XL. 2002. Loess in Kunlun Mountains and its implications on desert development and Tibetan Plateau uplift in west China. Science in China Series D  45: 289– 299. Google Scholar CrossRef Search ADS   Farris JS, Källersjö M, Kluge AG, Bult C. 1995. Testing the significance of incongruence. Cladistics  10: 315– 319. Google Scholar CrossRef Search ADS   Favre A, Päckert M, Pauls SU, Jähnig SC, Uhl D, Michalak I, Muellner-Riehl AN. 2015. The role of the uplift of the Qinghai-Tibetan Plateau for the evolution of Tibetan biotas. Biological Reviews of the Cambridge Philosophical Society  90: 236– 253. Google Scholar CrossRef Search ADS PubMed  Felsenstein J. 1985. Confidence limits on phylogenies: an approach using the bootstrap. Evolution  39: 783– 791. Google Scholar CrossRef Search ADS PubMed  Feng YX, Wang XQ, Pan KY, Hong DY. 1998. A revaluation of the systematic positions of the Cercidiphyllaceae and Daphniphyllaceae based on rbcL gene sequence analysis, with reference to the relationship in the “lower” Hamamelidae. Acta Phytotaxonomica Sinica  36: 411– 422. Freeland JR, Kirk H, Petersen SD. 2011. Molecular ecology, 2nd edn . Chichester: John Wiley & Sons. Google Scholar CrossRef Search ADS   Fuertes-Aguilar J, Nieto-Feliner G. 2003. Additive polymorphisms and reticulation in an ITS phylogeny of thrifts (Armeria, Plumbaginaceae). Molecular Phylogenetics and Evolution  28: 430– 447. Google Scholar CrossRef Search ADS PubMed  Fuertes-Aguilar J, Rosselló JA, Nieto-Feliner G. 1999. Nuclear ribosomal DNA (nrDNA) concerted evolution in natural and artificial hybrids of Armeria (Plumbaginaceae). Molecular Ecology  8: 1341– 1346. Google Scholar CrossRef Search ADS PubMed  Fu YX. 1997. Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics  147: 915– 925. Google Scholar PubMed  Funk DJ, Omland KE. 2003. Species level paraphyly and polyphyly: frequency, causes, and consequences, with insights from animal mitochondrial DNA. Annual Review of Ecology, Evolution, and Systematics  34: 397– 423. Google Scholar CrossRef Search ADS   Gao YD, Harris AJ, Zhou SD, He XJ. 2013. Evolutionary events in Lilium (including Nomocharis, Liliaceae) are temporally correlated with orogenies of the Q-T plateau and the Hengduan Mountains. Molecular Phylogenetics and Evolution  68: 443– 460. Google Scholar CrossRef Search ADS PubMed  Gao QB, Zhang DJ, Duan YZ, Zhang FQ, Li YH, Fu PC, Chen SL. 2012. Intraspecific divergences of Rhodiola alsia (Crassulaceae) based on plastid DNA and internal transcribed spacer fragments. Botanical Journal of the Linnean Society  168: 204– 215. Google Scholar CrossRef Search ADS   Garnhart NJ. 2001. Binthere V1.0, a program to bin AFLP data . Durham: University of New Hampshire. Guo ZT, Ruddiman WF, Hao QZ, Wu HB, Qiao YS, Zhu RX, Peng SZ, Wei JJ, Yuan BY, Liu TS. 2002. Onset of Asian desertification by 22 Myr ago inferred from loess deposits in China. Nature  416: 159– 163. Google Scholar CrossRef Search ADS PubMed  Hamilton MB. 2009. Population genetics . Chichester: John Wiley & Sons Ltd. Harrison TM, Copeland P, Kidd WS, Yin A. 1992. Raising Tibet. Science  255: 1663– 1670. Google Scholar CrossRef Search ADS PubMed  Harvey PH, May RM, Nee S. 1994. Phylogenies without fossils. Evolution  48: 523– 529. Google Scholar CrossRef Search ADS PubMed  Hauquier F, Leliaert F, Rigaux A, Derycke S, Vanreusel A. 2017. Distinct genetic differentiation and species diversification within two marine nematodes with different habitat preference in Antarctic sediments. BMC Evolutionary Biology  17: 120. Google Scholar CrossRef Search ADS PubMed  Hewitt GM. 1996. Some genetic consequences of ice ages, and their role in divergence and speciation. Biological Journal of the Linnean Society  58: 247– 276. Google Scholar CrossRef Search ADS   Hewitt GM. 1999. Post-glacial re-colonization of European Biota. Biological Journal of the Linnean Society  68: 87– 112. Google Scholar CrossRef Search ADS   Hewitt GM. 2000. The genetic legacy of the Quaternary ice ages. Nature  405: 907– 913. Google Scholar CrossRef Search ADS PubMed  Hewitt GM. 2004. Genetic consequences of climatic oscillations in the Quaternary. Philosophical Transactions of the Royal Society of London B: Biological Sciences  359: 183– 195. Google Scholar CrossRef Search ADS   Hewitt GM. 2011. Quaternary phylogeography: the roots of hybrid zones. Genetica  139: 617– 638. Google Scholar CrossRef Search ADS PubMed  Hudson RR. 1990. Gene genealogies and the coalescent process. In: Futuyma D, Antonovics J, eds. Oxford surveys in evolutionary biology . Oxford: Oxford University Press, 1– 44. Jia DR, Abbott RJ, Liu TL, Mao KS, Bartish IV, Liu JQ. 2012. Out of the Qinghai-Tibet Plateau: evidence for the origin and dispersal of Eurasian temperate plants from a phylogeographic study of Hippophaë rhamnoides (Elaeagnaceae). The New Phytologist  194: 1123– 1133. Google Scholar CrossRef Search ADS PubMed  Jia DR, Liu TL, Wang LY, Zhou DW, Liu JQ. 2011. Evolutionary history of an alpine shrub Hippophae tibetana (Elaeagnaceae): allopatric divergence and regional expansion. Botanical Journal of the Linnean Society  102: 37– 50. Google Scholar CrossRef Search ADS   Jiang DB, Lang XM, Tian ZP, Guo DL. 2011. Last glacial maximum climate over China from PMIP simulations. Palaeogeography, Palaeoclimatology, Palaeoecology  309: 347– 357. Google Scholar CrossRef Search ADS   Kuritzin A, Kischka T, Schmitz J, Churakov G. 2016. Incomplete lineage sorting and hybridization statistics for large-scale retroposon insertion data. PLoS Computational Biology  12: e1004812. Google Scholar CrossRef Search ADS PubMed  de Lafontaine G, Amasifuen Guerra CA, Ducousso A, Sánchez-Goñi MF, Petit RJ. 2014. Beyond skepticism: uncovering cryptic refugia using multiple lines of evidence. The New Phytologist  204: 450– 454. Google Scholar CrossRef Search ADS PubMed  Li JJ, Shi YF, Li BY. 1995. Uplift of the Qinghai-Xizang (Tibet) Plateau and global change . Lanzhou: Lanzhou University Press. Li ZH, Zou JB, Mao KS, Lin K, Li HP, Liu JQ, Källman T. 2012. Population genetic evidence for complex evolutionary histories of four high altitude juniper species in the Qinghai-Tibetan Plateau. Evolution  66: 831– 845. Google Scholar CrossRef Search ADS PubMed  Librado P, Rozas J. 2009. DnaSP v5: a software for comprehensive analysis of DNA polymorphism data. Bioinformatics  25: 1451– 1452. Google Scholar CrossRef Search ADS PubMed  Liu JQ. 2004. Uniformity of karyotypes in Ligularia (Asteraceae: Senecioneae), a highly-diversified genus of the eastern Qinghai-Tibet Plateau highlands and adjacent areas. Botanical Journal of the Linnean Society  144: 329– 342. Google Scholar CrossRef Search ADS   Liu B, Abbott RJ, Lu Z, Tian B, Liu J. 2014a. Diploid hybrid origin of Ostryopsis intermedia (Betulaceae) in the Qinghai-Tibet Plateau triggered by Quaternary climate change. Molecular Ecology  23: 3013– 3027. Google Scholar CrossRef Search ADS   Liu JQ, Duan YW, Hao G, Ge XJ, Sun H. 2014b. Evolutionary history and underlying adaptation of alpine plants on the Qinghai-Tibet Plateau. Journal of Systematics and Evolution  52: 241– 249. Google Scholar CrossRef Search ADS   Liu JQ, Gao TG, Chen ZD, Lu AM. 2002. Molecular phylogeny and biogeography of the Qinghai-Tibet Plateau endemic Nannoglottis (Asteraceae). Molecular Phylogenetics and Evolution  23: 307– 325. Google Scholar CrossRef Search ADS PubMed  Liu YP, Su X, He YH, Han LM, Huang YY, Wang ZZ. 2015. Evolutionary history of Orinus thoroldii (Poaceae), endemic to the western Qinghai-Tibetan Plateau in China. Biochemical Systematics and Ecology  59: 159– 167. Google Scholar CrossRef Search ADS   Liu JQ, Sun YS, Ge XJ, Gao LM, Qiu YX. 2012. Phylogeographic studies of plants in China: advances in the past and directions in the future. Journal of Systematics and Evolution  50: 267– 275. Google Scholar CrossRef Search ADS   Liu JQ, Wang YJ, Wang AL, Hideaki O, Abbott RJ. 2006. Radiation and diversification within the Ligularia-Cremanthodium-Parasenecio complex (Asteraceae) triggered by uplift of the Qinghai-Tibetan Plateau. Molecular Phylogenetics and Evolution  38: 31– 49. Google Scholar CrossRef Search ADS PubMed  Ma YZ, Li ZH, Wang X, Shang BL, Wu GL, Wang YJ. 2014. Phylogeography of the genus Dasiphora (Rosaceae) in the Qinghai-Tibetan Plateau: divergence blurred by expansion. Botanical Journal of the Linnean Society  111: 777– 788. Google Scholar CrossRef Search ADS   Mao K, Hao G, Liu J, Adams RP, Milne RI. 2010. Diversification and biogeography of Juniperus (Cupressaceae): variable diversification rates and multiple intercontinental dispersals. The New Phytologist  188: 254– 272. Google Scholar CrossRef Search ADS PubMed  Mayr E. 1942. Systematics and the origin of species, from the viewpoint of a zoologist . New York: Columbia University Press. Meng L, Chen G, Li Z, Yang Y, Wang Z, Wang L. 2015. Refugial isolation and range expansions drive the genetic structure of Oxyria sinensis (Polygonaceae) in the Himalaya-Hengduan Mountains. Scientific Reports  5: 10396. Google Scholar CrossRef Search ADS PubMed  Meng HH, Gao XY, Huang JF, Zhang ML. 2014. Plant phylogeography in arid northwest China: retrospectives and perspectives. Journal of Systematics and Evolution  52: 1– 16. Google Scholar CrossRef Search ADS   Meng L, Yang R, Abbott RJ, Miehe G, Hu T, Liu J. 2007. Mitochondrial and chloroplast phylogeography of Picea crassifolia Kom. (Pinaceae) in the Qinghai-Tibetan Plateau and adjacent highlands. Molecular Ecology  16: 4128– 4137. Google Scholar CrossRef Search ADS PubMed  Miraldo A, Hewitt GM, Paulo OS, Emerson BC. 2011. Phylogeography and demographic history of Lacerta lepida in the Iberian Peninsula: multiple refugia, range expansions and secondary contact zones. BMC Evolutionary Biology  11: 170. Google Scholar CrossRef Search ADS PubMed  Mittermeier RA, Gil PR, Hoffman M, Pilgrim J, Brooks T, Mittermeier CG, Lamoreux J, da Fonseca GAB. 2005. Hotspots revisited: earth’s biologically richest and most endangered terrestrial ecoregions . Chicago: University of Chicago Press. Nason JD, Hamrick JL, Fleming TH. 2002. Historical vicariance and postglacial colonization effects on the evolution of genetic structure in Lophocereus, a Sonoran Desert columnar cactus. Evolution  56: 2214– 2226. Google Scholar CrossRef Search ADS PubMed  Nei M, Li WH. 1979. Mathematical model for studying genetical variation in terms of restriction endonucleases. Proceedings of the National Academy of Sciences of the United States of America  76: 5269– 5273. Google Scholar CrossRef Search ADS PubMed  Nosil P. 2012. Ecological speciation . New York: Oxford University Press. Google Scholar CrossRef Search ADS   Opgenoorth L, Vendramin GG, Mao K, Miehe G, Miehe S, Liepelt S, Liu J, Ziegenhagen B. 2010. Tree endurance on the Tibetan Plateau marks the world’s highest known tree line of the Last Glacial Maximum. The New Phytologist  185: 332– 342. Google Scholar CrossRef Search ADS PubMed  Osborne CP, Beerling DJ. 2006. Nature’s green revolution: the remarkable evolutionary rise of C4 plants. Philosophical Transactions of the Royal Society of London B: Biological Sciences  361: 173– 194. Google Scholar CrossRef Search ADS   Parchman TL, Benkman CW, Britch SC. 2006. Patterns of genetic variation in the adaptive radiation of New World crossbills (Aves: Loxia). Molecular Ecology  15: 1873– 1887. Google Scholar CrossRef Search ADS PubMed  Peterson PM, Romaschenko K, Arrieta YH. 2016. A molecular phylogeny and classification of the Cynodonteae (Poaceae: Chloridoideae) with four new genera: Orthacanthus, Triplasiella, Tripogonella, and Zaqiqah; three new subtribes: Dactylocteniinae, Orininae, and Zaqiqahinae; and a subgeneric classification of Distichlis. Taxon  65: 1263– 1287. Google Scholar CrossRef Search ADS   Peterson PM, Romaschenko K, Snow N, Johnson G. 2012. A molecular phylogeny and classification of Leptochloa (Poaceae: Chloridoideae: Chlorideae) sensu lato and related genera. Annals of Botany  109: 1317– 1330. Google Scholar CrossRef Search ADS PubMed  Petit R, Aguinagalde I, de Beaulieu JL, Bittkau C, Brewer S, Cheddadi R, Ennos R, Fineschi S, Grivet D, Lascoux M, Mohanty A, Müller-Starck G, Demesure-Musch B, Palmé A, Martín JP, Rendell S, Vendramin GG. 2003. Glacial refugia: hotspots but not melting pots of genetic diversity. Science  300: 1563– 1565. Google Scholar CrossRef Search ADS PubMed  Petit RJ, Excoffier L. 2009. Gene flow and species delimitation. Trends in Ecology & Evolution  24: 386– 393. Google Scholar CrossRef Search ADS PubMed  Pons O, Chaouche K. 1995. Estimation, variance and optimal sampling of gene diversity II. Diploid locus. Theoretical and Applied Genetics  91: 122– 130. Google Scholar CrossRef Search ADS PubMed  Pons O, Petit RJ. 1996. Measuring and testing genetic differentiation with ordered versus unordered alleles. Genetics  144: 1237– 1245. Google Scholar PubMed  Posada D. 1998. jModelTest: phylogenetic model averaging. Molecular Biology and Evolution  25: 1235– 1256. Pritchard JK, Stephens M, Donnelly P. 2000. Inference of population structure using multilocus genotype data. Genetics  155: 945– 959. Google Scholar PubMed  Qiu YX, Fu CX, Comes HP. 2011. Plant molecular phylogeography in China and adjacent regions: tracing the genetic imprints of Quaternary climate and environmental change in the world’s most diverse temperate flora. Molecular Phylogenetics and Evolution  59: 225– 244. Google Scholar CrossRef Search ADS PubMed  Qiu YX, Liu YF, Kang M, Yi GM, Huang HW. 2013. Spatial and temporal population genetic variation and structure of Nothotsuga longibracteata (Pinaceae), a relic conifer species endemic to subtropical China. Genetics and Molecular Biology  36: 598– 607. Google Scholar CrossRef Search ADS PubMed  Rambaut A. 2009. FIGTREE 1.3.1. Available at: http://tree.bio.ed.ac.uk/software/figtree (accessed 7 March 2014). Rambaut A, Drummond AJ. 2007. Tracer. Available at: http://beast.bio.ed.ac.uk/tracer/ (accessed 11 December 2013). Rambaut A, Drummond AJ. 2009. TRACER v1.5. Available at: http://tree.bio.ed.ac.uk/software/tracer/ (accessed 7 March 2014). Rieseberg LH, Church SA, Morjan CL. 2004. Integration of populations and differentiation of species. The New Phytologist  161: 59– 69. Google Scholar CrossRef Search ADS PubMed  Rogers AR, Harpending H. 1992. Population growth makes waves in the distribution of pairwise genetic differences. Molecular Biology and Evolution  9: 552– 569. Google Scholar PubMed  Rohlf FJ. 2000. NTSYS-pc version 2.1: numerical taxonomy and multivariance analysis system . Setauket: Applied Biostatistics Inc., Exeter Software. Ronquist F, Huelsenbeck JP. 2003. MrBayes 3: Bayesian phylogenetic inference under mixed models. Bioinformatics  19: 1572– 1574. Google Scholar CrossRef Search ADS PubMed  Ruddiman W. 1998. Geology: early uplift in Tibet? Nature  394: 723– 725. Google Scholar CrossRef Search ADS   Sakaguchi S, Takeuchi Y, Yamasaki M, Sakurai S, Isagi Y. 2011. Lineage admixture during postglacial range expansion is responsible for the increased gene diversity of Kalopanax septemlobus in a recently colonised territory. Heredity  107: 338– 348. Google Scholar CrossRef Search ADS PubMed  Sang T, Crawford D, Stuessy T. 1997. Chloroplast DNA phylogeny, reticulate evolution, and biogeography of Paeonia (Paeoniaceae). American Journal of Botany  84: 1120. Google Scholar CrossRef Search ADS PubMed  Sardell JM, Uy JAC. 2016. Hybridization following recent secondary contact results in asymmetric genotypic and phenotypic introgression between island species of Myzomela honeyeaters. Evolution  70: 257– 269. Google Scholar CrossRef Search ADS PubMed  Schluter D. 2009. Evidence for ecological speciation and its alternative. Science  323: 737– 741. Google Scholar CrossRef Search ADS PubMed  Segovia RA, Pérez MF, Hinojosa LF. 2012. Genetic evidence for glacial refugia of the temperate tree Eucryphia cordifolia (Cunoniaceae) in southern South America. American Journal of Botany  99: 121– 129. Google Scholar CrossRef Search ADS PubMed  Shackleton NJ, Opdyke ND. 1973. Oxygen isotope and palaeomagnetic stratigraphy of Equatorial Pacific core V28-238: oxygen isotope temperatures and ice volumes on a 105 year and 106 year scale. Quaternary Research  3: 39– 55. Google Scholar CrossRef Search ADS   Shaw J, Shafer HL, Leonard OR, Margaret JK, Schorr M, Morris AB. 2017. Chloroplast DNA sequence utility for the lowest phylogenetic and phylogeographic inferences in angiosperms: the tortoise and the hare IV. American Journal of Botany  101: 1987– 2004. Google Scholar CrossRef Search ADS   Shi YF, Li JJ, Li BY. 1998. Uplift and environmental changes of Qinghai-Tibetan Plateau in the late Cenozoic . Guangzhou: Guangdong Science and Technology Press. Slatkin M, Hudson RR. 1991. Pairwise comparisons of mitochondrial DNA sequences in stable and exponentially growing populations. Genetics  129: 555– 562. Google Scholar PubMed  Soreng RJ, Peterson PM, Romaschenko K, Davidse G, Teisher JK, Clark LG, Barbéra P, Gillespie LJ, Zuloaga FO. 2017. A worldwide phylogenetic classification of the Poaceae (Gramineae) II: an update and comparison of two 2015 classifications. Journal of Systematics and Evolution  55: 259– 290. Google Scholar CrossRef Search ADS   Sousa V, Hey J. 2013. Understanding the origin of species with genome-scale data: modelling gene flow. Nature Reviews Genetics  14: 404– 414. Google Scholar CrossRef Search ADS PubMed  Stephens M, Smith NJ, Donnelly P. 2011. A new statistical method for haplotype reconstruction from population data. American Journal of Human Genetics  68: 978– 989. Google Scholar CrossRef Search ADS   Stewart JR, Lister AM, Barnes I, Dalen L. 2010. Refugia revisited: individualistic responses of species in space and time. Proceedings of the Royal Society of London B: Biological Sciences  277: 661– 671. Google Scholar CrossRef Search ADS   Su X, Cai LB. 2009. Orinus longiglumis (Poaceae: Chloridoideae), a new species from Xizang (Tibet), China. Annales Botanici Fennici  46: 143– 147. Google Scholar CrossRef Search ADS   Su X, Liu YP, Wu GL, Luo WC, Liu JQ. 2017. A taxonomic revision of Orinus (Poaceae) with a new species, O. intermedius, from the Qinghai-Tibet Plateau. Novon  25: 206– 213. Google Scholar CrossRef Search ADS   Su X, Wu G, Li L, Liu J. 2015. Species delimitation in plants using the Qinghai-Tibet Plateau endemic Orinus (Poaceae: Tridentinae) as an example. Annals of Botany  116: 35– 48. Google Scholar CrossRef Search ADS PubMed  Sun Y, Abbott RJ, Li L, Li L, Zou J, Liu J. 2014. Evolutionary history of Purple cone spruce (Picea purpurea) in the Qinghai-Tibet Plateau: homoploid hybrid origin and Pleistocene expansion. Molecular Ecology  23: 343– 359. Google Scholar CrossRef Search ADS PubMed  Sun YB, An ZS. 2001. History and variability of Asian interior aridity recorded by eolian flux in the Chinese Loess Plateau during the past 7 Ma. Science in China Series D: Earth Sciences  45: 420– 429. Google Scholar CrossRef Search ADS   Sun H, McLewin W, Fay MF. 2001. Molecular phylogeny of Helleborus (Ranunculaceae), with an emphasis on the East Asia-Mediterranean disjunction. Taxon  50: 1001– 1018. Google Scholar CrossRef Search ADS   Swofford DL. 2002. PAUP*: Phylogenetic analysis using parsimony (*and other methods), version 4 . Sunderland: Sinauer Associates. Tajima F. 1989. Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics  123: 585– 595. Google Scholar PubMed  Tamura K, Peterson D, Peterson N, Stecher G, Nei M, Kumar S. 2011. MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods. Molecular Biology and Evolution  28: 2731– 2739. Google Scholar CrossRef Search ADS PubMed  Tang LY, Shen CM. 1996. Late Cenozoic vegetational history and climatic characteristics of Qinghai-Xizang Plateau. Acta Micropalaeontologica Sinica  13: 321– 337. Tian X, Luo J, Wang A, Mao K, Liu J. 2011. On the origin of the woody buckwheat Fagopyrum tibeticum (=Parapteropyrum tibeticum) in the Qinghai-Tibetan Plateau. Molecular Phylogenetics and Evolution  61: 515– 520. Google Scholar CrossRef Search ADS PubMed  Tzedakis PC, Emerson BC, Hewitt GM. 2013. Cryptic or mystic? Glacial tree refugia in northern Europe. Trends in Ecology & Evolution  28: 696– 704. Google Scholar CrossRef Search ADS PubMed  Tzedakis PC, Lawson IT, Forgley MR, Hewitt GM. 2002. Buffered tree population changes in Quaternary refugium: evolutionary implications. Science  292: 267– 269. Vicentini A, Barber JC, Aliscioni SS, Giussani LM, Kellogg EA. 2008. The age of the grasses and clusters of origins of C4 photosynthesis. Global Change Biology  14: 2963– 2977. Google Scholar CrossRef Search ADS   de Villiers MJ, Pirie MD, Hughes M, Möller M, Edwards TJ, Bellstedt DU. 2013. An approach to identify putative hybrids in the ‘coalescent stochasticity zone’, as exemplified in the African plant genus Streptocarpus (Gesneriaceae). The New Phytologist  198: 284– 300. Google Scholar CrossRef Search ADS PubMed  Vos P, Hogers R, Bleeker M, Reijans M, van de Lee T, Hornes M, Frijters A, Pot J, Peleman J, Kuiper M. 1995. AFLP: a new technique for DNA fingerprinting. Nucleic Acids Research  23: 4407– 4414. Google Scholar CrossRef Search ADS PubMed  Voss N, Eckstein RL, Durka W. 2012. Range expansion of a selfing polyploid plant despite widespread genetic uniformity. Annals of Botany  110: 585– 593. Google Scholar CrossRef Search ADS PubMed  Wan DS, Feng JJ, Jiang DC, Mao KS, Duan YW, Miehe G, Opgenoorth L. 2016. The Quaternary evolutionary history, potential distribution dynamics, and conservation implications for a Qinghai-Tibet Plateau endemic herbaceous perennial, Anisodus tanguticus (Solanaceae). Ecology and Evolution  6: 1977– 1995. Google Scholar CrossRef Search ADS PubMed  Wang Q, Abbott RJ, Yu QS, Lin K, Liu JQ. 2013. Pleistocene climate change and the origin of two desert plant species, Pugionium cornutum and Pugionium dolabratum (Brassicaceae), in northwest China. The New Phytologist  199: 277– 287. Google Scholar CrossRef Search ADS PubMed  Wang L, Abbott RJ, Zheng W, Chen P, Wang Y, Liu J. 2009a. History and evolution of alpine plants endemic to the Qinghai-Tibetan Plateau: Aconitum gymnandrum (Ranunculaceae). Molecular Ecology  18: 709– 721. Google Scholar CrossRef Search ADS   Wang GN, He XY, Miehe G, Mao KS. 2014. Phylogeography of the Qinghai-Tibet Plateau endemic alpine herb Pomatosace filicula (Primulaceae). Journal of Systematics and Evolution  52: 1– 14. Google Scholar CrossRef Search ADS   Wang LY, Ikeda H, Liu TL, Wang YJ, Liu JQ. 2009b. Repeated range expansion and glacial endurance of Potentilla glabra (Rosaceae) in the Qinghai-Tibetan Plateau. Journal of Integrative Plant Biology  51: 698– 706. Google Scholar CrossRef Search ADS   Wang YJ, Susanna A, Von Raab-Straube E, Milne R, Liu JQ. 2009c. Island-like radiation of Saussurea (Asteraceae: Cardueae) triggered by uplifts of the Qinghai-Tibetan Plateau. Botanical Journal of the Linnean Society  97: 893– 903. Google Scholar CrossRef Search ADS   Weir BS. 1996. Genetic data analysis II . Sunderland: Sinauer Associates. Wen J, Nie ZL, Ickert-Bond SM. 2016. Intercontinental disjunctions between eastern Asia and western North America in vascular plants highlight the biogeographic importance of the Bering land bridge from late Cretaceous to Neogene. Journal of Systematics and Evolution  54: 469– 490. Google Scholar CrossRef Search ADS   Wen J, Zhang JQ, Nie ZL, Zhong Y, Sun H. 2014. Evolutionary diversifications of plants on the Qinghai-Tibetan Plateau. Frontiers in Genetics  5: 4. Google Scholar PubMed  White TJ, Bruns T, Lee S, Taylor J. 1990. Amplification and direct sequencing of fungal ribosomal RNA genes for phylogenetics. In: Innis M, Gelfand D, Sninsky J, White T, eds. PCR protocols . San Diego: Academic Press, 315– 322. Google Scholar CrossRef Search ADS   Wolfe KH, Li WH, Sharp PM. 1987. Rates of nucleotide substitution vary greatly among plant mitochondrial, chloroplast, and nuclear DNAs. Proceedings of the National Academy of Sciences of the United States of America  84: 9054– 9058. Google Scholar CrossRef Search ADS PubMed  Wu CY. 1988. Hengduan mountains flora and their significance. Journal of Japanese Botany  63: 297– 311. Wu ZH, Barosh PJ, Wu ZH, Hu DG, Zhao X, Ye PS. 2008. Vast early Miocene lakes of the central Tibetan Plateau. Geological Society of America Bulletin  120: 1326– 1337. Google Scholar CrossRef Search ADS   Wu LL, Cui XK, Milne RI, Sun YS, Liu JQ. 2010. Multiple autopolyploidizations and range expansion of Allium przewalskianum Regel. (Alliaceae) in the Qinghai-Tibetan Plateau. Molecular Ecology  19: 1691– 1704. Google Scholar CrossRef Search ADS PubMed  Xing YW, Ree RH. 2017. Uplift-driven diversification in the Hengduan Mountains, a temperate biodiversity hotspot. Proceedings of the National Academy of Sciences of the United States of America  114: E3444– E3451. Google Scholar CrossRef Search ADS PubMed  Xu T, Abbott RJ, Milne RI, Mao K, Du FK, Wu G, Ciren Z, Miehe G, Liu J. 2010. Phylogeography and allopatric divergence of cypress species (Cupressus L.) in the Qinghai-Tibetan Plateau and adjacent regions. BMC Evolutionary Biology  10: 194. Google Scholar CrossRef Search ADS PubMed  Yang FS, Li YF, Ding X, Wang XQ. 2008. Extensive population expansion of Pedicularis longiflora (Orobanchaceae) on the Qinghai-Tibetan Plateau and its correlation with the Quaternary climate change. Molecular Ecology  17: 5135– 5145. Google Scholar CrossRef Search ADS PubMed  Yang FS, Qin AL, Li YF, Wang XQ. 2012. Great genetic differentiation among populations of Meconopsis integrifolia and its implication for plant speciation in the Qinghai-Tibetan Plateau. PLoS One  7: e37196. Google Scholar CrossRef Search ADS PubMed  Yap IV, Nelson RJ. 1996. Winboot: a program for performing bootstrap analysis of binary data to determine the confidence limits of UPGMA-based dendrograms . IRRI Discussion Paper Series No. 14. Manila: International Rice Research Institute. Yasuda N, Taquet C, Nagai S, Fortes M, Fan TY, Harii S, Yoshida T, Sito Y, Nadaoka K. 2015. Genetic diversity, paraphyly and incomplete lineage sorting of mtDNA, ITS2 and microsatellite flanking region in closely related Heliopora species (Octocorallia). Molecular Phylogenetics and Evolution  93: 161– 171. Google Scholar CrossRef Search ADS PubMed  Yu QS, Wang Q, Wu GL, Ma YZ, He XY, Wang X, Xie PH, Hu LH, Liu JQ. 2013. Genetic differentiation and delimitation of Pugionium dolabratum and Pugionium cornutum (Brassicaceae). Plant Systematics and Evolution  299: 1355– 1365. Google Scholar CrossRef Search ADS   Zhang Q, Chiang TY, George M, Liu JQ, Abbott RJ. 2005. Phylogeography of the Qinghai-Tibetan Plateau endemic Juniperus przewalskii (Cupressaceae) inferred from chloroplast DNA sequence variation. Molecular Ecology  14: 3513– 3524. Google Scholar CrossRef Search ADS PubMed  Zhang JQ, Meng SY, Wen J, Rao GY. 2015. DNA barcoding of Rhodiola (Crassulaceae): a case study on a group of recently diversified medicinal plants from the Qinghai-Tibetan Plateau. PLoS One  10: 1– 15. Zhao JL, Gugger PF, Xia YM, Li QJ. 2016. Ecological divergence of two closely related Roscoea species associated with late Quaternary climate change. Journal of Biogeography  43: 1990– 2001. Google Scholar CrossRef Search ADS   Zheng D. 1996. The system of physico-geographical regions of the Qinghai-Xizang (Tibet) Plateau. Science China: Earth Sciences  39: 410– 417. Zheng HB, Powell CM, Rea DK, Wang JL, Wang PX. 2004. Late Miocene and mid-Pliocene enhancement of the East Asian monsoon as viewed from the land and sea. Global and Planetary Change  41: 147– 155. Google Scholar CrossRef Search ADS   Zheng BX, Xu QQ, Shen YP. 2002. The relationship between climate change and Quaternary glacial cycles on the Qinghai-Tibetan Plateau: review and speculation. Quaternary International  97–98: 93– 101. Google Scholar CrossRef Search ADS   Zhou Y, Duvaux L, Ren G, Zhang L, Savolainen O, Liu J. 2017. Importance of incomplete lineage sorting and introgression in the origin of shared genetic variation between two closely related pines with overlapping distributions. Heredity  118: 211– 220. Google Scholar CrossRef Search ADS PubMed  Zhou SL, Qian P. 2003. Matrix Generator (MG): a program for creating 0/1 matrix from DNA fragments. Acta Botanica Sinica  45: 766. Zimmer EA, Wen J. 2012. Using nuclear gene data for plant phylogenetics: progress and prospects. Molecular Phylogenetics and Evolution  65: 774– 785. Google Scholar CrossRef Search ADS PubMed  Zimmer EA, Wen J. 2015. Using nuclear gene data for plant phylogenetics: progress and prospects II. Next-gen approaches. Journal of Systematics and Evolution  53: 371– 379. Google Scholar CrossRef Search ADS   Zou JB, Peng XL, Li L, Liu JQ, Miehe G, Opgenoorth L. 2012. Molecular phylogeography and evolutionary history of Picea likiangensis in the Qinghai-Tibetan Plateau inferred from mitochondrial and chloroplast DNA sequence variation. Journal of Systematics and Evolution  50: 341– 350. Google Scholar CrossRef Search ADS   Zuo YJ, Wen J, Ma JS, Zhou SL. 2015. Evolutionary radiation of the Panax bipinnatifidus species complex (Araliaceae) in the Sino-Himalayan region of eastern Asia as inferred from AFLP analysis. Journal of Systematics and Evolution  53: 210– 220. Google Scholar CrossRef Search ADS   © 2018 The Linnean Society of London, Botanical Journal of the Linnean Society

Journal

Botanical Journal of the Linnean SocietyOxford University Press

Published: Feb 1, 2018

There are no references for this article.

You’re reading a free preview. Subscribe to read the entire article.


DeepDyve is your
personal research library

It’s your single place to instantly
discover and read the research
that matters to you.

Enjoy affordable access to
over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month

Explore the DeepDyve Library

Search

Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly

Organize

Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place.

Access

Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals.

Your journals are on DeepDyve

Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more.

All the latest content is available, no embargo periods.

See the journals in your area

DeepDyve

Freelancer

DeepDyve

Pro

Price

FREE

$49/month
$360/year

Save searches from
Google Scholar,
PubMed

Create lists to
organize your research

Export lists, citations

Read DeepDyve articles

Abstract access only

Unlimited access to over
18 million full-text articles

Print

20 pages / month

PDF Discount

20% off