Landscape Genetic Structure of a Streamside Tree Species Euptelea pleiospermum (Eupteleaceae): Contrasting Roles of River Valley and Mountain Ridge

We used landscape genetics and statistical models to test how landscape features influence connectivity or create barriers to dispersal for a mountain riparian tree species, Euptelea pleiospermum. Young leaves from 1078 individuals belonging to 36 populations at elevations of 900–2000 m along upper reaches of four rivers were genotyped using eight nuclear microsatellite markers. We found no evidence for the unidirectional dispersal hypothesis in E. pleiospermum within each river. The linear dispersal pattern along each river valley is mostly consistent with the “classical metapopulaton” model. Mountain ridges separating rivers were genetic barriers for this wind-pollinated tree species with anemochorous seeds, whereas river valleys provided important corridors for dispersal. Gene flow among populations along elevational gradients within each river prevails over gene flow among populations at similar elevations but from different rivers. This pattern of gene flow is likely to promote elevational range shifts of plant populations and to hinder local adaptation along elevational gradients. This study provides a paradigm to determine which of the two strategies (migration or adaptation) will be adopted by mountain riparian plants under climate warming.


Introduction
Migration and adaptation are two main responses of mountain plants to climate warming [1]. The predominant direction of gene flow (along elevational gradients or among similar elevations) is of importance in determining which of the two strategies will be adopted [2,3,4,5]. Here, we propose two models of gene flow ( Fig. 1) for mountain plants inhabiting restricted linear habitats (e.g. ridge, streamside, valley or adjacent mountains etc.) that transverse elevation gradients. First, gene flow among populations at similar elevations but along adjacent linear habitats is higher than that within the same linear habitat at different elevations ( Fig. 1A) because of similar microclimate and flowering times (phenology effect). Support for this model of gene flow has been shown for Poa hiemata [2] and Castanopsis eyrei [3]. Second, gene flow is more likely to occur within each linear habitat than between them (Fig. 1B) because of barriers (e.g. valley or ridge) between linear habitats along elevational gradients (geographic barrier effect). Support for this model has been found in Ainsliaea faurieana [6] and Silene ciliata [5]. The first model of gene flow (Fig. 1A) could promote local adaptation within an elevational zone but also potentially hinder elevational range shifts of plant populations, whereas the second model (Fig. 1B) would render the opposite trend.
For mountain riparian plants, gene flow is influenced by mountain landscape characteristics, such as elevation, stream flow, river valley and mountain ridge [6,7,8]. As dispersal is an important determinant of the level and direction of gene flow, it is essential for understanding how landscape features influence connectivity or create barriers to dispersal [9,10,11].
River valleys are long-recognized corridors for riparian or aquatic plants, which are characterized by linear distribution patterns [12,13]. Previous studies have proposed six hypothetical models about the patterns of dispersal and connectivity between linearly arranged populations [12,14]. The two extreme models are ''spatially extended population'' and ''regional ensemble'' [15] with free and without any ongoing gene flow among populations, respectively. For species characterized predominantly by unidirectional water-mediated dispersal, the linear asymmetrical adjacent flow model and the linear asymmetrical non-adjacent flow model were proposed [16]. The two models predict that downstream populations may harbor higher genetic diversity [17], and a few studies on riparian plants have confirmed this prediction [6,18,19]. However, wind-or animal-mediated upstream seed and/or pollen dispersal, which are likely to be promoted by valleys (corridors) along rivers, may weaken or counteract this cumulative effect [8,12,20,21,22]. Therefore, the two bidirectional dispersal models, the ''classical stepping-stone'' model and the ''classical metapopulation'' model, seem to be more common in studied linear populations, especially the latter.
Mountain ridge is the predominant factor causing limited gene flow and genetic differentiation in mountain plants [6,23]. As streamside in mountainous areas is always at the bottom of a deep valley, linear distribution patterns constrained by highlands on both sides render populations with sharp boundaries and strong natural fragmentation; such spatial structure hampers gene flow among rivers. Therefore, populations at similar elevation along different rivers may be genetically separated by ridges. The linear distribution pattern can, in turn, facilitate gene flow within river valleys, especially for wind-dispersed plants. Therefore, we expected that gene flow of mountain riparian plants would be consistent with the second model (Fig. 1B).
Euptelea pleiospermum Hook. f. et Thoms (Eupteleaceae) is a diploid (2n = 28) [24], deciduous, broad-leaved Tertiary-relict tree species endemic to China, India and Burma [25]. In China, this species distributes across wide altitudinal ranges along streamsides (720-3600 m a.s.l.) [25]. As other riparian plant species, E. pleiospermum is restricted to linear suitable habitats. It is windpollinated and flowers in early spring before leaf-out. Adult trees produce numerous indehiscent fruits. The small and light samaras are dispersed first by gravity and secondarily by wind and/or water [26]. For this species, Li [27] has observed water-mediated seed dispersal in a downstream direction by simulation experiments with dyed real and artificial fruits.
The purpose of the present study is to examine the landscape genetic structure of Euptelea pleiospermum populations located along four main rivers (elevation gradients) within one mountain. Specially, the following questions were addressed: (1) does the recognized water-mediated seed dispersal result in genetic diversity accumulation in downstream populations? If not, which of the other linear models is suitable for populations within rivers? (2) Do river valleys and mountain ridges act as migration corridors and genetic barriers, respectively?

Ethics Statement
All necessary permits were obtained for the described field studies. Mr. Jingyuan Yang from the Administration Bureau of the Shennongjia National Nature Reserve issued the permission for each location.

Study Area and Sample Collection
The study area is located in mountain riparian forests of the Shennongjia National Nature Reserve (SNNR) (31u219200-31u369200 N, 110u039050-110u339500 E), central China (Fig. 2). Its summit Shennongding (3105.4 m) is the highest peak in central China [28]. SNNR is an important part of the south-central China biodiversity hot-spot and is rich of Tertiary relicts and endemic plants [29,30]. There are four main river systems in the Shennongjia Mountains: Yandu River and Xiangxi River on the south-facing slope, which are tributaries of the Yangtze River, and Nan River and Du River on the north-facing slope, which are the tributaries of the Han River (Fig. 2). Euptelea pleiospermum is a common tree species in riparian forests at elevations between ,1200 and ,1900 m in this area [31].
Fieldwork was carried out in April and May 2008. Young leaves were collected and immediately dried in a 10 cm 65 cm plastic bag containing silica gel. This was repeated every ,100 m from high to low altitude along the four rivers. The sampling range depended on the altitudinal range of E. pleiospermum along each river. In this study, a total of 36 populations were sampled with six to 12 populations for each river (Fig. 2). Because of the sporadic distribution at some altitudes, sample sizes were sometimes relatively low, ranging from six to 61 per population (Table 1 and Fig. 2). A total of 1078 individuals were analyzed in this study.

DNA Extraction and Microsatellite Genotyping
All samples were stored at 4uC pending DNA extraction. Total DNA was extracted from approximately 0.5 g dried leaves using a modified CTAB method [32]. The DNA sample was diluted to 5-10 ng mL 21 . Fourteen published primer pairs [33] were tested using 36 individuals (one individual randomly chosen from each population) and eight nuclear microsatellite loci (EP021, EP036, EP059, EP081, EP087, EP091, EP278, and EP294) that detected a suitable level of polymorphism were selected to determine genotypes of the 1078 samples.
The microsatellite genotyping procedure was performed as described by Wei and Jiang [34]. PCR amplifications were carried out in a total volume of 10 mL consisting of 10-50 ng of template DNA, 1 mL 106 reaction buffers, 1.5 mM of MgCl 2 , 0.2 mM of dNTPs, 0.25 mM of each primer, and 0.5 U of Taq polymerase (Tiangen Biotech Co., LTD., Beijing, China ). PCR reaction was conducted in a PTC-100 thermal cycler (MJ Research, Inc., Cambridge, MA, USA). PCR reactions were performed with the following cycle profile: an initial denaturation at 94uC for 5 min; then 35 cycles of denaturation at 94uC for 50 s, annealing at 56uC for 50 s and extension at 72uC for 1 min; and final extension at 72uC for 10 min. After amplification, identification of the alleles was based on the position of the fragments in relation to a 25-bp marker ladder on a polyacrylamide gel of high resolution with silver stain.
Data Analysis GENETIX 4.05 [35] was used to calculate the following population diversity indices: mean number of alleles per locus (A), expected (H E ) and observed (H O ) heterozygosity, and inbreeding coefficient (F IS ) across all loci and locus-population combinations. FSTAT 2.9.3.2 [36] was used to estimate allelic richness (A R ) for each population and to test for linkage disequilibrium (LD). Departure from Hardy-Weinberg equilibrium (HWE) was tested by Fisher's exact test in GENEPOP version 4.0 [37]. The significant values for the LD were corrected for multiple comparisons by Bonferroni correction [38]. We used the program MICRO-CHECKER [39] to test for null alleles.
To reveal the pattern of genetic diversity along rivers, regression analyses between genetic diversity parameters (H E and A R ) and the position of populations along the course of the river were conducted using the SPSS statistics program (SPSS 15.0 for Windows; SPSS Inc., Chicago, IL, USA).
To provide direct evidence for the dispersal pattern along each river, contemporary and historical migration rates between populations along each river were calculated by using the programs BayesAss version 1.3 [40] and MIGRATE version 3.2.2.windows [41], respectively. BayesAss was run using a Markov chain Monte Carlo (MCMC) length of 3,000,000 with a burn-in period of 1,000,000 (initial conditions of Dp, Dm, and DF = 0.15). MIGRATE uses Bayesian inference with long chains (500,000 steps sampled, 5,000 steps recorded) and 1,000 burn-in per chain. We then compared mean values of downstream and upstream migration rates along each river valley, and significance was tested with a paired-samples t-test.
Pairwise population differentiation was estimated at two scales: among populations and among rivers. F ST [42] and D EST [43] were calculated using FSTAT 2.9.3.2 [36] and SMOGD 1.2.5 [44], respectively. Significance tests for F ST were based on 1000 randomizations. We compare genetic differentiation among populations from different rivers and that among populations within the same river (inter-F ST vs intra-F ST and inter-D EST vs intra-D EST , respectively), and significance was tested by performing an independent-samples t-test.
To analyze the level and partitioning of genetic variation among geographical groups, among populations within geographical groups, and within populations, analysis of molecular variance (AMOVA) [45] was performed with 1000 permutations in ARLEQUIN version 3.1 [46]. There were two kinds of geographical groups: populations along each river (Yandu, Xiangxi, Nan, and Du) and populations in each altitudinal zone (12, from 900 m to 2000 m).
To test whether the individuals from the same river valley or from similar altitudes can be clustered into genetically homogeneous groups, nonspatial and spatially sensitive Bayesian clustering were implemented by using STRUCTURE version 2.3.3 [47,48,49,50] and TESS version 2.3.1 [51]. Previous STRUC-TURE analysis employed a Markov chain Monte Carlo (MCMC) approach to cluster individuals into K panmictic groups without a priori assignment of individuals to sampling locations and by minimizing deviations from Hardy-Weinberg equilibrium and linkage equilibrium. The LOCPRIOR model, recently developed by Hubisz et al. [52], allowed us to infer weak population structure with the assistance of sample group information. To guide an empirical estimate of the number of genetic populations, we performed five replicates of each simulation from K = 1 to K = N P +3 (N P , the number of sampling locations), with a burn-in period of 100,000, and MCMC iterations of 1,000,000, based on the LOCPRIOR model and the admixture model with correlated allele frequencies. We applied Evanno et al.'s [53] ad hoc statistic, DK, to calculate the uppermost hierarchical level of structure (K). For the spatially sensitive TESS analysis (admixture model: BYM; spatial interaction parameter: 0.6), we set 10 runs per K max (maximal number of clusters), for K max ranging from two to six, with 300,000 total sweeps and a burn-in of 100,000.We used the program DISTRUCT version 1.1 [54] to create bar charts for the output data derived from STRUCTURE and TESS analyses.
To test for isolation by distance (IBD) [55], the correlation between genetic distances (F ST /(1-F ST )) and altitudinal distances (logarithmically transformed) among populations was executed using the Mantel test implemented in GenAlEx 6.41 [56]. We also performed partial Mantel test to check for the effect of altitudinal distance after accounting for geographical distance by using FSTAT 2.9.3.2 [36]. The IBD tests were performed at two spatial scales: each of the four rivers and the whole Shennongjia Mountains. Latitude and longitude coordinates for each population were used to calculate pairwise geographical distances between populations.

Microsatellite Polymorphism, Hardy-Weinberg Equilibrium and Linkage Disequilibrium
A total of 102 alleles at eight nuclear microsatellite loci were revealed across 1078 individuals of 36 sampled populations of E. pleiospermum. At the species level, all eight nuclear microsatellite markers were highly polymorphic, having 6-29 alleles per locus with an average of 12.8. At the intra-population level, mean number of alleles per locus (A) ranged from 3.88 to 6.88, and allelic richness (A R ) varied from 3.64 to 4.58 (Table 1 and Fig. 2A). The expected heterozygosity (H E ) ranged from 0.555 to 0.687 (Fig. 2B),  Table 1). The level of inbreeding for each population ranged from -0.104 to 0.335 and significant inbreeding coefficients (F IS ) were detected in 34 of 36 populations ( Table 1).
The test for HWE found that 135 of 288 locus-population combinations were significant, including 32 at locus EP021, 21 at EP036, three at EP059, 17 at EP081, 22 at EP087, nine at EP091, eight at EP278 and 23 at EP294. However, no locus displayed consistent deviation from HWE across all populations. Furthermore, MICROCHECKER revealed potential null alleles in five of the 36 populations, with three at locus EP021 and two at EP297. This suggested that departure from HWE is potentially resulted from both technical (e.g. null allele) and biological (e.g. inbreeding) reasons. The test for genotypic disequilibrium found that only one (EP0366EP294) of 28 locus pairs showed significant genotypic disequilibrium in 12 populations. Because no consistent genotypic disequilibrium was found between any locus pair across all populations, all loci were used for further analyses.

Linear Population Model along Rivers
Linear regression analyses did not reveal any significant negative relationship between genetic diversity and the position of populations along the course of the rivers ( Table 2 and Fig. 2), suggesting that there was no accumulation of genetic diversity in downstream populations. In accordance with this result, we found no significant difference between downstream and upstream migration rates, although the former was a little higher than the latter along all the four rivers ( Table 3, Table S1 and Table S2). Moderate levels of pairwise differentiation (F ST and D EST ) were found between populations within each river valley. Significance test revealed that most values of F ST were significantly different from zero (Table S3).
Mantel test detected significant IBD within the Yandu River and Du River (Table 4). However, after controlling for geographical distance, the partial Mantel test detected no significant relationships between genetic and altitudinal distances within any of the four rivers ( Table 4).

Roles of Mountain Ridges and River Valley
When individuals from the same river were seen as one population, all pairwise F ST values between rivers were significantly different from zero, and the pairwise D EST values also showed similar levels of genetic differentiation (Table 5).
We found moderate levels of genetic differentiation among most pair of populations from different rivers (Table S3). Furthermore, average F ST values were higher for inter-river population pairs than for intra-river population pairs (F ST inter: 0.06960.001; F ST intra: 0.06160.003, P = 0.007), and D EST values showed the same trend (D EST inter: 0.09160.002; D EST intra: 0.08260.004, P = 0.032) (See Table S3 for details).
The results of the AMOVA were shown at Table 6. There was significant differentiation among different rivers (0.99%; P,0.001). However, no significant differentiation was found among populations located at different altitudinal zones (0.20%; P = 0.206).
With the four rivers as sample group information, STRUC-TURE analysis revealed a clear peak of Evanno et al.'s [53] ad hoc statistic DK at K = 4. The bar plot assuming K = 4 showed that individuals from different rivers were clustered into distinct groups (Fig. 3A). However, the same analysis did not clearly detect genetic structure when the 12 altitudinal zones were used in the LOCPRIOR model (data not shown). Results from TESS paralleled those obtained from STRUCTURE, although there are some minor differences (Fig. 3).
Both Mantel and partial Mantel tests found significant correlations between F ST -based genetic distances and altitudinal distances within the whole Shennongjia Mountains (Table 4).

Linear Population Model within Rivers
We found no evidence for the unidirectional dispersal hypothesis in E. pleiospermum (Tables 2 and 3, Fig. 2), although watermediated seed dispersal in a downstream direction has been observed [27]. In other words, upstream dispersal of pollen and seed mediated by wind are efficient mechanisms to eliminate the downstream accumulation of genetic diversity caused by watermediated seed dispersal. Besides upstream dispersal, Honnay et al. [22] recently pointed out that higher seed recruitment opportunities in upstream habitats due to density dependence of recruitment is another explanation. However, this seems not the case for E. pleiospermum in the study area, because the pattern of genetic diversity along the Yandu River supports the centralmarginal hypothesis [57].   The ''regional ensemble'' and ''classical stepping-stone'' models were rejected because IBD was not observed within any of the four rivers ( Table 4), indicating that gene flow along each river was across a wide elevational range. The ''spatially extended population'', a single genetically uniform panmictic unit, was also rejected because most of the genetic differentiation values among populations were significant. As inferred above, the bidirectional long-distance dispersal within each river valley and the significant genetic differentiation among populations indicated that the linear dispersal pattern for E. pleiospermum is most consistent with the ''classical metapopulaton'' model. This model was also confirmed in another riparian tree species, Populus nigra [20]; anemochory is the most important dispersal mode for this species as well.

Contrasting Roles of Mountain Ridges and River Valleys
Our results suggested that mountain ridges and river valleys played contrasting roles in shaping the genetic structure of E. pleiospermum at the landscape scale. Mountain ridges between rivers can be recognized as genetic barriers, while river valleys act as corridors for gene flow. In other words, gene flow of E. pleiospermum is consistent with the second model (Fig. 1B).
Four results lead to the conclusion that mountain ridges act as genetic barriers. First, individuals from each river were clustered into distinct groups (Fig. 3). Furthermore, the 12 altitudinal zones were not clustered into different groups by STRUCTURE analysis and this result suggested that populations from different rivers at similar elevations could not be clustered into the same group. Second, the F ST and D EST values for inter-river comparisons were higher than intra-river comparisons. This has been found in another mountain riparian plant, Ainsliaea faurieana [6]. Third, we found significant IBD in the whole Shennongjia Mountains but lack of IBD within all four rivers (Table 4). Fourth, when each river was seen as a population, we detected significant genetic differentiation among them (Table 5). This was confirmed by the results of AMOVA, which showed that there was significant differentiation among rivers, whereas no significant differentiation was found among altitudinal zones (Table 6). Mountain ridges have been found as genetic barriers for another wind dispersed tree species, Betula maximowicziana [23].
However, unambiguous evidence for gene flow among rivers is that the X6 (1520 m) and X8 (1389 m) populations from Xiangxi River were clustered into the Nan River group, as revealed by STRUCTURE and TESS analyses (Fig. 3). This can be partly explained by the ex situ conservation performed with plant materials collected from the Nan River. There is one national highway connecting the two rivers, and there was an ex situ conservation site at ,1400 m along the Xiangxi River in the past. Furthermore, it has been reported that roads can serve as dispersal corridors for plants [58,59]. However, it is not yet clear why X1 (2003 m) and X9 (1270 m) populations were clustered into the Yandu River group (Fig. 3B).
Our results also support that linear riparian zones in deep valleys provide important corridors for gene flow. First, as mentioned above, lack of IBD in each river indicated that gene flow along the river valley was not restricted (Table 4). This has previously been reported for riparian and aquatic plants [6,14]. Second, as mentation above, the F ST and D EST values for interriver comparisons were higher as compared to intra-river comparisons. Third, the Bayesian clustering approach divided individuals from the four rivers into four distinct groups (Fig. 3), indirectly indicating that gene flow within each river was relatively high.

Conclusions and Implications
We found no accumulation of genetic diversity at the downstream populations of E. pleiospermum within each of the four rivers. This was supported by the estimates of contemporary and historical migration rates, which revealed that the downstream migration rate was not significantly higher than that in upstream direction. The linear dispersal pattern along each river valley was mostly consistent with the ''classical metapopulaton'' model. Mountain ridges among rivers were genetic barriers for this wind-pollinated tree species with anemochorous seeds, whereas river valleys provided important corridors for dispersal. In summary, the streamside population of E. pleiospermum is consistent with the second model of gene flow (Fig. 1B) for mountain plants restricted to linear habitats along elevational gradients. Here, it should be mentioned that the variation of sample size, especially small number of sample in some of the populations (less than 20), may influence in our conclusions. Gene flow among E. pleiospermum populations along elevational gradients within each river prevailed over gene flow among populations at similar elevations but from different rivers. This pattern of gene flow could promote elevational range shifts of plant populations but potentially hinder local adaptation along elevational gradients. Therefore, range shift along elevations is the most likely strategy for E. pleiospermum in the face of climate warming. As mountain riparian zones are crucial refuges for tertiary-relict plants [60,61,62,63], this study provides a paradigm to determine which of the two strategies (migration or adaptation) will be adopted by these plants under climate warming.

Supporting Information
Table S1 Contemporary migration rates between Euptelea pleiospermum populations along each river.   1 (B). Each individual is represented by a thin vertical line, which is partitioned into four colored segments that represent the individual's estimated membership fractions in the four clusters. Black vertical lines separate individuals of different populations. Labels below the plot provide population codes, which are the same as in Table 1. Labels above the plot (Du River, Nan River, Xiangxi River and Yandu River) are sampling information. doi:10.1371/journal.pone.0066928.g003