Phylogeography and genetic effects of habitat fragmentation on endemic Urophysa (Ranunculaceae) in Yungui Plateau and adjacent regions
Phylogeography and genetic effects of habitat fragmentation on endemic Urophysa (Ranunculaceae) in Yungui Plateau and adjacent regions
Deng-Feng Xie 0 1
Min-Jie Li 0 1
Jin-Bo Tan 0 1
Megan Price 1
Qun-Ying Xiao 0 1
Song- Dong Zhou 0 1
Yan Yu 0 1
Xing-Jin He 0 1
0 Key Laboratory of Bio-Resources and Eco-Environment of Ministry of Education, College of Life Sciences, Sichuan University , Chengdu , China , 2 Sichuan Key Laboratory of Conservation Biology on Endangered Wildlife, College of Life Sciences, Sichuan University , Chengdu, Sichuan , China
1 Editor: Tzen-Yuh Chiang, National Cheng Kung University , TAIWAN
Urophysa is a Chinese endemic genus with only two species (U. rockii and U. henryi) distributed in Yungui Plateau (Guizhou Province) and adjacent regions (i.e., Provinces of Hunan, Hubei, Chongqing and Sichuan). The aim of this study was to determine the genetic diversity and population differentiation within Urophysa and investigate the effect of the Yungui Plateau uplift and climate oscillations on evolution of Urophysa. In this study, micro-morphological characteristics, nine microsatellite loci (SSR), two nuclear loci (ITS and ETS) and two chloroplast fragments (psbA-trnH and trnL-trnF) were used to analyze the phylogenetic relationships and assess genetic and phylogeographical structure of Urophysa. Isolation by distance (IBD) was performed to research the effects of geographical isolation. We detected high genetic diversity at the species level but low genetic diversity within populations. Striking genetic differentiation (AMOVA) among populations and a significant phylogeographical structure (NST > GST, p < 0.01) were detected among U. henryi populations, along with significant effects of isolation by distance (IBD). Molecular clock estimation using calibration strategy and cpDNA substitution rate indicated that the divergence of U. henryi occurred during late Miocene to early Quaternary, when the orogeny of Yungui Plateau was violent. U. rockii originated at the early Quaternary and further differentiated at early Pleistocene. Our results suggested that habitat fragmentation played an important role in the genetic diversity and population differentiation of U. rockii and U. henryi. Heterogenous geomorphological configuration and complicated environment resulted from rapid uplift of the Yungui Plateau were inferred as important incentives for the modern phylogeograhpical pattern and species divergence of Urophysa. The geographical isolation, limited gene flow, specialized morphologies and the Pleistocene climatic oscillation greatly contributed to the allopatric divergence of U. rockii. Significant genetic drift and inbreeding were detected in these two species, in situ measures should be implemented to protect them.
Geological history and climate oscillations are important drivers in the evolution and genetic
structure of plant species [
]. Previous researches have indicated that landscape
heterogeneity caused by the uplift of the Qinghai-Tibetan Plateau (QTP) greatly contributed to the
evolution of many plants (i.e., Aconitum gymnandrum, Taxus wallichiana and Rhodiola kirilowii) in
the QTP [3±6]. However, it is uncertain how landscape heterogeneity from uplift in adjacent
areas (to the QTP) affected plant evolution [
], such as the Yungui Plateau. The Indian plate
has kept moving northward since the Cenozoic, colliding with the Eurasian Plate between 40
and 50 million years ago (Mya). As a result, the orogeny of the QTP is violent, which induced
the continuing uplift of the Yungui Plateau, contributing to the formation of its unique
geomorphological characteristics and heterogeneous landscape environment . The
heterogeneous landscape was typically characterized by low soil water content, periodic water
deficiency, and poor nutrient availability, which exert strong selective forces on plant
evolution, resulting in remarkably high species richness and endemism in the Yungui plateau [
Underlying this species differentiation of the plateau, there was a massive divergence in
population genetics and a promotion of ecological diversity. For example, the Yungui Plateau and
its adjacent regions have been regarded as an important center of origin for the East Asiatic
]. In addition, the unique characteristics of the plateau have given rise and refuge to a
variety of endemic plants [
], and now is a refuge for plants that are threatened elsewhere.
Orogeny always resulted in regions of most plants shifted and divided populations into
different spatial-temporal scales, along with breaking up of one patch of habitat into several
smaller patches [
]. This phenomenon was called habitat fragmentation, which plays a
major role in threatening plant species survival. Habitat fragmentation is often considered to
be the main driving force for local species extinction [
] as large populations split into
smaller populations with increasing geographic isolation [
]. The spatial isolation of
populations may restrict connectivity, resulting in low levels of gene flow between fragments, with
subsequently lower genetic diversity and higher genetic differentiation in/among remnant
1, 18, 19
]. Habitat fragmentation can reduce the fitness of remnant plants by
affecting population genetics and dynamics in several ways [
], including genetic erosion
] population divergence and random genetic drift [
]. Smaller isolated populations
often experience greater inbreeding, reduced reproduction rate and offspring survival .
Conversely, genetic variation could be maintained or even increased in fragmented
], and habitat fragmentation may play an important role in plant divergence and
allopatric speciation .
Urophysa Ulbr. (Ranunculaceae) is an Chinese endemic genus with two species, U. rockii
and U. henryi. These two species are morphologically distinct: the former has sacs near the
base of its petals while the latter has no sac (Fig 1). Additionally, U. rockii has a more restricted
distribution (only located in Jiangyou county, Sichuan Province) than U. henryi (distributed in
Yungui Plateau, south Chongqing, north Hunan and west Hubei). The two species' natural
populations are restricted to small and isolated areas separated by high mountains and deep
valleys, and grow in steep and karstic cliffs. The Yungui Plateau is bordered to the east by the
Mountains (Mts.) Daba, Wuling and Xuefeng, separating U. henryi Yungui Plateau
populations (populations in Guizhou Province) from neighboring populations (populations in
Hunan, Hubei and Chongqing Provinces). Moreover, these two species have strict habitat
requirements and are sensitive to environmental change and human activities (scenic spots,
power station). By field observations and laboratory experiments, we found that U. rockii and
U. henryi can not survive once leaving the karst limestone, and populations that located in
hydroelectric dams and tourist attractions possess less individuals than that have not been
2 / 24
Fig 1. The morphological characters of Urophysa rockii (a1±a4) and Urophysa henryi (b1±b4). 1:
Habitat feature; 2: Flowers; 3: Anatomical characteristics of petals; 4: Seed dispersal. The red arrows
represent spurred petals and saccate petals in a3 and b3, respectively.
effected by human activities. The plants are collected for Chinese traditional medicine for the
treatment of contusions and bruises. However, the population of the two species are in decline
]. There is ongoing research on the endangered U. rockii, including its growing
environment and conservation strategies , biological and ecological characteristics [
reproductive biology . However, there have been no studies on the evolutionary history
and distribution of these two species. Similarly, species differentiation and the effect of habitat
fragmentation in the Yungui Plateau have been little studied. Few studies that have been
conducted have found that the occurrence of physical barriers in the Yungui Plateau and its
adjacent regions have affect the genetic structure of species such as Eurycorymbus cavaleriei,
Saruma henryi and Dipentodon [30±32]. Therefore, we determined the phylogenetic
relationship between the two species, their lineage and cause of speciation and why they exhibit such
special distribution and distinct morphologies, particularly for the endangered U. rockii.
Phylogeographical approaches are useful for reconstructing evolutionary histories of species
and their close relatives [
], and are reliable for examining interspecific divergence and
]. Here, we used nine nuclear microsatellites (nSSR), two nuclear gene fragments (ITS
and ETS) and two cpDNA regions (psbA-trnH and trnL-trnF) to investigate the
phylogeography structure of U. rockii and U. henryi, and determine their divergence and evolution. It is
commonly believed that chloroplast genomes are susceptible to the influence of fragmentation
because they are maternally inherited, which can only be dispersed via seeds in sexual
reproduction. This is in contrast to nuclear regions, that are diploid and can be dispersed via seeds
or pollen. Therefore, we considered that the cpDNA markers were good indicators of
population historical dynamics and genetic effects [
]. Our specific aims were to: (1) Evaluate the
3 / 24
genetic diversity and population differentiation of U. rockii and U. henry. (2) Investigate
phylogeographical patterns of Urophysa and (3) Investigate the effect of the QTP uplift and climate
oscillations on the evolution of Urophysa, especially U. rockii. We believe that this study will be
useful for developing conservation strategies and restoration of these two endangered U. rockii
and endemic U. henry.
Materials and data
Plant sampling and morphological observation
A total of 190 individual plants from 14 populations were collected, covering almost the entire
geographic range of both species (Table 1), Nine to 16 individual plants were sampled from
the 14 populations in the Yungui Plateau (Guizhou Province) and its adjacent regions (i.e.,
Provinces of Hunan, Hubei, Chongqing and Sichuan Provinces). Healthy leaves were sampled
and dried in silica gel until total DNA was extracted. Voucher specimens were deposited in the
Herbarium of Sichuan University (SZ). The detailed population information of sampled
individuals was measured using a handheld GPS unit. Morphological characteristics were
identified using herbarium specimens and fresh materials. The mature pollen grains and leaves
dehydrated by graded ethanol were directly mounted on aluminum stubs using conducting
carbon adhesive tab, sputter-coated with gold, and then observed using the JSM-7500F
scanning electron microscope (SeM, Japan).
DNA extraction, sequencing and microsatellite genotyping
Total genomic DNA was isolated from silica-gel dried leaves using plant genomic DNA kit
(Tiangen Biotech, Beijing, China). Internal transcribed spacer (ITS) sequences were amplified
using the primers ITS4 and ITS5 [
] and ETS fragments were amplified based on the primers
ETS9bp and ETS18s [
]. The primers trnH and psbA [
]were used to amplify the psbA-trnH
sequences, and the trnL-trnF sequences were amplified by trnL and trnF [
chain reaction (PCR) applications were carried out in 30μL reaction volumes. For the ITS and
UR:Urophysa rockii; UH: Urophysa henryi; SC = Sichuan province; HN = Hunan province; HB = Hubei province; GZ = Guizhou province; CQ = Chongqing
4 / 24
ETS regions, reactions were conducted with the following program: an initial 4-min
denaturation at 94ÊC followed by 30 cycles of 45-sec denaturation at 94ÊC, 45-sec annealing at 56ÊC
for ITS (53ÊC for ETS) and 90-sec extension at 72ÊC with a final 5-min extension at 72ÊC [
For the psbA-trnH and trnL-trnF intergenic spacers, the PCR program began with 4-min
initial denaturing at 94ÊC followed by 30 cycles of 1-min denaturation at 94ÊC, 1-min annealing
at 54ÊC for psbA-trnH (45 sec at 64ÊC for trnL-trnF), and 2-min extension at 72ÊC. A final
extension was run for 7 min at 72ÊC for psbA-trnH (5 min at 72ÊC for trnL-trnF) [
PCR products were separated in 1.5% (w/v) agarose TAE gel and purified using Wizard PCR
preps DNA Purification System (Promega, Madison, WI, USA) following the manufacturer's
instructions. The purified PCR products were sequenced in an ABI Genetic Analyzer (Applied
Biosystems Inc., Foster City, CA, USA) in both directions using the PCR primers. The
chloroplast DNA sequences and nuclear fragments were edited and assembled using seqMan
], then aligned employing CLUSTAL W in MEGA 5 [
] and adjusted manually. All
haplotype sequences were deposited in the GenBank database under accession numbers
KR820593 to KR820702 (S1 Table).
According to the primers and amplification protocols developed for Li et al. [
selected nine pairs of SSR primers for our population genetic data analysis (S2 Table). A
preexperiment was performed using 20 pairs of EST primers, of which, nine pairs of primers that
showed significant polymorphism and high homology with relatives Aquilegia Li et al. [
were employed in following analyses (S3 Table). PCR products were separated on 3.5% of
agarose gel followed by staining with ethidium bromide. Alleles were sized using PeakScanner v.
1.0 software (Applied Biosystems).
Genetic diversity and divergence. The DnaSP version 5.0 [
] was used to identify
different haplotypes and to calculate haplotype diversity (Hd) and nucleotide diversity (π). To assess
the level of genetic variation and population differentiation, the average within populations
gene diversity (HS), total gene diversity (HT) and two population differentiation parameters,
] and NST [
] were estimated following the methods described by Pons and Petit [
using PERMUT [
]. The analysis molecular variance (AMOVA) in Arlequin version 3.5 [
was used to further quantify genetic differentiation between groups or subgroups, as well as
between populations within groups and among individuals within populations. In addition,
maximum parsimony median-joining method in NETWORK 18.104.22.168[
] was applied to
construct haplotype networks.
Fst is a measure of checking population differentiation and inferring the effects of gene flow
and drift [
]. Through constructing the regional scatterplots of the genetic differentiation
(quantified as Fst/(1 ±Fst)) on geographical distances and calculating the correlation
coefficients between them, we can evaluate the relative historical influences of gene flow and genetic
drift on regional population structure [
]. Based on this method, the Fst values between
pairwise populations were calculated using Arlequin version 3.5 [
]. Scatterplots of the Fst/(1 ±
Fst) on geographical distances were constructed, and correlation coefficients were calculated
along with the significance of correlation in GenAlEx version 6.0 [
], The linearized Fst
statistic [Fst/(1 ±Fst)] was compared with the matrix of geographical distance by means of a simple
Mantel test to detect isolation by distance and to evaluate the relative influences of gene flow
and genetic drift on the regional population structure. We used 999 random permutations to
test for the Mantel statistical significance.
For SSR analysis, the total number of detected alleles (NA), allelic richness (AR) and
inbreeding coefficient (FIS) were calculated using the software fstat 22.214.171.124 [
5 / 24
MICRO-CHECKER v. 2.2.3[
] was used to test for evidence of deviations from the Hardy±
Weinberg equilibrium, genotyping errors and null alleles, and compute the expected
heterozygosity (HE) and observed heterozygosity (HO). A cluster analysis was performed to characterize
population structure in STRUCTURE version 2.2 [
], along with 100 000 Markov chain
Monte Carlo (MCMC) cycles following 10 000 burn-in cycles, using the admixture model with
independent allele frequencies. Eleven replications were performed for each K, in the range
K = 2±10, and the optimal K value was determined according to Evanno et al. [
]. In addition,
the analysis molecular variance (AMOVA) in Arlequin version 3.5 [
] was used to measure
the partitioning of genetic variability within and among populations. In addition, the IBD
(isolation by distance) analysis was performed in GenAlEx version 6.0 [
Phylogenetic analysis. Maximum parsimony (MP) and Bayesian inference (BI) analysis
were conducted using the program PAUP version 4.0 beta 10 [
] and MrBayes v3.1.2 [
respectively. For parsimony analyses, heuristic searches were carried out with 1000 random
sequence replicates, with the tree bisection-reconnection (TBR) branch swapping and the
Mul-Trees options selected. All characters were weighted and unordered equally and gaps
were treated as missing data. Branch supporting values were estimated with 1000 replicates by
bootstrap analysis. For Bayesian inference analyses, the best-fit evolutionary model was
selected by the Akaike information criterion (AIC) using MrModeltest 2.2 [
]. The Bayesian
Markov chain Monte Carlo (MCMC) searches was performed for 2 × 10−7 generations with
four chains, sampling trees every 1000 generations and discarding the first 20% of sampled
trees as burn-in sample. The remaining trees were used for constructing a 50% majority-rule
consensus tree and calculating posterior probability (PP). In addition, maximum likelihood
(ML) analyses were performed by RAXML v7.2.8 [
] with 1000 bootstraps under the GTR
+ G substitution model, Nodes with a bootstrap value of 95% were considered
well-supported in this analysis. Additionally, we made a model test for each locus, and performed the
phylogenetic analysis again both in RAxML and MrBayes.
Demographic history analysis. Mismatch distribution analysis was undertaken to assess
whether genealogy experienced historical population expansions. It is assumed that if
populations experienced a sudden demographic expansion should display a unimodal and smooth
distribution, and the multimodal distributions are related to demographic equilibrium or
]. Gene flow between populations (NM) and Neutrality test including the Tajima's
D test [
] and Fu and Li's D and F statistics [
] were calculated by Arlequin version 3.5
Taking advantage of fossils from Prototinomiscium vangerowii (Menispermaceae) [
calibration points used in previous studies of Aquilegia [
], we estimated the divergence times
of Urophysa. Bayesian searches for tree topologies and node ages of cpDNA and nrDNA
dataset were conducted in BEAST using a GTR + G substitution model selected by JMODELTEST
] and an uncorrelated lognormal relaxed clock [
]. A Yule process was specified as tree
prior and the calibration priors was modelled as normal distributions with a mean time.
Unfortunately, no fossils are known for Urophysa or any closely related lineages. Therefore, we
chose the age of the fossil Prototinomiscium vangerowii (Menispermaceae). The
Menispermaceae are represented in the Cretaceous of Europe by endocarps assigned to the fossil genus
Prototinomiscium. Based on the oldest record, from the Turonian of Central Europe, dated to 91.0
], sets the minimum age of the split between Menispermaceae and Ranunculaceae
]. In addition, we employed the age interval 51±66 Mya reported by WikstroÈm et al. [
for the split of Ranunculus (subfamily Ranunculoideae) and Xanthorhiza (subfamily
Coptidoideae), the former genus being in the sister subfamily to Thalictroideae, which includes the
genus Urophysa [
]. Based on this information, we set priors of 91.0 ± 7.5 and 58.0 ± 2.5 Mya
respectively for each calibration point. A Yule speciation model was selected as the tree prior.
6 / 24
This is a simple model of speciation that is more appropriate when considering sequences
from different species. Considering Bastida et al. [
], we dated the origin of Thalictroideae as
a mean age of 27.61 Mya with a normal distribution and with a specified 95% confidence
interval (CI) of 26.59±28.56 Mya. The stem node of Aquilegia was dated to 10.18 Mya with a normal
distribution and with the specified 95% confidence interval (CI) of 9.21±11.14 Mya [
program BEAST version 1.8.0 [
] was employed to estimate the divergence of major lineages
of Urophysa, MCMC runs were performed, each of 2 × 107 generations with sampling every
2000 generations, following a burn-in of the initial 10% cycles. MCMC samples were inspected
in Tracer to confirm sampling adequacy and convergence of the chains to a stationary
distribution. All of the outgroup sequences used for time calibration were listed in S4 Table.
The substitution rates were also used to estimate time of divergence of Urophysa. Aquilegia
incurvata and Semiaquilegia adoxoides were chosen as outgroups [
]. Considering the
dominant geographical distribution of cpDNA haplotypes and poorly resolved relationships
among nrDNA haplotypes, the cpDNA dataset was used to estimate the divergence time. The
best-fit model was GTR + G for cpDNA data inferred from the Akaike information criterion
(AIC) with MrModeltest 2.2 . An uncorrelated lognormal model was selected to describe
the relaxed clock. The constant cpDNA substitution rates for most angiosperm species have
been estimated to be in the range 1±3 × 10−9 s/s/y [
]. Given Urophysa are perennial plants,
the general substitution rates of the plastid sequence (u = 1.52 × 10−9 s/s/y) was more suitable.
The Markov chain Monte Carlo (MCMC) analyses were run for 2 × 107 generations and trees
were sampled every 2000 generations. The first 10% of sampled trees was discarded as the
burn-in sample as checked with Tracer 1.5 [
]. The Tree Annotator version1.4.8 [
used to summarize the samples in the maximum clade credibility tree with the posterior
probability limit set to 0.5. The results were displayed in Figtree version1.3.1 [
The scanning electron micrographs of leaf epidermis and pollen are shown in S1 Fig. 1) U.
henryi had spindly and numerous epidermal hairs (a1-c1, b2-c2 and a3-c3), while U. rockii had
few hairs that were short and had swollen-bases (d1-d3). 2) Distinctive papillary (a2 and c2)
and stelliform (b2) surface ornamentation patterns of leaves were found in U. henryi, distinctly
different to the sinuous surface of U. rockii (d2). 3) The leaf epidermis in U. henryi was
ornamented by bar-shaped appendages (i.e., b5 and c5 except a5, which was smooth), which was
noticeably different from U. rockii with lineate appendages (d5). 4) The stomatas of U. henryi
can be divided into sunken (a5) and flat types (b5 and c5), but the stomatas of U. rockii were
raised (d5). In the polar view of pollen grains, it is slightly pointed in U. henryi (e1) but flat in
U. rockii (f1), and the pollen grains pores were sparse in U. henryi (e2) but numerous in U.
cpDNA sequence analysis
The two cpDNA regions were aligned along a total length of 1567 bp (psbA-trnH, 572 bp;
trnL-trnF, 995 bp) and 17 chlorotypes were generated, including 71 polymorphisms (S5
Table). Chlorotypes H1±H14 were found only in U. henryi and H15±H17 were found in U.
rockii (Fig 2 and Table 2). Only three chlorotypes (H1, H3, H15) were shared by two or more
populations and no chlorotype was shared by these two species. In addition, we detected no
haplotype located in the center of the chlorotype network, while many private cpDNA
haplotypes were detected in each population. At the species level, for U. rockii, Hd = 0.742 and π =
0.00117, for U. henryi, Hd = 0.917 and π = 0.01398. Population JY2 of U. rockii and SM of U.
7 / 24
Fig 2. Geographical distribution of cpDNA and NETWORK-derived genealogical relationship. (a) Distribution ranges of U. henryi and U.
rockii. (b) The distributions of cpDNA haplotypes (H1±H17) detected in U. rockii and U. henryi (population codes see Table 1). (c) The statistical
parsimony network of cpDNA haplotypes (H1±H17). The topographic map is from the website https://en.wikipedia.org/wiki/File:East_Asia_
henryi possessed the highest haplotype diversity and the maximum number of haplotypes. The
total genetic diversity (HT) was higher than the diversity within populations (HS) for both
species, and HT and HS were higher in U. henryi than in U. rockii. Additionally, NST was
significantly higher than GST (NST = 0.921, GST = 0.716, P < 0.01) but only in U. henryi (Table 3),
indicating a significant phylogeographical structure existing between populations of U. henryi.
AMOVA indicated that 41.32% of the total cpDNA variation was partitioned between
species and 57.06% of the variation could be attributed to variation between populations within
species (Table 4). For each species, statistically significant variation was detected between
populations (88.24% for U. rockii and 97.30% for U. henryi). The observed value in the mismatch
distribution analysis fitted multimodal curves (S2A±S2E Fig) and the results of Tajima's D and
Fu and Li's D and F statistics were not significantly negative (S6 Table), which indicated that
the populations did not undergo expansion. Mantel test results indicated a significant effect of
isolation by distance (IBD) at the genus range scale (r = 0.322, p = 0.001), as well as the species
range of U. henryi (r = 0.246, p = 0.039), while negative in U. rockii (r = ±0.150, p = 0.561)
The topological structure of haplotypes derived from Parsimony analyses was similar to
that from Bayesian tree and maximum likelihood tree analyses, therefore only the Bayesian
tree with maximum parsimony and maximum likelihood bootstrap support values was shown
in S4 Fig. In the Bayesian tree, three major clades (Clade I, II and III) corresponded to
geographical distributions of the populations. Clade I consisted of haplotypes H7±H13 that were
all from the populations of Yungui Plateau, Clade II included the haplotypes H1±H6 mainly
located in adjacent regions of Yungui Plateau (Hunan and Hubei Provinces) H14±H17
haplotypes formed the clade III, in which H15±H17 of U. rockii situated at the edge of Sichuan
Basin, were at the crown of phylogenetic tree. The Median-joining network of cpDNA
8 / 24
* Exclusive haplotype or sequence types
Hd, Haplotype diversity; π, Nucleotide diversity.
H30*(11), H31*(2), H32*(2)
H1*(2), H2*(2), H3*(3), H4*(2), H5*(2), H6*(2), H7*(1)
H10*(2), H11*(1), H12*(2), H13*(1), H14*(2), H15*(3), H16*(1), H17*(2), H18*
Average gene diversity within populations (HS), total gene diversity (HT), interpopulation differentiation (GST), and number of substitution types (NST).
9 / 24
FCT = differentiation among groups; FST = differentiation among populations; FSC = differentiation among populations within groups.
** P < 0.001, 1000 permutations.
haplotypes was consistent with the strict consensus tree produced by Bayesian analysis (Fig
2C). In addition, after testing the substitution model of each locus and combining to perform
phylogenetic analysis, we found that the topology of RAxML and MrBayes trees were similar,
and was consistent with topology which was obtained from cpDNA haplotypes (S7 Table and
From two time-calculated ways (calibration strategy and substitution rate), we estimated
the divergence time of Urophysa, and found that the calibrated time estimated by cpDNA data
is more accurate and reliable compared with time obtained by substitution rate. The origin
time of Urophysa was dated to approximately 10.29 Mya (95% HPD = 8.99±11.69 Mya) (Fig 3).
10 / 24
Fig 3. BEAST-derived chronograms for Urophysa cpDNA haplotypes. Numbers on the branches indicate the Bayesian posterior probabilities values.
Ages of the main clades are shown below the branches. Haplotypes of U. henryi are H1±H14 and H15±H17 belong to U. rockii.
U. henryi began to diverge at approximately 8.86 Mya (95% HPD = 6.24±11.00 Mya). However,
the origin of U. rockii did not occur until the early Quaternary, which was dated to
approximately 3.16 Mya (95% HPD = 0.78±6.33 Mya), and it diverged at the early Pleistocene,
approximately 1.48 Mya (95% HPD = 0.28±3.84 Mya). We obtained a roughly similar time from the
substitution rate method, the divergence time of U. henryi was dated to approximately 7.86
Mya (95% HPD = 4.29±13.07 Mya) and the divergence of U. rockii was dated to 1.62 Mya
(95% HPD = 0.42±3.48 Mya) (S5 Fig).
nrDNA sequence analysis
The length of aligned sequences was 665 bp for ITS, and 521 bp for ETS. Based on 69 nuclear
polymorphic sites combined the two data sets, 34 haplotypes (H1±H34) were detected (Fig 4A
and S8 Table). All of the haplotypes were unique to a population except H33 (Table 2), which
was shared by population JY3 and JY4. The haplotype diversity (Hd) ranged from 0 to 0.378
and the nucleotide diversity (π) from 0 to 0.00051 within U. rockii. The haplotype diversity
(Hd) within U. henryi ranged from 0 to 0.978 and the nucleotide diversity (π) from 0 to
0.00765. Diversity at the species level for U. rockii was Hd = 0.756 and π = 0.00287, and for U.
henryi was Hd = 0.943 and π = 0.00930. The highest haplotype diversity and the maximum
number of haplotypes were observed in population SM of U. henryi and JY2 of U. rockii.
AMOVA showed that 62.14% of the genetic variation occurred between species and 31.98%
of the variation was partitioned among all populations within species. Within both species,
11 / 24
Fig 4. The distributions of nrDNA genotype and Network relationships. (a) Geographical distribution of nrDNA haplotypes (H1±H34) detected
in U. rockii and U. henryi (population codes see Table 1). (b) Median-joining network of nrDNA haplotypes (H1±H34). Sizes of the circles in the
network are proportional to the observed frequencies of the haplotypes. The small black bars represent mutation steps and the black dots (MV)
represent missing haplotypes. The topographic map is from the website https://en.wikipedia.org/wiki/File:East_Asia_topographic_map.png.
variation between populations was significant and accounted for 82.06% and 53.11% of total
variation in U. henryi and U. rockii, respectively (Table 4). Similar to the findings for cpDNA,
mismatch distribution analysis was clearly multimodal (S2F±S2H Fig), as well as a positive and
non-significant result of neutrality test (p > 0.10) (S6 Table), which indicated the populations
did not undergo expansion.
A high total genetic diversity (HT = 0.989) and low average within-population genetic
diversity (HS = 0.397) was detected in U. henryi and in U. rockii (HT = 0.900, HS = 0.116).
Similar to cpDNA, the results of nrDNA (NST > GST, P < 0.01, Table 3) indicated a
significant phylogeographical structure existed between populations of U. henryi. Mantel test
results indicated a significant effect of isolation by distance (IBD) at the genus range scale
(r = 0.301, p = 0.004), whereas no significant IBD pattern presented in the two species
Only the Bayesian tree with parsimony bootstrap and maximum likelihood support values
was shown in S6 Fig because it has the same topology as the maximum parsimony tree and
maximum likelihood tree. U. henryi and U. rockii were clustered into a single clade with high
bootstrap support and posterior probability values (S6 Fig). The network of nrDNA haplotypes
was roughly consistent with the phylogenetic tree, but no ancestral haplotypes were detected
Divergence time estimation by nrDNA data (ITS) indicated that Urophysa began to diverge
at approximately 9.65 Mya (95% HPD = 7.48±11.48 Mya) (S7 Fig), which was slightly older
than that by cpDNA data. Estimations of divergence times of species within Urophysa are not
reliable because of low branches posterior probabilities in the nrDNA tree, and it is better to
use the time estimated by cpDNA data.
12 / 24
NA: Observed number of alleles; AR: Allelic richness; HO: observed heterozygosity; HE: expected heterozygosity; Y/N: yes/no; FIS: inbreeding coefficient at
the population level
*significant at P < 0.05; NM: Gene flow estimated from Nm = 0.25*(1-FST)/FST.
Genetic differentiation at the microsatellite loci
A few null alleles were detected by Micro-Checker v.2.2.3 test (Table 5). To test the impact of
these null alleles, we removed the microsatellite markers (A41, EST2, EST9) with null alleles in
both species and used the remaining six markers to perform the population structure analysis.
We found that excluding microsatellite markers with null alleles had little impact on our
results (S8 Fig), and U. rockii and U. henryi are distinctly separate. Therefore, all microsatellite
markers were used for relevant analyses. The results of SSRn indicated low mean allelic
richness and gene diversity in U. rockii (AR = 3.386; HO = 0.165; HE = 0.341) and U. henryi (AR =
5.464; HO = 0.249; HE = 0.438; see Table 5). All pairwise population genetic differentiations
were significant (p < 0.001), except between adjacent U. rockii populations JY3 and JY4
(pairwise FST = 0.054, p > 0.05). Pairwise FST estimates are listed in S9 Table. FIS values ranged
from 0.103 to 0.930 across the U. rockii with an average value of 0.476. FIS values ranged from
±0.148 to 1.000 for U. henryi, with an average value of 0.384.
The STUCTURE analysis, using the ΔK method, showed that the optimal K value was K = 2
(Fig 5A), which strongly supported two genetic clusters among our samples and generally
corresponded to the two respective species. When K = 5, we detected a small peak, which further
divided populations of U. henryi into four groups (Fig 5). AMOVA indicated that 31.7% of the
genetic variation occurred between species (FCT = 0.317, p < 0.001). Within each species, most
of the genetic variation was partitioned among populations (Table 4). A significant effect of
isolation by distance (IBD) was detected between the 14 populations of Urophysa (r = 0.646,
13 / 24
Fig 5. Bayesian clustering analysis for the species of Urophysa inferred by STRUCTURE. (a) Bayesian
inference of the cluster number (K). (b) Results for clusters (K = 2 and 5) as detected by STRUCTURE. K was
estimated using the distribution of ΔK (second order rate of change of the likelihood distribution). The bars on
the figure represent these individuals were sampled from the same species. Bar plots showing Bayesian
assignment probabilities. Each vertical bar corresponds to one individual. Populations are separated by black
bars and identified at the bottom.
p = 0.001) and between populations of U. henryi (r = 0.356, p = 0.020), while negative in U.
rockii (r = ±0.175, p = 0.484) (S3G±S3I Fig).
Genetic diversity and significant population differentiation
Based on cpDNA and nrDNA data sets, high haplotype diversity (Hd) and nucleotide diversity
(π) were observed in U. henryi and U. rockii at the species level (Table 2). We detected a high
level of total genetic diversity in U. henryi (HT = 0.959 for cpDNA, HT = 0.989 for nrDNA) and
in U. rockii (HT = 0.891 for cpDNA, HT = 0.900 for nrDNA), which was higher than other
14 / 24
species of Ranunculaceae, such as Aconitum gymnandrum (HT = 0.739 for cpDNA) and
Clematis sibirica (HT = 0.496 for nrDNA) [
]. Genetic diversity of species is closely related to
environmental factors (e.g., climate, topography) and life history characteristics (i.e., life cycle,
breeding system). U. henryi and U. rockii possess a relatively long life cycle and sexual
reproductive period and high offspring mortality, which may lead to low levels of genetic diversity [
It has been found that habitat fragmentation can reduce genetic diversity due to restricted gene
flow, genetic erosion and random genetic drift in isolated populations [
15, 19, 75, 76
other studies have shown that genetic variation could be maintained or even increased in
fragmented populations [
]. Individuals of U. henryi and U. rockii are confined to steep and
karstic limestones in ravines, which are naturally fragmented. These conditions are found
especially in the Yungui Plateau, as high mountains, deep valleys and fast-flowing rivers acting as
natural barriers enhancing isolation, drift and mutation [
]. Isolation of U. henryi and U. rockii
individuals was also supported by the IBD model and Neutrality test (S3 Fig, S6 and S9 Tables),
indicating a significant isolation-by-distance pattern between populations of these two species
and populations of U. henryi. Thus, the high genetic diversity observed in these two species may
be closely related to their fragmented habitats.
Despite a high level of cpDNA genetic diversity detected at the species level, the average
genetic diversity within-population was low (HS). We found a high coefficient of genetic
differentiation (GST) in U. henryi and in U. rockii, as well as significant phylogeographical structure (NST
> GST, P < 0.01) and variation (FSTcpDNA = 0.973, FSTnrDNA = 0.821 for U. henryi; FSTcpDNA =
0.882, FSTnrDNA = 0.531 for U. rockii) between populations. Additionally, a low level of genetic
variation was indicated by microsatellite markers (mean HO = 0.249 for U. henryi and HO = 0.165
for U. rockii) at the population level. While a high degree of genetic differentiation in
microsatellites was found in U. henryi (FST = 0.663) and in U. rockii (FST = 0.669) (Table 4). These results
indicate high genetic differentiation exists between populations. Akin to our study, the
combination of low genetic diversity and high genetic differentiation has been reported for other plants in
Yungui Plateau and adjacent regions, such as Cercidiphyllum japonicum, Tetrastigma
hemsleyanum and Cardiocrinum giganteum [77±79].
Limited gene flow may be a crucial factor resulting in low genetic diversity and high genetic
differentiation. It is believed that the breeding system of plant (i.e., flowering and seed
proliferation) is an important biological characteristic, which strongly influences the spatial-temporal
distribution of genetic variation [
] and affects processes that lead to speciation or extinction
]. Although the reproduction system and pollination mechanism of U. henryi has been little
studied, previous research of U. rockii reproduction has indicated the seeds are extremely
small and are dispersed into rock gaps due to mechanical strain when the follicle cracks
spontaneously after seeds have matured [
]. However, most seeds are washed away by rainwater,
and only a few seeds fall to the base of cliff where they may germinate. Petit et al. [
that species with seeds dispersed by gravity tended to show higher differentiation between
populations than species with wind-dispersed seeds. Interestingly, U. rockii and U. henryi are
entomophilous plants and bloom in the winter when lower temperatures reduce the visiting
frequency of pollinators, which may result in weak pollen flow. The poor gene flow mediated
by seeds and pollen was also supported by our results (NM) (Table 5 and S6 Table). Thus,
limited gene flow in such fragmented habitats could lead to significant genetic differentiation
between populations [
Human activities have significantly influenced the distribution of these two species. Over
the past few decades, many hydroelectric dams and tourist attractions have been built in areas
where the wild populations of U. henryi and U. rockii are located leading to the extensive
removal of their habitat. The number of individuals has also sharply declined due to excessive
15 / 24
collection for their medicinal value. All of these human activities have reduced the population
size, increased fragmentation and isolation, and enhanced population differentiation.
Demographic history of U. henryi
The time of origin obtained from the molecular clock estimation is generally congruent with
that from points calibration. The molecular clock estimated that the populations of U. henryi
in the Yungui Plateau (Clade I) diverged approximately at 8.86 (6.24±11.0) Mya, which is in
accordance with the rapid uplift time of the QTP [
]. It is believed that the QTP reached an
elevation similar to the present at about 8.0 Mya, but decreased following extensive faulting.
The most recent rapid uplift of the OTP occurred at around 3.6 Mya [
]. Given the geologic
close relationship between the Yungui Plateau and the QTP, the orogenic events of the Yungui
Plateau were similarly violent during the continuous uplift of the QTP [
]. Therefore, we
suggest that rapid QTP uplift and subsequent Yungui Plateau contortions significantly
contributed to the differentiation of U. henryi. Besides, global climate also fluctuated dramatically at
that time and the prevailed East Asian monsoon brought plenty of rainfall [
In the cpDNA tree, two clades of U. henryi were revealed (Clade I and clade II) (Fig 3 and
S4 Fig), which were consistent with geographical distribution of U. henryi (Fig 2). This
grouping was reflected by a striking differentiation between populations. U. henryi populations were
further divided into four groups when K = 5 in structure analysis, and this roughly
corresponded with geographical regions. Mountains such as Wuling, Dalou and Xuefeng. extend
from northeast to southwest with an average altitude of 2,000 m [
], and acted as geographical
barriers between Yungui Plateau (Clade I) and its adjacent populations (clade II) of U. henryi.
Previous phylogeographical studies have identified that Mt. Wuling and Mt. Xuefeng as the
major barriers to gene flow [
]. In addition, the significant increase in geological and
ecological variabilities during rapid uplift of the Yungui Plateau has promoted rapid divergence in
other small and isolated populations such as Ligularia±Cremanthodium±P arasenecio complex,
Babina pleuraden, Dipentodon and Eurycorymbus cavaleriei [
30, 31, 89, 90
]. Therefore, the
habitat of U. henryi was fragmented due to geographical barriers and fluctuating climate
conditions (i.e., unstable rainfall), and populations of U. henryi were separated causing high
The allopatric divergence of U. rockii
To investigate the phylogenetic relationship between U. rockii and U. henryi, we identified
distinct morphological traits of petals, leaf epidermis and pollen grains (Fig 1 and S1 Fig). Results
of SSRn and nrDNA (Figs 4 and 5, S6 Fig) showed that the U. rockii and U. henryi are distinctly
separate species and the genus Urophysa was a monophyly, which were consistent with
previous research [
]. The cpDNA phylogeny demonstrated that U. rockii was at the end of the
cpDNA phylogenetic tree with high bootstrap support values and exhibited significant
differentiation with U. henryi.
Habitat fragmentation is likely to significantly influence the divergence and allopatric
speciation of plants [
]. We estimated the origin time of U. rockii was 3.16 Mya. At that time,
orogeny of the Yungui Plateau triggered habitat fragmentation and significant geographic isolation
between populations. U. rockii grows exclusively on cliffs or fissures of rocks in China and
only a few populations have been found. Connectivity with U. henryi populations was hindered
by geographic barriers, the Sichuan Basin and Yangzi River, and thus there was no gene flow
(Fig 2 and S6 Table). The effect of barriers on gene flow in this region has been documented
for Myotis pilosus [
] and Tapiscia sinensis [
]. The significant isolation caused U. rockii to
undergo allopatric divergence. Additionally, We found considerable genetic differentiation
16 / 24
and limited gene flow between U. rockii and U. henryi (Table 5 and S6 Table). This high genetic
differentiation between populations is likely due to isolated populations with restricted gene
]. This is also supported by the Mantel test between the pairwise FST/(1- FST) and
geographic distance (S3 Fig). Moreover, local adaptation may have played an important role in
driving allopatric speciation of U. rockii. As expected, we identified extensive footprints of
local adaptation from its specialized morphologies. These specialized morphologies included
unusual floral organs such as petaloid sepals that can display various colors in different flower
phases, the petals with a nectar spur that could attract more pollinators, and a mass of small
seeds (thousand seed weight is 0.6684 ± 0.0038g) [
]. These specialized morphologies likely
contributed enhanced reproductive efficiency.
The divergence of the U. rockii was dated to approximately 1.48 Mya, coinciding with the
frequent climatic oscillation during the Pleistocene, which is considered as one of the most
important periods for genetic diversification [
]. It is believed that dramatic climate
fluctuation during the Pleistocene resulted in massive changes in plant population distributions, with
many species shifting to more suitable habitat and other spatial-temporal adjustments [
Mismatch analysis indicated that expansion did not occur in populations of U. rockii (S2E and
S2H Fig), and that U. rockii populations were presumably impacted by the climatic oscillation
of the Pleistocene. However, the mountain system located in north of Sichuan Basin was in
favor of preserving plant species [
] The mountains of Qinling, Daba and Micang could also
have acted as barriers to reduce the effects of the Pleistocene climate [
]. These mountains
have also been regarded as a key glacial refuge for other plants including Pinus massoniana,
Liriodendron chinense and Rhinolophus ferrumequinum [
79, 95, 97
]. Overall, we believe that
long-term geographical isolation, limited gene flow caused by habitat fragmentation and
specialized morphologies probably contributed to the allopatric divergence of U. rockii. The
climatic oscillation of the Pleistocene further promoted population divergence and resulted in
the current distribution of U. rockii.
Implications for conservation
U. rockii is an endangered species of China that has a small geographic range, few distinct
populations and extremely low germination rate (less than 2%) [
], and each population has a
limited number of individuals (all five populations >2,000) [
]. Genetic drift is likely to have
occurred, which could have led to the observed genetic diversity decline and increased
population differentiation. The correlation estimation of the genetic differentiation [Fst/(1 ±Fst)] and
geographical distances (isolation by distance) showed that genetic drift was much more
influential than gene flow on the distribution of genetic variability. Consequently, populations of U.
rockii are not at equilibrium, which may have resulted from its strict habitat requirement or
narrow distribution. A high value of FST indicates high genetic drift load, while a low value of
HS signifies high inbreeding [
]. In our study, genetic drift and inbreeding were both high
(Table 3 and Table 4), which was mirrored by the FIS from SSR data (Table 5), indicating a risk
Although U. rockii and U. henryi demonstrated high genetic diversity at the species level,
high genetic drift and inbreeding can be deleterious to species survival. Similarly, species with
small population sizes and fragmented distributions are vulnerable to extinction especially
with high genetic drift and inbreeding [
]. Protection of in situ populations and their
habitat, including removing threats (i.e., human disturbance) is essential for U. rockii’s survival in
the wild. Ex situ breeding programs can also be implemented using seeds of local provenance
to propagate different genotypes to ensure diversity and provide a source of seedling for
planting in the wild. Although U. henryi is not (yet) endangered, frequent human activities have
17 / 24
fragmented and degraded its habitat. It is possible that U. henryi will be threatened in the
future due to persistent human disturbance, restricted range and isolated populations. It is
necessary to adopt clear and practicable measures to restrict anthropogenic disturbances now,
as genetic diversity of wild populations can be maintained and conserved, rather than
implement measures after future declines.
S1 Fig. Scanning electron micrographs of leaf epidermis and pollen grains features. Leaf
epidermis: Urophysa henryi are shwed in a±c, Urophysa rockii are showed in d. In a±d, the
number after each letter indicate: 1±2: upper epidermis; 3±5: lower epidermis. Pollen grains:
U. henryi are showed in e1±e2, U. rockii are showed in f1±f2.).
S2 Fig. Mismacth distribution analysis for chloroplast DNA haplotypes (a±e) and nrDNA
haplotypes (f±h): (a) Urophysa; (b) U. henryi; (c) Clade I; (d) Clade II; (e) U. rockii (Clade III);
(f) Urophysa; (g) U. henryi; (h) U. rockii. The solid line represents expected (Exp) values under
a sudden population expansion, the dashed line shows observed (Obs) values.
S3 Fig. Scatterplots representing relationships between Plots of genetic distance [Fst/(1 ±
Fst)] and geographic distance (Km). At genus and species levels based on cpDNA (a±c),
nrDNA (d±f) and SSR (g±i) data. (a, d, g) Urophysa; (b, e, h) U. henryi; (d, f, i) U. rockii.
S4 Fig. Phylogenetic relationships based on the 17 cpDNA haplotypes. Numbers on the
branches indicate the maximum Parsimony bootstrap, maximum likelihood support value
and Bayesian posterior probabilities, respectively.
S5 Fig. BEAST-derived chronograms for Urophysa cpDNA haplotype based on constant
cpDNA substitution rate (u = 1.52 × 10−9 s/s/y). Numbers on the branches indicate the
Bayesian posterior probabilities values. Ages of the main clades are shown below the branches.
Different colors represent different haplotypes of species: blue, the haplotypes of U. rockii; red,
the haplotypes of U. henryi.
S6 Fig. Phylogenetic relationships of the 34 nrDNA haplotypes. Numbers on the branches
indicate the maximum Parsimony bootstrap, maximum likelihood support value and Bayesian
posterior probabilities, respectively.
S7 Fig. Divergence time calibration for Urophysa nrDNA haplotypes. Numbers on the
branches indicate the Bayesian posterior probabilities values. Ages of the main clades are
shown below the branches. Haplotypes H1±H27 and H28±H34 are possessed by U. henryi and
U. rockii, respectively.
S8 Fig. Histogram of the STRUCTURE analysis for the six microsatellite markers that did
not present null alleles. (a) Bayesian inference of the cluster number (K). (b) Results for
clusters (K = 2) as detected by STRUCTURE. The bars on the figure represent these individuals
that were retrieved from the same species. Bar plots showing Bayesian assignment
probabilities. Each vertical bar corresponds to one individual. Populations are separated by black bars
18 / 24
and identified at the bottom.
S1 Table. Information of our haplotype sequences deposited in the GenBank. Numbers in
the brackets is the sequence number of each haplotype included. : represent the specific
haplotype of each population.
S2 Table. Microsatellite markers used in this study. For each primer pair, forward (F) and
reverse (R) primer sequence, repeat motif (Repeat), size of cloned allele (bp), optimal PCR
annealing temperature (Ta).
S3 Table. Homology test of each SSR locus between Urophysa and Aquilegia.
S4 Table. psbA-trnH and trnl-trnF sequences of outgroups from Genbank for time
S5 Table. Variable sites of the aligned two chloroplast DNA fragments (psbA-trnH and
trnL-trnF) among Urophysa.
S6 Table. Parameters of Neutrality test and gene flow among populations. NM: gene flow.
S7 Table. Model test for each of locus and phylogenetic analysis. Numbers on the branches
indicate the maximum likelihood support value and Bayesian posterior probabilities,
respectively. Different colors represent different populations of species: red, the populations of U.
henryi; blue, the populations of U. rockii.
S8 Table. Variable sites of the aligned ITS and ETS fragments among Urophysa.
S9 Table. Pairwise FST values among the 14 populations of U. henryi and U. rockii based on
SSRn data. Note: Values in bold were not significantly different from zero after sequential
Bonferroni correction. : significant at p < 0.001.
We are grateful to De-Qing Huang, Yun-Dong Gao and Xiang-Guang Ma for their guidance
in paper revision and to Mr. Zheng-Yu Liu from the Chongqing Institute Of Medicinal Plant
Cultivation for providing the materials information of the CQ population. This work was
supported by the National Natural Science Foundation of China (Grant Nos. 31270241, 31470009,
31570198, 31500188), and the Specimen Platform of China, Teaching Specimen's sub-platform
(Available website: http://mnh.scu.edu.cn/), the Science and Technology Basic Work (Grant
Conceptualization: Deng-Feng Xie, Yan Yu, Xing-Jin He.
19 / 24
Data curation: Deng-Feng Xie, Min-Jie Li, Jin-Bo Tan, Yan Yu.
Formal analysis: Deng-Feng Xie, Megan Price, Qun-Ying Xiao.
Funding acquisition: Yan Yu, Xing-Jin He.
Investigation: Deng-Feng Xie, Jin-Bo Tan, Xing-Jin He.
Methodology: Deng-Feng Xie, Song-Dong Zhou, Yan Yu.
Project administration: Song-Dong Zhou, Yan Yu, Xing-Jin He.
Software: Deng-Feng Xie, Min-Jie Li, Jin-Bo Tan.
Supervision: Deng-Feng Xie, Xing-Jin He.
Writing ± original draft: Deng-Feng Xie, Min-Jie Li, Jin-Bo Tan, Qun-Ying Xiao.
Writing ± review & editing: Deng-Feng Xie, Megan Price, Qun-Ying Xiao, Yan Yu.
20 / 24
21 / 24
22 / 24
23 / 24
1. Young A , Boyle T , Brown T. The population genetic consequences of habitat fragmentation for plants . Trends in Ecology & Evolution . 1996 ; 11 ( 10 ): 413 ± 418 .
2. Hewitt G. The genetic legacy of the Quaternary ice ages . Nature . 2000 ; 405 ( 6789 ): 907 ± 913 . https:// doi.org/10.1038/35016000 PMID: 10879524
3. Cheng J , Liu XQ , Gao ZJ , Tang DY , Yue JW . Effect of the Tibetan Plateau uplifting on the geological environment of the Yunnan Plateau . Geoscience. 2001 ; 15 ( 3 ): 290 ± 296 .
Wang LY , Abbott RJ , Zheng W , Chen P , Wang YJ , Liu JQ . History and evolution of alpine plants endemic to the Qinghai-Tibetan Plateau: Aconitum gymnandrum (Ranunculaceae) . Molecular Ecology . 2009b ; 18 ( 4 ): 709 ± 721 .
5. Liu J , MoÈller M , Provan J , Gao LM , Poudel RC , Li DZ . Geological and ecological factors drive cryptic speciation of yews in a biodiversity hotspot . New Phytologist . 2013 ; 199 ( 4 ): 1093 ± 1108 . https://doi.org/ 10.1111/nph.12336 PMID: 23718262
6. Zhang JQ , Meng SY , Rao GY . Phylogeography of Rhodiola kirilowii (Crassulaceae): A Story of Miocene Divergence and Quaternary Expansion . PloS One . 2014 ; 9 ( 11 ):e112923. https://doi.org/10. 1371/journal.pone. 0112923 PMID: 25389750
7. Meng LH , Yang R , Abbott RJ , Miehe G , Hu TH , Liu JQ . Mitochondrial and chloroplast phylogeography of Picea crassifolia (Pinaceae) in the Qinghai-Tibetan Plateau and adjacent highlands . Molecular Ecology . 2007 ; 16 ( 19 ): 4128 ± 4137 . https://doi.org/10.1111/j. 1365 - 294X . 2007 . 03459 . x PMID : 17894761
8. Li Y , Stocks M , HemmilaÈ S , KaÈllman T , Zhu H , Zhou Y , et al. Demographic histories of four spruce (Picea) species of the Qinghai-Tibetan Plateau and neighboring areas inferred from multiple nuclear loci . Molecular Biology and Evolution . 2010 ; 27 ( 5 ): 1001 ± 1014 . https://doi.org/10.1093/molbev/ msp301 PMID: 20031927
Wu FY , Huang BC , Ye K , Fang AM . Collapsed Himalayan-Tibetan orogen and the rising . Tibetan Plateau Acta Petrologica Sinica . 2008 ; 24 :1± 30 .
10. Kang M , Tao J , Wang J , Ren C , Qi Q , Xiang QY , Huang H . Adaptive and nonadaptive genome size evolution in Karst endemic flora of China . New Phytologist. 2014 ; 202 ( 4 ): 1371 ± 1381 . https://doi.org/ 10.1111/nph.12726 PMID: 24533910
11. Wang WT . On some distribution patterns and some migration routes found in the eastern Asiatic region . Acta Phytotaxonomica Sinica . 1992 ; 30 :1± 24 .
12. Li XW , Li J. A preliminary floristic study on the seed plants from the region of Hengduan Mountain . Acta Botanica Yunnanica . 1993 ; 15 ( 3 ): 217 ± 231 .
13. Zhang DF , Fengqu L , Jianmin B . Eco-environmental effects of the Qinghai-Tibet Plateau uplift during the Quaternary in China . Environmental Geology . 2000 ; 39 ( 12 ): 1352 ± 1358 .
14. Sun H , Li ZM . Qinghai-Tibet Plateau uplift and its impact on Tethys flora . Advances In Earth Sciences . 2003 ; 18 : 852 ± 862 .
15. Aguilar R , Quesada M , Ashworth L , Herrerias-Diego Y , Lobo J. Genetic consequences of habitat fragmentation in plant populations: susceptible signals in plant traits and methodological approaches . Molecular Ecology . 2008 ; 17 ( 24 ): 5177 ± 5188 . https://doi.org/10.1111/j. 1365 - 294X . 2008 . 03971 . x PMID : 19120995
16. Honnay O , Jacquemyn H . Susceptibility of common and plant species to the genetic consequences of habitat fragmentation . Conservation Biology . 2007 ; 21 ( 3 ): 823 ± 831 . https://doi.org/10.1111/j.1523- 1739 . 2006 . 00646 . x PMID : 17531059
17. Lienert J . Habitat fragmentation effects on fitness of plant populations-a review . Journal for Nature Conservation . 2004 ; 12 ( 1 ): 53 ± 72 .
18. Sork VL , Smouse PE . Genetic analysis of landscape connectivity in tree populations . Landscape Ecology . 2006 ; 21 ( 6 ): 821 ± 836 .
19. Heinken T , Weber E. Consequences of habitat fragmentation for plant species: Do we know enough? Perspectives in Plant Ecology , Evolution and Systematics . 2013 ; 15 ( 4 ): 205 ± 216 .
20. Dubreuil M , Riba M , GonzaÂlez-MartÂõnez SC , Vendramin GG , Sebastiani F , Mayol M. Genetic effects of chronic habitat fragmentation revisited: strong genetic structure in a temperate tree, Taxus baccata (Taxaceae), with great dispersal capability . American Journal of Botany . 2010 ; 97 ( 2 ): 303 ± 310 . https:// doi.org/10.3732/ajb.0900148 PMID: 21622391
21. Ellstrand NC , Elam DR . Population genetic consequences of small population size: implications for plant conservation . Annual Review of Ecology and Systematics . 1993 ; 24 ( 1 ): 217 ± 242 .
22. Keller LF , Waller DM . Inbreeding effects in wild populations . Trends in Ecology & Evolution . 2002 ; 17 ( 5 ): 230 ± 241 .
23. Carson HL . Increased genetic variance after a population bottleneck . Trends in Ecology & Evolution . 1990 ; 5 ( 7 ): 228 ± 230 .
24. Xu TT , Abbott RJ , Milne RI , Mao K , Du FK , Wu G , et al. Phylogeography and allopatric divergence of cypress species (Cupressus L.) in the Qinghai-Tibetan Plateau and adjacent regions . BMC Evolutionary Biology . 2010 ; 10 ( 1 ): 194 .
25. Du BG , Zhu DY , Yang YJ , Shen J , Yang FL , Su ZY . Living situation and protection strategies of endangered Urophysa rockii . Jiangsu Journal of Agricutural Sciences . 2010 ; 1 : 324 ± 325 .
26. Xie DF , Zhang L , Hu HY , Guo XL , He XJ . Fragmented habitat drives significant genetic divergence in the Chinese endemic plant, Urophysa henryi (Ranuculaceae) . Biochemical Systematics and Ecology 2016 ; 69 : 76 ± 82 .
27. Wang JX , He XJ , Xu W , Meng WK , Su ZY . Preliminary study on Urophysa rockii. II. Biological characteristics, ecological characteristics and community analysis . Journal Sichuan Forestry Sciecne Technology. 2011a; 32 ( 4 ): 28 ± 39 .
28. Zhang YX , Hu HY , He XJ . Genetic diversity of Urophysa rockii Ulbrich, an endangered and rare species, detected by ISSR . Acta Botanica Boreali-Occidentalia Sinica. 2013b; 33 ( 6 ): 1098 ± 1105 .
29. Zhang YX , Hu HY , Yang LJ , Wang CB , He XJ . Seed dispersal and germination of an endangered and rare species Urophysa rockii (Ranunculaceae) . Acta Botanica Boreali-Occidentalia Sinica. 2013c; 35 ( 3 ): 303 ± 309 .
30. Yuan QJ , Zhang ZY , Peng H , Ge S . Chloroplast phylogeography of Dipentodon (Dipentodontaceae) in southwest China and northern Vietnam . Molecular Ecology . 2008 ; 17 ( 4 ): 1054 ± 1065 . https://doi.org/ 10.1111/j. 1365 - 294X . 2007 . 03628 . x PMID : 18208489
31. Wang J , Gao PX , Kang M , Lowe AJ , Huang HW . Refugia within refugia: the case study of a canopy tree (Eurycorymbus cavaleriei) in subtropical China . Journal of Biogeography . 2009a; 36 ( 11 ): 2156 ± 2164 .
32. Zhou TH , Li S , Qian ZQ , Su HL , Huang ZH , Guo ZG , Zhao GF . Strong phylogeographic pattern of cpDNA variation reveals multiple glacial refugia for Saruma henryi Oliv. (Aristolochiaceae), an endangered herb endemic to China . Molecular Phylogenetics and Evolution . 2010 ; 57 ( 1 ): 176 ± 188 . https:// doi.org/10.1016/j.ympev. 2010 . 07 .001 PMID: 20637294
33. Hickerson MJ , Carstens BC , Cavender-Bares J , Crandall KA , Graham CH , Johnson JB , et al. Phylogeography's past, present and future: 10 years after Avise . Molecular Phylogenetics and Evolution . 2010 ; 54 ( 1 ): 291 ± 301 . https://doi.org/10.1016/j.ympev. 2009 . 09 .016 PMID: 19755165
34. Petit RJ , Pineau E , Demesure B , Bacilieri R , Ducousso A , Kremer A . Chloroplast DNA footprints of postglacial recolonization by oaks . Proceedings of the National Academy of Sciences USA . 1997 ; 94 ( 18 ): 9996 ± 10001 .
35. White TJ , Bruns T , Lee S , Taylor J. Amplification and direct sequencing of fungal ribosomal RNA genes for phylogenetics. PCR protocols: A guide to methods and applications . 1990 ; 18 ( 1 ): 315 ± 322 .
36. Wright SD , Yong CG , Wichman SR , Dawson JW , Gardner RC . Stepping stones to Hawaii: a transequatorial dispersal pathway for Metrosideros (Myrtaceae) inferred from nrDNA (ITS + ETS) . Journal of Biogeography . 2001 ; 28 ( 6 ): 769 ± 774 .
37. Hamilton MB . Four primer pairs for the amplification of chloroplast intergenic regions with intraspecific variation . Molecular Ecology . 1999 ; 8 : 521 ± 523 . PMID: 10199016
38. Taberlet P , Gielly L , Pautou G , Bouvet J . Universal primers for amplification of three non-coding regions of chloroplast DNA . Plant Molecular Biology . 1991 ; 17 ( 5 ): 1105 ± 1109 . PMID: 1932684
39. Burland TG . DNASTAR's Lasergene sequence analysis software . Bioinformatics Methods and Protocols . 2000 ; 132 : 71 ± 91 .
40. Tamura K , Peterson D , Peterson N , Stecher G , Nei M , Kumar S. MEGA5: molecular evolutionary genetics analysis using maximum likelihood, evolutionary distance, and maximum parsimony methods . Molecular Biology and Evolution . 2011 ; 28 ( 10 ): 2731 ± 2739 . https://doi.org/10.1093/molbev/ msr121 PMID: 21546353
41. Li LF , Ping D , Liao QL , Xiao HX . Genomic and EST microsatellite markers for Aquilegia flabellata and cross-amplification in A. oxysepala (Ranunculaceae) . American Journal of Botany . 2011 ; 98 ( 8 ): 213 ± 215 .
42. Librado P , Rozas J. DnaSP v5: a software for comprehensive analysis of DNA polymorphism data . Bioinformatics . 2009 ; 25 ( 11 ): 1451 ± 1452 . https://doi.org/10.1093/bioinformatics/btp187 PMID: 19346325
43. Nei M. Molecular evolutionary genetics . Columbia: Columbia University Press. 1987 .
44. Grivet D , Petit RJ . Phylogeography of the common ivy (Hedera sp.) in Europe: genetic differentiation through space and time . Molecular Ecology . 2002 ; 11 ( 8 ): 1351 ± 1362 . PMID: 12144657
45. Pons O , Petit RJ . Measuring and testing genetic differentiation with ordered versus unordered alleles . Genetics . 1996 ; 144 ( 3 ): 1237 ± 1245 . PMID: 8913764
46. Excoffier L , Lischer HEL . Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows . Molecular Ecology Resources . 2010 ; 10 ( 3 ): 564 ± 567 . https:// doi.org/10.1111/j.1755- 0998 . 2010 . 02847 . x PMID : 21565059
47. Bandelt HJ , Forster P , Rohl A . Median-joining networks for inferring intraspecific phylogenies . Molecular Biology and Evolution . 1999 ; 16 ( 1 ): 37 ± 48 . PMID: 10331250
48. Hutchison DW , Templeton AR . Correlation of pairwise genetic and geographic distance measures: inferring the relative influences of gene flow and drift on the distribution of genetic variability . Evolution . 1999 ; 1898 ± 1914 . https://doi.org/10.1111/j.1558- 5646 . 1999 .tb04571. x PMID: 28565459
49. Rousset F. Genetic differentiation and estimation of gene flow from F-statistics under isolation by distance . Genetics . 1997 ; 145 ( 4 ): 1219 ± 1228 . PMID: 9093870
50. Peakall R , Smouse PE . GENALEX 6: genetic analysis in Excel. Population genetic software for teaching and research . Molecular Ecology Notes . 2006 ; 6 ( 1 ): 288 ± 295 .
51. Kimura M. A simple method for estimating evolutionary rates of base substitutions through comparative studies of nucleotide sequences . Journal of Molecular Evolution . 1980 ; 16 ( 2 ): 111 ± 120 . PMID: 7463489
52. Van Oosterhout C , Hutchinson WF , Wills DP , Shipley P . MICRO-CHECKER: software for identifying and correcting genotyping errors in microsatellite data . Molecular Ecology Resources . 2004 ; 4 ( 3 ): 535 ± 538 .
53. Pritchard JK , Stephens M , Donnelly P . Inference of population structure using multilocus genotype data . Genetics . 2000 ; 155 ( 2 ): 945 ± 959 . PMID: 10835412
54. Evanno G , Regnaut S , Goudet J . Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study . Molecular Ecology . 2005 ; 14 ( 8 ): 2611 ± 2620 . https://doi.org/10. 1111/j. 1365 - 294X . 2005 . 02553 . x PMID : 15969739
55. Swofford DL . PAUP*: Phylogenetic analysis using parsimony (*and other methods), Version 4 . Sinauer Associates, Sunderland; 2003 .
56. Ronquist F , Huelsenbeck JP . MrBayes 3: Bayesian phylogenetic inference under mixed models . Bioinformatics . 2003 ; 19 ( 12 ): 1572 ± 1574 . PMID: 12912839
57. Nylander JAA . MrModeltest 2.2. Computer program and documentation distributed by the author . Evolutionary Biology Centre , Uppsala: Uppsala University; 2004 .
58. Stamatakis A , Hoover P , Rougemont J. A rapid bootstrap algorithm for the RAxML web servers . Systematic Biology . 2008 ; 57 ( 5 ): 758 ± 771 . https://doi.org/10.1080/10635150802429642 PMID: 18853362
59. Rogers AR , Harpending H . Population growth makes waves in the distribution of pairwise genetic differences . Molecular Biology and Evolution . 1992 ; 9 ( 3 ): 552 ± 569 . PMID: 1316531
60. Tajima F. Statistical method for testing the neutral mutation hypothesis by DNA polymorphism . Genetics . 1989 ; 123 ( 3 ): 585 ± 595 . PMID: 2513255
61. Fu YX . Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection . Genetics . 1997 ; 147 ( 2 ): 915 ± 925 . PMID: 9335623
62. Knobloch E , Mai DH . Monographie der FruÈchte und Samen in der Kreide von Mitteleuropa . Rozpravy uÂstrednÂõho uÂstavu geologickeÂnho Praha . 1986 ; 47 :1± 219 .
63. Bastida JM , AlcaÂntara JM , Rey PJ , Vargas P , Herrera CM . Extended phylogeny of Aquilegia: the biogeographical and ecological patterns of two simultaneous but contrasting radiations . Plant Systematics and Evolution . 2010 ; 284 ( 3 ±4): 171 ± 185 .
64. Posada D. jModelTest: phylogenetic model averaging . Molecular Biology and Evolution . 2008 ; 25 ( 7 ): 1253 ± 1256 . https://doi.org/10.1093/molbev/msn083 PMID: 18397919
65. Drummond AJ , Nicholls GK , Rodrigo AG , Solomon W. Estimating mutation parameters, population history and genealogy simultaneously from temporally spaced sequence data . Genetics . 2002 ; 161 ( 3 ): 1307 ± 1320 . PMID: 12136032
66. Anderson CL , Bremer K , Friis EM . Dating phylogenetically basal eucots using rbcl sequences and multiple fossil reference points . American Journal of Botany . 2005 ; 92 : 1737 ± 1748 . https://doi.org/10. 3732/ajb.92.10.1737 PMID: 21646091
67. WikstroÈm N , Savolainen V , Chase M. Evolution of angiosperms: calibrating the family tree . Proceedings of the Royal Society of London B: Biological Sciences . 2001 ; 268 ( 1482 ): 2211 ± 2222 .
68. Ro K-E , Keener CS , McPheron BA . Molecular phylogenetic study of the Ranunculaceae: utility of the nuclear 26S ribosomal DNA in inferring intrafamilial relationships . Molecular Phylogenetics and Evolution . 1997 ; 8 : 117 ± 127 . https://doi.org/10.1006/mpev. 1997 .0413 PMID: 9299218
69. Drummond AJ , Rambaut A . BEAST: Bayesian evolutionary analysis by sampling trees . BMC Evolutionary Biology . 2007 ; 7 ( 1 ): 214 .
70. Li CY . Classification and Systematics of the Aquilegiinae Tamura . Beijing: The Chinese Academy of Science . 2006 ; 7 : 3 .
71. Wolfe KH , Li WH , Sharp PM . Rates of nucleotide substitution vary greatly among plant mitochondrial, chloroplast, and nuclear DNAs . Proceedings of the National Academy of Sciences USA . 1987 ; 84 ( 24 ): 9054 ± 9058 .
72. Rambaut A. 2009 . FIGTREE1.3 .1 [online]. Available from http://tree.bio.ed.ac.uk/software/figtree/; 2009 .
73. Zhang HX , Zhang ML , Sanderson SC . Retreating or standing: Responses of forest species and steppe species to climate change in arid Eastern Central Asia . PloS One. 2013a; 8 ( 4 ): e61954 .
74. Hamrick JL , Godt MJW . Effects of life history traits on genetic diversity in plant species . Philosophical Transactions Of The Royal Society Of London Series B-Biological Sciences . 1996 ; 351 ( 1345 ): 1291 ± 1298 .
75. Dauber J , Biesmeijer JC , Gabriel D , Kunin WE , Lamborn E , Meyer B , et al. Effects of patch size and density on flower visitation and seed set of wild plants: a pan-European approach . Journal of Ecology . 2010 ; 98 ( 1 ): 188 ± 196 .
76. Vranckx GUY , Jacquemyn H , Muys B , Honnay O . Meta-analysis of susceptibility of woody plants to loss of genetic diversity through habitat fragmentation . Conservation Biology . 2012 ; 26 ( 2 ): 228 ± 237 . https://doi.org/10.1111/j.1523- 1739 . 2011 . 01778 . x PMID : 22044646
77. Qi XS , Chen C , Comes HP , Sakaguchi S , Liu YH , Tanaka N , et al. Molecular data and ecological niche modelling reveal a highly dynamic evolutionary history of the East Asian Tertiary relict Cercidiphyllum (Cercidiphyllaceae) . New Phytologist . 2012 ; 196 ( 2 ): 617 ± 630 . https://doi.org/10.1111/j.1469- 8137 . 2012 . 04242 . x PMID : 22845876
78. Wang YH , Jiang WM , Comes HP , Hu FS , Qiu YX , Fu CX . Molecular phylogeography and ecological niche modelling of a widespread herbaceous climber, Tetrastigma hemsleyanum (Vitaceae): insights into Plio±Pleistocene range dynamics of evergreen forest in subtropical China . New Phytologist. 2015 ; 206 ( 2 ): 852 ± 867 . https://doi.org/10.1111/nph.13261 PMID: 25639152
79. Yang AH , Dick CW , Yao XH , Huang HW . Impacts of biogeographic history and marginal population genetics on species range limits: a case study of Liriodendron chinense . Scientific Reports . 2016 ; 6 .
80. Barrett SCH . Mating strategies in flowering plants: the outcrossing-selfing paradigm and beyond . Philosophical Transactions Of The Royal Society Of London Series B-Biological Sciences . 2003 ; 358 ( 1434 ): 991 ± 1004 .
81. Holsinger KE . Reproductive systems and evolution in vascular plants . Proceedings of the National Academy of Sciences USA . 2000 ; 97 ( 13 ): 7037 ± 7042 .
82. Petit RJ , Aguinagalde I , de Beaulieu JL , Bittkau C , Brewer S , Cheddadi R , et al. Glacial refugia: hotspots but not melting pots of genetic diversity . Science . 2003 ; 300 ( 5625 ): 1563 ± 1565 . https://doi.org/ 10.1126/science.1083264 PMID: 12791991
83. Harrison TM , Copeland P , Kidd WSF , Yin A. Raising Tibet . Science . 1992 ; 255 ( 5052 ): 1663 ± 1670 . https://doi.org/10.1126/science.255.5052.1663 PMID: 17749419
84. Qiu YX , Li JH , Liu HL , Chen YY , Fu CX . Population structure and genetic diversity of Dysosma versipellis (Berberidaceae), a rare endemic from China . Biochemical Systematics and Ecology . 2006 ; 34 ( 10 ): 745 ± 752 .
85. Wang CB , Wang T , Su YJ . Phylogeography of cephalotaxus oliveri (cephalotaxaceae) in relation to habitat heterogeneity, physical barriers and the uplift of the yungui plateau . Molecular Phylogenetics and Evolution . 2014 ; 80 : 205 ± 216 . https://doi.org/10.1016/j.ympev. 2014 . 08 .015 PMID: 25160902
86. An Z , Kutzbach JE , Prell WL , Porter SC . Evolution of Asian monsoons and phased uplift of the Himalaya-Tibetan plateau since Late Miocene times . Nature . 2001 ; 411 ( 6833 ): 62 ± 66 . https://doi.org/10. 1038/35075035 PMID: 11333976
87. Wang ZH , Yang LM , Yang CB . The Study about the Typical Continental Ecosystem of Yunnan-Guizhou Plateau . Beijing: Science Press. 2011b .
88. Dai PF . Phylogeography, phylogeny and genetic diversity of Chimonanthus (Calvcanthaceae) . Shanxi, Northwest University (Doctoral dissertation). 2012 .
89. Liu J , Wang Y , Wang A , Hideaki O , Abbott RJ . Radiation and diversification within the Ligularia±Cremanthodium±Parasenecio complex (Asteraceae) triggered by uplift of the Qinghai±Tibetan Plateau . Molecular Phylogenetics and Evolution . 2006 ; 38 : 31 ± 49 . https://doi.org/10.1016/j.ympev. 2005 . 09 .010 PMID: 16290033
90. Li ZJ , Yu GH , Rao DQ , Yang JX . Phylogeography and demographic history of Babina pleuraden (Anura, Ranidae) in southwestern China . PLoS ONE . 2012 ; 7:e34013 . https://doi.org/10.1371/journal. pone. 0034013 PMID: 22448286
91. Lu GJ , Lin AQ , Luo JH , Blondel DV , Meiklejohn KA , Sun K , et al. Phylogeography of the Rickett's bigfooted bat, Myotis pilosus (Chiroptera: Vespertilionidae): a novel pattern of genetic structure of bats in China . BMC Evolutionary Biology . 2013 ; 13 ( 1 ): 241 .
92. Zhang JJ , Li ZZ , Fritsch PW , Tian H , Yang AH , Yao XH . Phylogeography and genetic structure of a Tertiary relict tree species, Tapiscia sinensis (Tapisciaceae): implications for conservation . Annals of Botany . 2015 ; mcv112 .
93. Shimono Y , Watanabe M , Hirao AS , Wada N , Kudo G . Morphological and genetic variations of Potentilla matsumurae (Rosaceae) between fellfield and snowbed populations . American Journal of Botany . 2009 ; 96 ( 4 ): 728 ± 737 . https://doi.org/10.3732/ajb.0800242 PMID: 21628228
94. Hewitt GM . Genetic consequences of climatic oscillations in the Quaternary . Philosophical Transactions Of The Royal Society Of London Series B-Biological Sciences . 2004 ; 359 ( 1442 ): 183 ± 195 .
95. Qing GF . The study about geographic origin and evolution regularity of Pinus massoniana . Journal of Forestry Research . 2002 ; 15 : 406 ± 412 .
96. Chen DM , Zhang XX , Kang HZ , Sun X , Yin S , Du H , et al. Phylogeography of Quercus variabilis based on chloroplast DNA sequence in East Asia: multiple glacial refugia and mainland-migrated island populations . PloS One . 2012 ; 7 ( 10 ):e47268. https://doi.org/10.1371/journal.pone. 0047268 PMID: 23115642
97. Flanders J , Wei L , Rossiter SJ , Zhang S. Identifying the effects of the pleistocene on the greater horseshoe bat, Rhinolophus ferrumequinum, in east Asia using ecological niche modelling and phylogenetic analyses . Journal of Biogeography . 2011 ; 38 ( 3 ): 439 ± 452 .
98. Wang JX , He XJ . Preliminary study on Urophysa rockii. I. Literature research and biological characteristics of Urophysa rockii . Journal Sichuan Forestry Sciecne Technology . 2011 ; 32 : 69 ± 73 .
99. JaquieÂry J , Guillaume F , Perrin N. Predicting the deleterious effects of mutation load in fragmented populations . Conservation Biology . 2009 ; 23 ( 1 ): 207 ± 218 . https://doi.org/10.1111/j.1523- 1739 . 2008 . 01052 . x PMID : 18847439
100. Lande R . Anthropogenic, ecological and genetic factors in extinction and conservation . Researches on Population Ecology . 1998 ; 40 ( 3 ): 259 ± 269 .