BioMed Research International
Volume 2016 (2016), Article ID 5792708, 13 pages
http://dx.doi.org/10.1155/2016/5792708

Research Article

Phylogeographic Structure of a Tethyan Relict Capparis spinosa (Capparaceae) Traces Pleistocene Geologic and Climatic Changes in the Western Himalayas, Tianshan Mountains, and Adjacent Desert Regions

Qian Wang,1,2 Ming-Li Zhang,1,3 and Lin-Ke Yin¹

¹Key Laboratory of Biogeography and Bioresource in Arid Land, Xinjiang Institute of Ecology and Geography, Chinese Academy of Sciences, Urumqi, Xinjiang 830011, China
²University of Chinese Academy of Sciences, Beijing 100049, China
³Institute of Botany, Chinese Academy of Sciences, Beijing 100093, China

Received 27 February 2016; Accepted 27 April 2016

Academic Editor: William H. Piel

Copyright © 2016 Qian Wang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Complex geological movements more or less affected or changed floristic structures, while the alternation of glacials and interglacials is presumed to have further shaped the present discontinuous genetic pattern of temperate plants. Here we consider Capparis spinosa, a xeromorphic Tethyan relict, to discuss its divergence pattern and explore how it responded in a stepwise fashion to Pleistocene geologic and climatic changes. 267 individuals from 31 populations were sampled and 24 haplotypes were identified, based on three cpDNA fragments (trnL-trnF, rps12-rpl20, and ndhF). SAMOVA clustered the 31 populations into 5 major clades. AMOVA suggests that gene flow between them might be restricted by vicariance. Molecular clock dating indicates that intraspecific divergence began in early Pleistocene, consistent with a time of intense uplift of the Himalaya and Tianshan Mountains, and intensified in mid-Pleistocene. Species distribution modeling suggests range reduction in the high mountains during the Last Glacial Maximum (LGM) as a result of cold climates when glacier advanced, while gorges at midelevations in Tianshan appear to have served as refugia. Populations of low-altitude desert regions, on the other hand, probably experienced only marginal impacts from glaciation, according to the high levels of genetic diversity.

1. Introduction

An essential purpose of biogeography is to comprehend the mechanisms resulting in the current distribution patterns of organisms during their evolutionary histories [1, 2]. Events since the late Paleogene, especially the Quaternary geologic and climatic changes, likely had a far-reaching influence on the spatial geographic structure and genetic diversification of temperate species [3, 4]. After the Oligocene, the inland Tethys or Paratethys Sea retreated, accompanied by an increasingly droughty climate, resulting in elaboration of the xerophytic Tethyan Flora [5, 6]. From the Neogene to Quaternary, with the increasing aridification and uplift of the Himalayas, Tibet Plateau, and Tianshan Mountains [7, 8], the originally continuous arid Tethyan Flora gradually evolved into the North Temperate, Mediterranean, West Asian-Central Asian, Central Asian, and Himalayan elements of today [9, 10]. Complex geological movements more or less affected or changed floristic structures, while rapid climatic shifts likewise impacted plant life worldwide during this time. Significantly, the alternation of glacials and interglacials during the Pleistocene is presumed to have played a vital role in further shaping the present discontinuous genetic pattern for a majority of temperate plant species in the Northern Hemisphere [6, 11, 12].

Phylogeography, a discipline which began more than 20 years ago, has mostly considered the response of geographic distribution and genetic structure of current species to historical events and the cycles of glaciation [13–16]. Chloroplast DNA (cpDNA), suitable for reconstructing phylogenetic patterns of plant species, is effective in exploring and estimating phylogeographical history of temperate plants caused by climatic change [17–20]. Recent scholars have concentrated on the refuge during the Last Glacial Maximum (LGM) and range expansion during inter- or postglacial periods [21–24]. With the most focus on mountain plants in Europe, North America, and East Asia [25–27], phylogeographic patterns and demographical history of eremophytes in arid Northwestern China are also receiving increasing attention [28–30].

We use Capparis spinosa, a typical xerophyte, to verify the hypothesis of Pleistocene historical influence through phylogeographic analysis and discuss how its spatial genetic structure responded to paleogeologic and paleoclimatic changes. The perennial creeping subshrub is a type member of a relatively large genus of the Capparaceae, with a natural distribution in the Mediterranean, West, and Central Asia [31, 32]. In arid Western China, there is only one single genus (Capparis Tourn. ex L.) with one single species (C. spinosa L.) belonging exclusively to this family [33]. C. spinosa ranges over the Western Himalayas, Eastern Pamir Plateau, Tianshan Mountains, Turpan-Hami Basin, and the Hexi Corridor. Previous studies on Capparis have paid more attention to autecology, phenology, quantitative morphology, and reproductive ecology [34–37]. Molecular studies related to the genetic variation of the species are little explored [38, 39] and have not considered genetic divergence from a phylogeographic point of view.

Here we employed three maternally inherited cpDNA spacers, trnL-trnF, rps12-rpl20, and ndhF, to trace historical divergence patterns of C. spinosa. We design to address the following concrete issues. () According to phylogenetic relationships yielded by cpDNA sequence variation, what is the genetic pattern of the species? () When and how have allopatric divergence and genetic differentiation of this xerophyte been influenced by glacial-interglacial cycles associated with geologic and climatic changes?

2. Materials and Methods

2.1. Population Sampling

A total of 267 individuals of C. spinosa from 31 natural populations were investigated in the present study (Table 1), throughout almost the entire geographic distribution of the species in Eastern Central Asia, including Xinjiang, western Gansu, and western Tibet. For each population, leaf samples were collected from 5 to 12 mature and large perennial individuals separated by at least 100 m, in order to avoid inbreeding and clones. The latitude, longitude, and altitude at each collection center were measured using a global positioning system (GPS) receiver. About 10 g of fresh leaf materials was gathered per individual, dried in silica gel in the field, and stored at −20°C in the laboratory. Capparis bodinieri Lévl. was chosen as an outgroup for the phylogenic study, based on previous taxonomic studies on Capparis [31]. Voucher specimens for all individuals were deposited in the Herbarium of the Xinjiang Institute of Ecology and Geography, Chinese Academy of Sciences (XJBI).

Table 1: Details of geographical locations, sample sizes, cpDNA haplotypes, and diversity indices for 31 populations of C. spinosa.

Code
Population locality
Latitude/longitude
Altitude (m)
No.
Haplotype
(±SD)
(±SD)

Overall

267

0.8925 ± 0.0075
0.0037 ± 0.0019

(I) The Western Himalayas group

    (1) ZD
Zanda, Tibet
31°30′N/79°38′E
3,598
9
H22, H23, H24
0.4167 ± 0.1907
0.0021 ± 0.0013

    (2) QS
Qusum, Tibet
31°29′N/79°48′E
3,645
5
H22
0
0

(II) The Eastern Pamir group

    (3) YGS
Yengisar, Xinjiang
38°43′N/76°14′E
2,010
9
H6
0
0

    (4) TXG
Taxkorgan, Xinjiang
37°47′N/75°13′E
3,694
5
H6
0
0

    (5) AKT
Akto, Xinjiang
38°59′N/75°32′E
2,441
9
H6
0
0

    (6) WQ
Wuqia, Xinjiang
39°44′N/75°40′E
1,945
10
H6, H21
0.4667 ± 0.1318
0.0006 ± 0.0005

    (7) ATX
Artux, Xinjiang
39°52′N/76°47′E
1,303
10
H6, H7
0.4667 ± 0.1318
0.0006 ± 0.0005

    (8) JS
Jiashi, Xinjiang
39°56′N/78°01′E
1,162
9
H6
0
0

(III) The Ili Valley group

    (9) GL
Gongliu, Xinjiang
43°37′N/81°49′E
710
9
H7
0
0

    (10) YN
Yining, Xinjiang
43°37′N/82°08′E
750
9
H7
0
0

    (11) BL
Bole, Xinjiang
44°42′N/82°05′E
470
9
H7
0
0

    (12) WS
Wusu, Xinjiang
44°19′N/84°19′E
685
9
H7, H12
0.5000 ± 0.1283
0.0013 ± 0.0009

(IV) The north side of the Tianshan Mountains group

    (13) KRMY
Karamay, Xinjiang
44°56′N/84°46′E
296
5
H12
0
0

    (14) TL
Toli, Xinjiang
45°55′N/83°36′E
684
5
H17
0
0

    (15) SHZ
Shihezi, Xinjiang
44°16′N/85°57′E
356
5
H11
0
0

    (16) SW
Shawan, Xinjiang
44°19′N/85°35′E
584
9
H5
0
0

    (17) MNS
Manas, Xinjiang
44°11′N/86°08′E
642
8
H5, H10
0.4286 ± 0.1687
0.0002 ± 0.0002

    (18) URMQ
Urumqi, Xinjiang
43°44′N/86°57′E
1,060
8
H10, H11
0.4286 ± 0.1687
0.0010 ± 0.0007

    (19) JMS
Jimsar, Xinjiang
43°55′N/89°07′E
945
8
H5
0
0

(V) The south side of the Tianshan Mountains group

    (20) KP
Kalpin, Xinjiang
40°35′N/79°36′E
1,131
10
H1, H5
0.3556 ± 0.1591
0.0008 ± 0.0006

    (21) WNS
Wensu, Xinjiang
41°29′N/79°54′E
1,296
9
H1, H8, H20
0.7500 ± 0.0786
0.0016 ± 0.0010

    (22) BC
Baicheng, Xinjiang
41°59′N/83°03′E
1,369
8
H1, H8
0.5714 ± 0.0945
0.0003 ± 0.0003

    (23) LT
Luntai, Xinjiang
41°14′N/84°12′E
921
10
H8, H14
0.3556 ± 0.1591
0.0011 ± 0.0007

    (24) KRL
Korla, Xinjiang
41°50′N/86°03′E
1,072
8
H1, H8, H13
0.6071 ± 0.1640
0.0004 ± 0.0004

    (25) HJ
Hejing, Xinjiang
42°31′N/86°16′E
1,403
10
H5, H8
0.4667 ± 0.1318
0.0012 ± 0.0008

    (26) TKS
Toksun, Xinjiang
42°49′N/88°38′E
95
12
H1, H5
0.5455 ± 0.0615
0.0012 ± 0.0008

    (27) TRP
Turpan, Xinjiang
42°52′N/89°11′E
−84
12
H1, H5, H11, H18, H19
0.8485 ± 0.0586
0.0017 ± 0.0010

    (28) SS
Shanshan, Xinjiang
43°11′N/90°09′E
1,043
9
H10, H15, H16
0.7222 ± 0.0967
0.0015 ± 0.0009

    (29) HM
Hami, Xinjiang
42°31′N/94°09′E
773
9
H9, H10, H11
0.7500 ± 0.0786
0.0013 ± 0.0009

    (30) DH
Dunhuang, Gansu
40°25′N/94°44′E
1,051
10
H1, H2, H3
0.3778 ± 0.1813
0.0005 ± 0.0004

    (31) GZ
Guazhou, Gansu
40°22′N/95°34′E
1,133
10
H1, H2, H4, H5
0.7333 ± 0.1005
0.0011 ± 0.0007

No.: number of sample individuals; : haplotype diversity; : nucleotide diversity.

2.2. DNA Extraction, Amplification, and Sequencing

Silica-gel-dried leaf tissues were ground into powder in liquid nitrogen with a Fastprep-24 device (Sample Preparation System, MP Biomedicals, USA). Total genomic DNA was extracted from about 50 mg of powdered tissues using a modified 2x CTAB method [40]. Through preliminary screening from eighteen pairs of primers (see Table in Supplementary Material available online at http://dx.doi.org/10.1155/2016/5792708), three cpDNA regions, trnL-trnF, rps12-rpl20, and ndhF (329F, 927R), were found because they have more sequence variation in individuals among populations. Polymerase chain reaction (PCR) amplification and DNA sequencing were performed using primer pairs for these three regions. PCRs were carried out in 30 μL reaction volumes including 1.5 μL of 10x buffer, 2 μL of 25 mmol/L MgCl₂, 2.2 μL of 50 ng/μL each primer (Sangon, Shanghai, China), 3 μL of 2.5 mmol/L dNTP, 0.3 μL of 5 U/μL Taq DNA polymerase, and 1.0 μL of 5 ng/μL genomic DNA. PCR amplification was run on a Gene-Amp PCR system 9700 DNA Thermal Cycler (Applied Biosystems, Foster City, CA, USA) with parameters set as follows: an initial start at 94°C with 4 min; followed by 30 cycles of 94°C with 30 s, then 59°C, 52°C, and 52°C, respectively, for trnL-trnF, rps12-rpl20, and ndhF with 30 s and 72°C with 90 s; plus a subsequent additional extension at 72°C with 10 min. Sizes of PCR products were visualized on 1% TAE-agarose gel electrophoresis. All PCR products were purified with PCR product purification kit (Shanghai SBS, Biotech Ltd., China) and then sequenced in both directions for trnL-trnF, rps12-rpl20, and ndhF using an ABI Prism 3730 xl automated sequencer in Sangon Biological Engineering Technology & Service Co., Ltd., Shanghai, China.

2.3. Data Analysis

The chloroplast DNA sequences were edited by SeqMan (Lasergene, DNASTAR Inc., Madison, Wisconsin, USA). All sequences were aligned with Clustal X version 1.81 [41], then refined, and adjusted. All resulting cpDNA sequences identified by nucleotide variations, signifying different haplotypes, were deposited in GenBank. The corresponding accession numbers are KU940216–KU940218 for trnL-trnF, KU866560–KU866562 and KU940220 for rps12-rpl20, and KU866548–KU866558 and KU940221–KU940223 for ndhF (329F, 927R); and accessions numbers of C. bodinieri as the outgroup are KU940219 for trnL-trnF, KU866563 for rps12-rpl20, and KU866559 for ndhF (329F, 927R). Different cpDNA haplotypes were identified by DnaSP 5.0 [42]. Geographic distribution of detected haplotypes was mapped using ArcMap 9.3 (ESRI, Redlands, CA, USA). Phylogenetic relationships among cpDNA haplotypes were associated via a median-joining network (MJN) method [43], which was constructed using Network 4.6.1.0 (available at http://www.fluxus-engineering.com/sharenet_rn.htm).

To maximize the differences of haplotype composition between geographical locations, spatial analysis of molecular variance algorithm was conducted using SAMOVA 1.0 procedure. This method divided the sampled populations into inferred partition ( groups) based on genetic differentiation and geographical homogeneity [44]. We initially set the values ranged from 2 to 12 by repeated simulated annealing process until largest index was obtained. However, when , at least one of the groups contained only one population, implying that information on the group structure would have been disappearing [45–48]. To maintain the grouping pattern and exclude configurations with a single population in any of the groups, finally valid range of was set as 2–5.

Parameters of Nei’s diversity, containing haplotype diversity () and nucleotide diversity (), were computed in Arlequin 3.1 program [49, 50]. To evaluate genetic variation within populations, among populations within groups, and among groups, AMOVA (analysis of molecular variance) was performed in the program Arlequin version 3.1 [50], and significance tests were conducted using 1,000 permutations. The statistics of average gene diversity and population differentiation, that is, the average gene diversity within populations (), total gene diversity across all populations (), the number of substitution types (), and genetic differentiation over populations (), were calculated using Permut version 1.0 [51] (available at http://www.pierroton.inra.fr/genetics/labo/Software/PermutCpSSR/index.html). To test for the existence of phylogeographical structure, the comparison of versus was conducted using Permut based on 1,000 random permutations of haplotypes across populations. If is significantly higher than , then genealogically closely related haplotypes tend to occur together within populations, indicating the presence of phylogeographical structure [51]. Significance of values was examined by the Test using the Haplonst procedure [51].

To estimate the divergence times of cpDNA haplotypes, the BEAST version 1.6.1 [52] was carried out for lineage analysis, following a Bayesian relaxed molecular clock. A Bayesian Markov Chain Monte Carlo (MCMC) procedure was run for 20,000,000 generations with a coalescent-based tree priority rule and a HKY substitution model. Given the uncertainty of the cpDNA mutation rate, the normal distribution priority with a mean value of 2 × 10⁻⁹ s s⁻¹ yr⁻¹ (substitutions per site per year) and a SD of 6.08 × 10⁻⁹ for most angiosperm species were adopted as a criterion for the analytical operation [24, 53], to cover the rate variation within the 95% range of the distribution for the present calculation of molecular divergence times of genotypes. The combined parameters were validated in Tracer version 1.5 [52] to check whether the parameter value was greater than 1 and the effective sample size value was larger than 200, as an indication that these sequence data were suitable for a relaxed molecular clock model. Trees were finally compiled in FigTree 1.3.1. The statistical support of clades depended on Bayesian posterior probability.

To investigate spatial genetic patterns of geographical distance among populations for each defined regional group, genetic landscape shape analysis was performed using Alleles In Space [54]. Primarily, a connective Delaunay triangulation network was constructed based on the geographic coordinates of all sampled points, according to the Watson (1992) and Brouns (2003) methods [55, 56]. Secondly, the average genetic distances resulting from the cpDNA sequence data were calculated among populations in the network. Afterwards, the genetic distance matrix was interpolated onto spatial grids of the connective network, forming a landscape shape. A final three-dimensional surface plot was generated, where the abscissa and ordinate signified geographical coordinates and the vertical coordinate were equivalent to genetic distance.

Possible demographic historical expansions for defined groups were inferred by mismatch distribution analysis using Arlequin with 1,000 parametric bootstrap replicates [57]. Populations that have experienced recent expansion are expected to have a unimodal curve in the pairwise mismatch distribution, whereas stable populations should have a bi- or multimodal shape [58, 59]. Neutrality tests (Tajima’s and Fu’s statistics) [60, 61] were used to detect further evidence for probable historical processes [62, 63]. A positive or negative value (not 0) may imply that a population underwent expansion or bottlenecks [64], while a negative value probably indicates a recent population expansion [61].

To explore the potential distribution changes of the species between the present-day and the Last Glacial Maximum (LGM; approximately 18–21 ka), the ecological niche modeling was performed in MAXENT version 3.3.1 [65], following the maximum entropy algorithm. We input an available coordinate set of occurrence points from the species’ entire geographical distribution in Eastern Central Asia based on the Chinese Virtual Herbarium (http://www.cvh.org.cn/cms/) and our field survey locations (Table 1). The current climatic envelopes were predicted based on WorldClim dataset [66] and the LGM layers in accordance with the Model for Interdisciplinary Research on Climate (MIROC) [67] and the Community Climate System Model (CCSM) [68]. We acquired 19 BIOCLIM variables to model ecological niche at 2.5 arcmin resolution from the WorldClim database (available at http://www.worldclim.org/download/). To avoid the potential problems of overfitting [69], we removed the highly correlated bioclimatic variables, retaining nine least correlated ones. Model performance was identified by the area under the receiver operating characteristic curve (AUC). The default parameters to form ecological niche modeling were chosen as follows: regularization multiplier = 1.0, percentage of the data set for train the model = 75%, percentage for testing = 25%, and the threshold for available presence data = 10 percentiles. In addition, a Jackknife test was used to measure the significant contributions of the BIOCLIM variables on the occurrence prediction for each model. Finally, potential distribution ranges for the present-day and LGM were projected in ArcMap 9.3 (ESRI, Redlands, CA, USA).

3. Results

3.1. Chloroplast Variation and Haplotype Geographical Distribution

The aligned sequence length of the three cpDNA regions was 2254 bp (923 bp for trnL-trnF, 728 bp for rps12-rpl20, and 603 bp for ndhF (329F, 927R)). Based on the total sequences, thirty-six polymorphic sites (thirty nucleotide substitutions and six indels) were obtained (Table ). A total of 24 haplotypes (H1–H24) were identified from 267 samples among 31 populations across the entire geographic range of C. spinosa (Figure 1(a), Table 1).

(a)

(b)

(c)

(d)

Figure 1: Geographical distribution and phylogenetic network of haplotypes in sampled locations of Capparis spinosa. (a) Distribution of the cpDNA haplotypes. The abbreviated letters of population localities correspond to those shown in Table 1. Pie charts show the different haplotypes and their proportions in each population; different colors are consistent with matching haplotypes in the figure legend. (b) Median-joining network reflecting the cpDNA haplotype relationships. Small, open red circles represent potentially intermediate genotypes with other mutational steps between real haplotypes; circle sizes are proportional to haplotype frequencies. (c, d) Two different kinds of habitats of the species: (c) eroding mountain slopes and (d) gravel deserts.

According to the results of SAMOVA, the 31 studied populations were divided into five defined phylogeographical groups: (I) the Western Himalayas Group, (II) the Eastern Pamir Group, (III) the Ili Valley Group, (IV) the north side of the Tianshan Mountains Group, and (V) the south side of the Tianshan Mountains Group (Table 1). The main genealogical relationships of the network connection diagram (Figure 1(b)) as well as the BEAST tree (Figure 2) suggested the similar divergence trend.

Figure 2: Phylogenetic chronogram of 24 haplotypes in C. spinosa generated from BEAST analysis. Divergence time (million years ago) of nodes are shown; blue horizontal bars indicating the 95% ranges of highest posterior density.

Distribution of the haplotypes showed high geographic structure. Haplotypes H22, H6, H7, H1, H8, H5, and H10 had the most frequent and widespread distributions: H22 only was distributed in the populations of Group I; H6 occurred in all populations of Group II; H7 only occurred in the populations of Group III; H1 and H8 were dominantly distributed in the populations of Group V, while H5 and H10 were mainly present in the populations of Group IV (Figure 1(a)).

3.2. Genetic Diversity and Population Genetic Structure

Genetic diversity analysis revealed that total genetic diversity (, s.e. = 0.0172) was much higher than average gene diversity (, s.e. = 0.0533), which indicated that there was considerable population differentiation across the distributional range. The Test showed that total (0.843, s.e. = 0.033) was significantly higher than (0.655, s.e. = 0.057) (, ), indicating cpDNA variation with signal phylogeographic structure for the species.

AMOVA revealed that a larger proportion of the variation was interpopulational than intrapopulational (85.20% versus 14.80%). Differences among geographical groups were significant, and AMOVA results showed that 71.83% (, ) of total variation existed between the five groups (Table 2). For the two genetic diversity indexes, haplotype diversities () within the 31 populations ranged from 0 to 0.8485, and nucleotide diversities () ranged from 0 to 0.0021. Considering the whole species, high degrees of haplotype diversity () and nucleotide diversity () were found. At the level of population, six populations (TRP, HM, WNS, GZ, SS, and KRL) in Group V showed considerably high degrees of genetic diversity (Table 1).

Table 2: Results of AMOVA analysis of molecular variance for populations and geographical groups of C. spinosa.

Source of variation
d.f.
SS
VC
PV (%)
Fixation index

Among groups
4
738.961
3.723
71.83
: 0.71831

Among populations within groups
26
175.341
0.693
13.37
: 0.47458

Within populations
236
181.039
0.767
14.80
: 0.85199

Total
266
1095.341
5.183

d.f.: degrees of freedom; SS: sum of squares; VC: variance components; PV: percentage of variation; : correlation of chlorotypes within populations relative to groups; : correlation within populations relative to the total; : correlation within groups relative to the total; ; 1000 permutations.

3.3. Molecular Divergence Time Estimates and Demographic History Analyses

The phylogenetic tree yielded by BEAST analysis was strongly coupled with geography, demonstrating that the coalescent sequences were suitable to be employed in phylogenetic and phylogeographic analysis. The results indicate that the identified cpDNA haplotypes are clustered into six main clades. Divergence times estimated from BEAST analysis were Pleistocene, and the main divergence between the geographical groups began at 1.179 Ma, in early Pleistocene (Figure 2).

The three-dimensional surface plot produced by genetic landscape shape analyses showed that spatial genetic distances gradually increased from the western populations to eastern populations (Figure 3). For mismatch analysis, the observed mismatch curve was not unimodal but strongly discordant with a model of sudden range expansion (Figure 4), indicating the stable populations. Additionally, no supporting evidence was provided by Tajima’s and Fu’s values for recent range expansion events.

Figure 3: Spatial genetic landscape shapes constructed by interpolation analysis for C. spinosa. The abscissae and ordinates correspond to geographical coordinates covering the entire distributional populations, and the vertical axes represent genetic distances.

Figure 4: Pairwise mismatch distributions of cpDNA sequences of C. spinosa.

3.4. Potential Distribution Modeling

Under the selected model, the test AUC value for the ENM of C. spinosa in the present-day was 0.995; and, under the MIROC climate scenario, the value was 0.997 in the LGM projection. SDM did not yield suitable or credible distribution in LGM based on the CCSM model. The Jackknife analysis of regularized training gain indicated that there were five environmental variables contributing more highly to the distribution model in the present-day: precipitation of the wettest month (22.64%), mean temperature of the coldest quarter (22.97%), precipitation of the coldest quarter (15.48%), mean temperature of the driest quarter (8.25%), and precipitation seasonality (6.37%). The five most contributive bioclimatic variables of the MIROC model in the LGM were annual mean temperature (24.32%), precipitation of the coldest quarter (20.18%), mean temperature of the coldest quarter (13.87%), mean temperature of the driest quarter (12.18%), and precipitation of the wettest month (11.89%). The range differences under the distinct climatic conditions of the two models revealed a diminished amount of suitable habitat for the species during the LGM period. In comparison with the present potential distribution, the MIROC scenario indicated that the high-altitude ranges of the species in the Western Himalayas, eastern Pamir, and western Tianshan Mountains were diminished (Figure 5), according to the LGM modeling distribution. However, some localities in the middle Tianshan, such as WNS and KRL, were not affected by this in the LGM model (Figure 5), and thus midelevation gorges environments might have provided suitable survival habitats for the species in the glacial epoch.

Figure 5: Species distribution modeling of C. spinosa during present-day and LGM periods in arid Eastern Central Asia. Potential distribution in the LGM epoch is based on scenarios of the MIROC model. Yellow localities: populations sensitive to glaciation in the high mountains. Blue localities: populations at middle elevations. Green localities: gorges served as refugia in the mid-Tianshan. Red localities: populations of low-altitude desert region.

4. Discussion

4.1. Genetic Variation

At the species level, a high total gene diversity () was detected. Compared with other species where three intergenic spacers have been sequenced, gene diversity in C. spinosa was higher, for instance, Juniperus sabina () [70] and Gymnocarpos przewalskii () [29]. The results of genetic diversity analysis also give indication of a high level of total haplotype diversity () (Table 1). High values of the diversity indices might have been promoted by the accumulation of variation in C. spinosa, which, to a great extent, may have been associated with differences in geological and topographical conditions, complex habitat patterns, and gradually intensive aridification occurring in Eastern Central Asia during the Quaternary period.

Dramatic genetic variation among the total populations was shown, accounting for a far greater proportion than that within populations (Table 2). The genetic differences among populations are believed to be closely related to the degree of gene exchange [71], which depends on seed transfer for cpDNA in most angiosperms [72]. The seed dispersal modes of C. spinosa mainly include anemochory, myrmecochory, and endozoochory by rodent vector [73]. In the current survey regions, the mountains and sprawling desert zones isolated habitat into disjunctive geographical units. In this situation, gene flow between populations was potentially restricted in short distance. Especially in small sized populations, individuals might be susceptible to inbreeding. Therefore, gene flow was limited attributing to the incapacity of seed dispersal across impassable physical barriers by anemochory or zoochory, resulting in population differentiation.

4.2. Influence of Pleistocene Geological and Climatic Changes on Phylogeographical Structure

The timescales of intraspecific divergences are usually related to historical geologic and climatic events [3, 12]. The estimated genetic divergence time between the five geographical groups started at approximately 1.179 million years ago based on molecular dating results of the BEAST analysis (Figure 2), which falls into early Pleistocene, the phase of intense uplift of Himalayas [74]. From early to middle Pleistocene, most large-scale piedmont glaciers developed on the northern slope of the Himalayas [75], which may have been the impetus causing haplotypes H23, H24, and H22 belonging to the Western Himalayan populations to branch off the phylogenetic tree (0.808–0.621 Ma) (Figure 2). At the intersection of the Himalaya, Kunlun, and Tianshan Mountains, the Pamir uplift was structured from at least the end of the Paleogene. In the Quaternary stage, owing to at least three large-scale glaciations, piedmont moraines intermittently formed in this region [76]. In connection with these fluctuations, the Western Himalayan and Eastern Pamirian populations from high-altitude localities experienced large-scale differentiation, accompanied with adaptation to colder climates, resulting in different degrees of genetic diversity and conspicuous population differences. We found that haplotypes H22 and H6 were, respectively, distributed over every sampled location in Western Himalayan and Eastern Pamir region, while the rare haplotypes H23, H24, and H21 were randomly distributed in single populations. Although BEAST analysis and network diagram did not clearly show that haplotype H22 clustered into the same clade with H23 and H24, it occupies a key transitional position (Figures 1(b) and 2). Genetically H22 contained a haplotype approximating those in the Tianshan region, while geographically it is distributed in the Western Himalayan region, indicating a possible common Tethyan origin and then allopatric divergence and local evolution.

In comparison with the formation of the Himalayas and Pamirs, neotectonic movement of Tianshan Mountains was somewhat lagging, occurring in Pliocene and intensively uplifting in early Pleistocene [7, 8]. Since Quaternary diastrophism and climatic changes, together with more than three glacial-interglacial cycles [77], not only changed the topography but also resulted in temperature and precipitation differences between western and eastern ranges, north and south slopes, and high and low altitudes [74], rapid uplift of the Tianshan range as a barrier would, of necessity, have changed the local atmospheric circulation and screened the sources of moist airflow, intensifying the arid climate in adjacent regions. Later, since the mid-Pleistocene, enhanced drought and cold occurred in the Tianshan and the surroundings at 0.6–0.2 Ma [78], as evidenced by loess sediments on northern slope of the Tianshan Mountains [74]. During this period, a weakened southwest monsoon and strengthened plateau winter monsoon resulted in increasing aridification and desert expansion [79]. Compared with the geographic patterns under present-day climatic conditions, distribution ranges of C. spinosa during the LGM period were contracted (Figure 5). Although distribution areas of the species changed, the main distribution in the low-altitude eastern desert region was stable during glacial periods and only marginally influenced. The great spatial and genetic distances among locations in the eastern Tianshan and Hexi Corridor can be perceived in the 3D surface plots (Figure 3). This characteristic of the spatial genetic landscape analysis also demonstrates that the species has had relatively stable habitats in the low-altitude desert region during an extensive evolutionary history. At the stage of the mid-Pleistocene period, the continuously expanding Taklimakan [80] and intensive aridification accelerated the worsening drought and heat for the lower elevations of the eastern Tianshan, Turpan-Hami Basin, and Hexi Corridor populations, which must have influenced recent divergence among these areas. In general, genetic divergence among geographic groups in Pleistocene appears to have been driven mainly by intense uplift of the Himalayas and Tianshan Mountains. The cold-dry to warm-humid climatic cycles during the late Quaternary are inferred to have promoted genetic divergence within groups.

4.3. Potential Glacial Refugia and Demographic Dispersal in Past Scenarios

According to previous glaciological research, maximum glaciation was likely to have occurred in the early phase of the mid-Pleistocene (0.8–0.6 Ma) in Eastern Central Asia [74, 77]. During this period and the next two glacial stages, owing to large-scale glacial advances, glaciers in Tianshan at one time extended as far down as the piedmont belts [74]. In comparison with the present distribution range, the LGM localities of C. spinosa were reduced in the high-altitude ranges of the Western Himalayas, Pamirs, and western Tianshan Mountains (Figure 5), as a result of frigid climate in the wake of glacial advances. The increasingly severe climate in the glacial epoch would have damaged the ecological habitats of some xeric plant species [81]; however, it seems likely that C. spinosa could have survived in some stable gorges of Tianshan Mountains. These suitable refugia were less influenced by climatic oscillations and therefore could maintain richer degrees of genetic diversity among populations [25, 82]. Refugial zones are generally recognizable because of high levels of genetic diversity and unique genotypes in species populations [12, 22] and are usually located in the middle altitude of mountains, which could supply warmer habitats for species to resist cold glacial climates [83]. A number of glacial refugia have been revealed worldwide on the basis of molecular approaches [84–86]. However, refugia localities in the Tianshan and adjacent arid regions have remained rather cryptic, as there seems to be a lack of applicable studies [23]. As a representative xerophyte, C. spinosa selected relatively warm gorges (WNS and KRL) in the middle Tianshan Mountains as potential refugia during glacial epochs, in accordance with the high levels of genetic diversity within these populations.

Unlike the traditional northward expansions during warm interglacial periods [18, 87], in the current study, the dispersal of C. spinosa was profoundly affected not only by temperature variation but also by aridification in Eastern Central Asia; thus xeric species may well prefer relatively arid areas during the warm and moist interglacial stages. Considering the network diagram, genotypes H1 and H8 were found with high rates of incidence on the south side of Tianshan group, indicating that these two were predominant haplotypes in Group V (Figure 1). These characteristics imply that the species should have undergone large-scale interglacial expansions from glacial refugia in this region. However, multimodal mismatch distribution shape (Figure 4) and positive Fu’s value failed to provide any significant evidence for recent range expansion. According to Printzen et al. (2003), extensive intraspecific differentiation is probably connected with range fragmentation or long-distance dispersal [88]. Evidenced from our field observations, C. spinosa is a xeromorphic species that prefers to distribute itself in eroding rocky mountains and arid piedmont proluvial or alluvial gravelly fans (Figures 1(c) and 1(d)). In the end of early Pleistocene and the early mid-Pleistocene, intensified aridification and the formation of Taklimakan desert resulted in dry and rainless conditions in piedmont zones of Tianshan [74, 89]. However, in contrast with the cold climate in glacial phase, in the interglacial and postglacial period the climate became warm and humid. At this time, glacier-melt water increased the volume of river runoff on the edge of the Tarim Basin [74]. Thus, the movement of the wind and river runoff provided appropriate conditions for survival and seed dispersal of C. spinosa. These could harbor relatively stable habitats, preserving the species’ existing genetic variation and genetic differences among populations. There is considerable evidence, such as increasing values, palynological components, multiple alluvial sequences, and lake sediments, that reveals the alternation of glacials and interglacials occurring in the Tianshan Mountains during the Pleistocene [74, 90, 91]. This circularly fluctuant climate propelled demographic dispersal in the interglacial intervals. Interestingly, unlike most other species with decreased levels of genetic diversity related to desert expansion during interglacial and postglacial colonization [29, 30], as xeric Tethyan relic, C. spinosa showed successively increased levels of genetic diversity and population sizes, especially when it contacted the more droughty and hot gravelly deserts in Turpan-Hami Basin and Hexi Corridor. One reasonable explanation for this phenomenon is the probability that the suitable warm temperate continental climate in eastern gravel deserts is similar to the ancient hot and dry subtropical summer climate of the Tethyan Flora, preserving the characteristics of the original climate and environment. The other corollary is that the desert region could play the role of a shelter influenced only marginally by glaciation.

5. Conclusions

The current study focuses on how complicated Pleistocene geological and cyclical climatological events influenced the phylogeographic structure and genetic divergence of C. spinosa in arid Eastern Central Asia. Strongly spatial phylogeographic patterns were documented in the species. Intense uplift of the Himalaya and Tianshan Mountains around the early Pleistocene, coupled with a cold-dry climate and consequent aridification during the stage of Quaternary glaciation in this survey region, played important roles both in triggering and in shaping the current phylogeographic structure of the species group. There were at least three glacial-interglacial cycles in the Himalayas, as well as the Pamir and Tianshan Mountains during the Pleistocene. During glacial epochs, potential refugia were inferred for the gorges of the middle Tianshan Mountains, while, in interglacials and the postglacial period, the species experienced dispersal. The phylogeography of C. spinosa in the present research provides the basis for future studies on how xerophilous plant species across both arid high-altitude rocky mountains (Figure 1(c)) and low-altitude gravel deserts (Figure 1(d)) responded in a stepwise fashion to Quaternary geologic and climatic events.

Competing Interests

The authors declare that there are no competing interests regarding the publication of this paper.

Acknowledgments

The authors are grateful to Dr. Stewart C. Sanderson in the Shrub Sciences Laboratory, Rocky Mountain Research Station, US Department of Agriculture, Utah, USA, for his helpful suggestions and English improvement of the paper. This study was financially supported by China National Key Basic Research Program (2014CB954201), CAS Key Laboratory of Biogeography and Bioresource in Arid Land, and Xinjiang Institute of Ecology and Geography, Chinese Academy of Sciences.

References

C. B. Cox and P. D. Moore, Biogeography: an Ecological and Evolutionary Approach, Blackwell Science, Oxford, UK, 2000.

R. I. Milne and R. J. Abbott, “The origin and evolution of tertiary relict floras,” Advances in Botanical Research, vol. 38, pp. 281–314, 2002.

G. M. Hewitt, “The genetic legacy of the quaternary ice ages,” Nature, vol. 405, no. 6789, pp. 907–913, 2000.

G. M. Hewitt, “Genetic consequences of climatic oscillations in the Quaternary,” Philosophical Transactions of the Royal Society B: Biological Sciences, vol. 359, no. 1442, pp. 183–195, 2004.

E. Palamarev, “Paleobotanical evidences of the Tertiary history and origin of the Mediterranean sclerophyll dendroflora,” Plant Systematics and Evolution, vol. 162, no. 1, pp. 93–107, 1989.

H. Sun, “Tethys retreat and Himalayas-Hengduanshan Mountains uplift and their significance on the origin and development of the Sino-Himalayan elements and Alpine flora,” Acta Botanica Yunnanica, vol. 24, no. 3, pp. 273–288, 2002 (Chinese).

J.-M. Sun, R.-X. Zhu, and J. Bowler, “Timing of the Tianshan Mountains uplift constrained by magnetostratigraphic analysis of molasse deposits,” Earth and Planetary Science Letters, vol. 219, no. 3-4, pp. 239–253, 2004.

L.-N. Wang, J.-Q. Ji, D.-X. Sun et al., “The uplift history of south-western Tianshan-Implications from AFT analysis of detrital samples,” Chinese Journal of Geophysics, vol. 53, no. 4, pp. 931–945, 2010 (Chinese).

A. Takhtajan, T. J. Crovello, and A. Cronquist, Floristic Regions of the World, University of California Press, Berkeley, Calif, USA, 1986.

Z.-Y. Wu, H. Sun, Z.-K. Zhou, D.-Z. Li, and H. Peng, Floristics of Seed Plants from China, Science Press, Beijing, China, 2010 (Chinese).

L. L. Anderson, F. S. Hu, D. M. Nelson, R. J. Petit, and K. N. Paige, “Ice-age endurance: DNA evidence of a white spruce refugium in Alaska,” Proceedings of the National Academy of Sciences of the United States of America, vol. 103, no. 33, pp. 12447–12450, 2006.

R. J. Petit, I. Aguinagalde, J.-L. de Beaulieu et al., “Glacial refugia: hotspots but not melting pots of genetic diversity,” Science, vol. 300, no. 5625, pp. 1563–1565, 2003.

J. C. Avise, “The history and purview of phylogeography: a personal reflection,” Molecular Ecology, vol. 7, no. 4, pp. 371–379, 1998.

A. Terrab, P. Schönswetter, S. Talavera, E. Vela, and T. F. Stuessy, “Range-wide phylogeography of Juniperus thurifera L., a presumptive keystone species of western Mediterranean vegetation during cold stages of the Pleistocene,” Molecular Phylogenetics and Evolution, vol. 48, no. 1, pp. 94–102, 2008.

M. J. Hickerson, B. C. Carstens, J. Cavender-Bares et al., “Phylogeography's past, present, and future: 10 years after Avise, 2000,” Molecular Phylogenetics and Evolution, vol. 54, no. 1, pp. 291–301, 2010.

J.-Q. Liu, Y.-S. Sun, X.-J. Ge, L.-M. Gao, and Y.-X. Qiu, “Phylogeographic studies of plants in China: advances in the past and directions in the future,” Journal of Systematics and Evolution, vol. 50, no. 4, pp. 267–275, 2012.

A. C. Newton, T. R. Allnutt, A. C. M. Gillies, A. J. Lowe, and R. A. Ennos, “Molecular phylogeography, intraspecific variation and the conservation of tree species,” Trends in Ecology and Evolution, vol. 14, no. 4, pp. 140–145, 1999.

I. V. Bartish, J. W. Kadereit, and H. P. Comes, “Late Quaternary history of Hippophaë rhamnoides L. (Elaeagnaceae) inferred from chalcone synthase intron (Chsi) sequences and chloroplast DNA variation,” Molecular Ecology, vol. 15, no. 13, pp. 4065–4083, 2006.

X.-J. Ge, C.-C. Hwang, Z.-H. Liu et al., “Conservation genetics and phylogeography of endangered and endemic shrub Tetraena mongolica (Zygophyllaceae) in Inner Mongolia, China,” BMC Genetics, vol. 12, article 1, 2011.

Y.-Y. Liao, A. W. Gichira, Q.-F. Wang, and J.-M. Chen, “Molecular phylogeography of four endemic Sagittaria species (Alismataceae) in the Sino-Japanese Floristic Region of East Asia,” Botanical Journal of the Linnean Society, vol. 180, no. 1, pp. 6–20, 2016.

P. Schönswetter, I. Stehlik, R. Holderegger, and A. Tribsch, “Molecular evidence for glacial refugia of mountain plants in the European Alps,” Molecular Ecology, vol. 14, no. 11, pp. 3547–3555, 2005.

J. R. Stewart, A. M. Lister, I. Barnes, and L. Dalén, “Refugia revisited: individualistic responses of species in space and time,” Proceedings of the Royal Society B: Biological Sciences, vol. 277, no. 1682, pp. 661–671, 2010.

Y.-X. Qiu, C.-X. Fu, and H. P. Comes, “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, vol. 59, no. 1, pp. 225–244, 2011.

D.-R. Jia, R.-J. Abbott, T.-L. Liu, K.-S. Mao, I. V. Bartish, and J.-Q. Liu, “Out of the Qinghai-Tibet Plateau: evidence for the origin and dispersal of Eurasian temperate plants from a phylogeographic study of Hippophaë rhamnoides (Elaeagnaceae),” New Phytologist, vol. 194, no. 4, pp. 1123–1133, 2012.

H. P. Comes and J. W. Kadereit, “The effect of quaternary climatic changes on plant distribution and evolution,” Trends in Plant Science, vol. 3, no. 11, pp. 432–438, 1998.

O. Bettin, C. Cornejo, P. J. Edwards, and R. Holderegger, “Phylogeography of the high alpine plant Senecio halleri (Asteraceae) in the European Alps: in situ glacial survival with postglacial stepwise dispersal into peripheral areas,” Molecular Ecology, vol. 16, no. 12, pp. 2517–2524, 2007.

X.-J. Ge, Y.-C. Chiang, C.-H. Chou, and T.-Y. Chiang, “Nested clade analysis of Dunnia sinensis (Rubiaceae), a monotypic genus from China based on organelle DNA sequences,” Conservation Genetics, vol. 3, no. 4, pp. 351–362, 2002.

Z. Su, M. Zhang, and S. C. Sanderson, “Chloroplast phylogeography of Helianthemum songaricum (Cistaceae) from northwestern China: implications for preservation of genetic diversity,” Conservation Genetics, vol. 12, no. 6, pp. 1525–1537, 2011.

S. M. Ma, M. L. Zhang, and S. C. Sanderson, “Phylogeography of the rare Gymnocarpos przewalskii (Caryophyllaceae): indications of multiple glacial refugia in north-western China,” Australian Journal of Botany, vol. 60, no. 1, pp. 20–31, 2012.

H.-H. Meng and M.-L. Zhang, “Diversification of plant species in arid Northwest China: species-level phylogeographical history of Lagochilus Bunge ex Bentham (Lamiaceae),” Molecular Phylogenetics and Evolution, vol. 68, no. 3, pp. 398–409, 2013.

M. Jacobs, “The genus Capparis (Capparaceae) from the Indus to the Pacific,” Blumea, vol. 12, pp. 385–541, 1965.

M.-L. Zhang and G. C. Tucker, “Capparaceae,” in Flora of China, Z.-Y. Wu and P. H. Raven, Eds., pp. 433–450, Science Press/Missouri Botanical Garden Pres, Beijing, China, 2008.

Z.-X. An, “Capparis L,” in Flora Xinjiangensis, Z.-M. Mao, Ed., pp. 35–36, Xinjiang Science, Technology and Hygiene Publishing House, Ürümqi, China, 1995.

M. Zohary, “The species of Capparis in the Mediterranean and the near eastern countries,” Bulletin of the Research Council of Israel, vol. 8, pp. 49–64, 1960.

R. N. Higton and J. R. Akeroyd, “Variation in Capparis spinosa L. in Europe,” Flora Europaea, vol. 106, no. 2, pp. 104–112, 1991.

T. Zhang and D.-Y. Tan, “An examination of the function of male flowers in an andromonoecious shrub Capparis spinosa,” Journal of Integrative Plant Biology, vol. 51, no. 3, pp. 316–324, 2009.

E. Saadaoui, A. Guetat, N. Tlili, M. El Gazzah, and A. Khaldi, “Subspecific variability of Tunisian wild populations of Capparis Spinosa L,” Journal of Medicinal Plants Research, vol. 5, no. 17, pp. 4339–4348, 2011.

H. Nosrati, M. A. H. Feizi, M. Mazinani, and A. R. Haghighi, “Effect of population size on genetic variation levels in Capparis spinosa (Capparaceae) detected by RAPDs,” EurAsian Journal of BioSciences, vol. 6, pp. 70–75, 2012.

Ö. Özbek and A. Kara, “Genetic variation in natural populations of Capparis from Turkey, as revealed by RAPD analysis,” Plant Systematics and Evolution, vol. 299, no. 10, pp. 1911–1933, 2013.

S. O. Rogers and A. J. Bendich, “Extraction of DNA from milligram amounts of fresh, herbarium and mummified plant tissues,” Plant Molecular Biology, vol. 5, no. 2, pp. 69–76, 1985.

J. D. Thompson, T. J. Gibson, F. Plewniak, F. Jeanmougin, and D. G. Higgins, “The CLUSTAL X windows interface: flexible strategies for multiple sequence alignment aided by quality analysis tools,” Nucleic Acids Research, vol. 25, no. 24, pp. 4876–4882, 1997.

P. Librado and J. Rozas, “DnaSP v5: a software for comprehensive analysis of DNA polymorphism data,” Bioinformatics, vol. 25, no. 11, pp. 1451–1452, 2009.

H.-J. Bandelt, P. Forster, and A. Röhl, “Median-joining networks for inferring intraspecific phylogenies,” Molecular Biology and Evolution, vol. 16, no. 1, pp. 37–48, 1999.

I. Dupanloup, S. Schneider, and L. Excoffier, “A simulated annealing approach to define the genetic structure of populations,” Molecular Ecology, vol. 11, no. 12, pp. 2571–2581, 2002.

J. Godbout, J. P. Jaramillo-Correa, J. Beaulieu, and J. Bousquet, “A mitochondrial DNA minisatellite reveals the postglacial history of jack pine (Pinus banksiana), a broad-range North American conifer,” Molecular Ecology, vol. 14, no. 11, pp. 3497–3512, 2005.

D. Magri, G. G. Vendramin, B. Comps et al., “A new scenario for the Quaternary history of European beech populations: palaeobotanical evidence and genetic consequences,” New Phytologist, vol. 171, no. 1, pp. 199–221, 2006.

M. Heuertz, S. Fineschi, M. Anzidei et al., “Chloroplast DNA variation and postglacial recolonization of common ash (Fraxinus excelsior L.) in Europe,” Molecular Ecology, vol. 13, no. 11, pp. 3437–3452, 2004.

T. Iwasaki, K. Aoki, A. Seo, and N. Murakami, “Comparative phylogeography of four component species of deciduous broad-leaved forests in Japan based on chloroplast DNA variation,” Journal of Plant Research, vol. 125, no. 2, pp. 207–221, 2012.

M. Nei, Molecular Evolutionary Genetics, Columbia University Press, New York, NY, USA, 1987.

L. Excoffier, G. Laval, and S. Schneider, “Arlequin (version 3.0): an integrated software package for population genetics data analysis,” Evolutionary Bioinformatics, vol. 1, pp. 47–50, 2005.

O. Pons and R. J. Petit, “Measuring and testing genetic differentiation with ordered versus unordered alleles,” Genetics, vol. 144, no. 3, pp. 1237–1245, 1996.

A. J. Drummond and A. Rambaut, “BEAST: Bayesian evolutionary analysis by sampling trees,” BMC Evolutionary Biology, vol. 7, no. 1, article 214, 2007.

K. H. Wolfe, W. H. Li, and P. M. Sharp, “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, vol. 84, no. 24, pp. 9054–9058, 1987.

M. P. Miller, “Alleles In Space (AIS): computer software for the joint analysis of interindividual spatial and genetic information,” Journal of Heredity, vol. 96, no. 6, pp. 722–724, 2005.

D. F. Watson, Contouring: A Guide to the Analysis and Display of Spatial Data (with Programs on Diskette), Pergamon Press, Oxford, UK, 1992.

G. Brouns, A. De Wulf, and D. Constales, “Delaunay triangulation algorithms useful for multibeam echosounding,” Journal of Surveying Engineering, vol. 129, no. 2, pp. 79–84, 2003.

S. Schneider and L. Excoffier, “Estimation of past demographic parameters from the distribution of pairwise differences when the mutation rates vary among sites: application to human mitochondrial DNA,” Genetics, vol. 152, no. 3, pp. 1079–1089, 1999.

M. Slatkin and R. R. Hudson, “Pairwise comparisons of mitochondrial DNA sequences in stable and exponentially growing populations,” Genetics, vol. 129, no. 2, pp. 555–562, 1991.

A. R. Rogers and H. Harpending, “Population growth makes waves in the distribution of pairwise genetic differences,” Molecular Biology and Evolution, vol. 9, no. 3, pp. 552–569, 1992.

F. Tajima, “Statistical method for testing the neutral mutation hypothesis by DNA polymorphism,” Genetics, vol. 123, no. 3, pp. 585–595, 1989.

Y.-X. Fu, “Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection,” Genetics, vol. 147, no. 2, pp. 915–925, 1997.

L. A. Hurtado, T. Erez, S. Castrezana, and T. A. Markow, “Contrasting population genetic patterns and evolutionary histories among sympatric Sonoran Desert cactophilic Drosophila,” Molecular Ecology, vol. 13, no. 6, pp. 1365–1375, 2004.

C. I. Smith and B. D. Farrell, “Range expansions in the flightless longhorn cactus beetles, Moneilema gigas and Moneilema armatum, in response to Pleistocene climate changes,” Molecular Ecology, vol. 14, no. 4, pp. 1025–1044, 2005.

F. Tajima, “The amount of DNA polymorphism maintained in a finite population when the neutral mutation rate varies among sites,” Genetics, vol. 143, no. 3, pp. 1457–1465, 1996.

S. J. Phillips, R. P. Anderson, and R. E. Schapire, “Maximum entropy modeling of species geographic distributions,” Ecological Modelling, vol. 190, no. 3-4, pp. 231–259, 2006.

R. J. Hijmans, S. E. Cameron, J. L. Parra, P. G. Jones, and A. Jarvis, “Very high resolution interpolated climate surfaces for global land areas,” International Journal of Climatology, vol. 25, no. 15, pp. 1965–1978, 2005.

H. Hasumi and S. Emori, K-1 Coupled GCM (MIROC) Description, Center for Climate System Research, University of Tokyo, Tokyo, Japan, 2004.

W. D. Collins, C. M. Bitz, M. L. Blackmon et al., “The Community Climate System Model version 3 (CCSM3),” Journal of Climate, vol. 19, no. 11, pp. 2122–2143, 2006.

A. T. Peterson and Y. Nakazawa, “Environmental data sets matter in ecological niche modelling: an example with Solenopsis invicta and Solenopsis richteri,” Global Ecology and Biogeography, vol. 17, no. 1, pp. 135–144, 2008.

Y.-P. Guo, R. Zhang, C.-Y. Chen, D.-W. Zhou, and J.-Q. Liu, “Allopatric divergence and regional range expansion of juniperus sabina in China,” Journal of Systematics and Evolution, vol. 48, no. 3, pp. 153–160, 2010.

M. Slatkin, “Gene flow and the geographic structure of natural populations,” Science, vol. 236, no. 4803, pp. 787–792, 1987.

D. E. McCauley, “The use of chloroplast DNA polymorphism in studies of gene flow in plants,” Trends in Ecology & Evolution, vol. 10, no. 5, pp. 198–202, 1995.

T. Zhang, The Reproductive Ecology on Capparis spinosa L. (Capparaceae), Xinjiang Agricultural University, Urumqi, China, 2008 (Chinese).

Y.-F. Shi, Z.-J. Cui, and Z. Su, The Quaternary Glaciations and Environmental Variations in China, Hebei Science and Technology Press, Shijiazhuang, China, 2006 (Chinese).

S. C. Porter, “Quaternary glacial record in Swat Kohistan, West Pakistan,” Geological Society of America Bulletin, vol. 81, no. 5, pp. 1421–1446, 1970.

Xinjiang Integrated Survey Group of Chinese Academy of Sciences, Xinjiang Geomorphology, Science Press, Beijing, China, 1978 (Chinese).

X. Xu, A. Kleidon, L. Miller, S. Wang, L. Wang, and G. Dong, “Late Quaternary glaciation in the Tianshan and implications for palaeoclimatic change: a review,” Boreas, vol. 39, no. 2, pp. 215–232, 2010.

J.-J. Li, “The patterns of environmental changes since late Pleistocene in Northwestern China,” Quaternary Science Reviews, vol. 3, pp. 197–204, 1990.

X. Wu, F. Wang, and Z.-S. An, “The phase and height of tibetan plateau rise in late cenozoic,” in Loess Quaternary Geology Global Changes, D.-S. Liu and Z.-S. An, Eds., pp. 1–13, Scientific Publishing House, Beijing, China, 1992.

Z. Wu, “Approach to the genesis of the Taklamakan desert,” Acta Geographica Sinica, vol. 36, no. 3, pp. 280–291, 1981 (Chinese).

Y.-X. Liu, “A study on origin and formation of the Chinese desert floras,” Acta Phytotaxonmica Sinica, vol. 33, no. 2, pp. 131–143, 1995 (Chinese).

A. Widmer and C. Lexer, “Glacial refugia: sanctuaries for allelic richness, but not for gene diversity,” Trends in Ecology and Evolution, vol. 16, no. 6, pp. 267–269, 2001.

Y.-J. Zhang, M. Stöck, P. Zhang, X.-L. Wang, H. Zhou, and L.-H. Qu, “Phylogeography of a widespread terrestrial vertebrate in a barely-studied Palearctic region: green toads (Bufo viridis subgroup) indicate glacial refugia in Eastern Central Asia,” Genetica, vol. 134, no. 3, pp. 353–365, 2008.

R. J. Abbott, L. C. Smith, R. I. Milne, R. M. M. Crawford, K. Wolff, and J. Balfour, “Molecular analysis of plant migration and refugia in the Arctic,” Science, vol. 289, no. 5483, pp. 1343–1346, 2000.

T. Lacourse, R. W. Mathewes, and D. W. Fedje, “Late-glacial vegetation dynamics of the Queen Charlotte Islands and adjacent continental shelf, British Columbia, Canada,” Palaeogeography, Palaeoclimatology, Palaeoecology, vol. 226, no. 1-2, pp. 36–57, 2005.

L.-Y. Wang, R. J. Abbott, W. Zheng, P. Chen, Y. Wang, and J. Liu, “History and evolution of alpine plants endemic to the Qinghai-Tibetan Plateau: Aconitum gymnandrum (Ranunculaceae),” Molecular Ecology, vol. 18, no. 4, pp. 709–721, 2009.

G. M. Hewitt, “Quaternary phylogeography: the roots of hybrid zones,” Genetica, vol. 139, no. 5, pp. 617–638, 2011.

C. Printzen, S. Ekman, and T. Tønsberg, “Phylogeography of Cavernularia hultenii: evidence of slow genetic drift in a widely disjunct lichen,” Molecular Ecology, vol. 12, no. 6, pp. 1473–1486, 2003.

H.-B. Zheng, C. M. Powell, K. Butcher, and J. Cao, “Late Neogene loess deposition in southern Tarim Basin: tectonic and palaeoenvironmental implications,” Tectonophysics, vol. 375, no. 1–4, pp. 49–59, 2003.

Q. Wen and Y. Shi, “The quaternary climo-environment changes in Chaiwopu basin of Xinjiang region,” Chinese Geographical Science, vol. 3, no. 2, pp. 147–158, 1993.

H. Lu, D. W. Burbank, and Y. Li, “Alluvial sequence in the north piedmont of the Chinese Tian Shan over the past 550 kyr and its relationship to climate change,” Palaeogeography, Palaeoclimatology, Palaeoecology, vol. 285, no. 3-4, pp. 343–353, 2010.