Effect of Anthropogenic Landscape Features on Population Genetic Differentiation of Przewalski's Gazelle: Main Role of Human Settlement

Anthropogenic landscapes influence evolutionary processes such as population genetic differentiation, however, not every type of landscape features exert the same effect on a species, hence it is necessary to estimate their relative effect for species management and conservation. Przewalski's gazelle (Procapra przewalskii), which inhabits a human-altered area on Qinghai-Tibet Plateau, is one of the most endangered antelope species in the world. Here, we report a landscape genetic study on Przewalski's gazelle. We used skin and fecal samples of 169 wild gazelles collected from nine populations and thirteen microsatellite markers to assess the genetic effect of anthropogenic landscape features on this species. For comparison, the genetic effect of geographical distance and topography were also evaluated. We found significant genetic differentiation, six genetic groups and restricted dispersal pattern in Przewalski's gazelle. Topography, human settlement and road appear to be responsible for observed genetic differentiation as they were significantly correlated with both genetic distance measures [FST/(1−FST) and F′ST/(1−F′ST)] in Mantel tests. IBD (isolation by distance) was also inferred as a significant factor in Mantel tests when genetic distance was measured as FST/(1−FST). However, using partial Mantel tests, AICc calculations, causal modeling and AMOVA analysis, we found that human settlement was the main factor shaping current genetic differentiation among those tested. Altogether, our results reveal the relative influence of geographical distance, topography and three anthropogenic landscape-type on population genetic differentiation of Przewalski's gazelle and provide useful information for conservation measures on this endangered species.


Introduction
The potential effect of landscape features on dispersal and population differentiation of wild animals has been well recognized [1]. With the rapid expansion of human population and associated land use, concerns have been raised about the influence of anthropogenic landscape features because they could impede gene flow, lead to population isolation and genetic differentiation [2,3], and reduce genetic variation and evolutionary potential [4,5]. This research field is now a focus in population and conservation genetics [6], and many recent studies have made much progress including in invertebrates [e.g. 7], amphibians [e.g. 8], reptiles [e.g. 9] and mammals [e.g. 10]. However, despite some studies on the genetic effect of certain anthropogenic landscape type such as road [5,11], most previous researches in this field have tended to present the total effect of different anthropogenic landscapes. As not every type of landscape features exert the same effect on a species, it is necessary to assess their relative effect for species management and conservation. Moreover, given that landscape effects are a function of both environmental features [12,13] and the biological characters of a species [14,15], studies across more diverse taxa and regions are needed to fully understand the influences of man-made landscape features on wild animals. In central Asia, many ungulates are threatened or endangered [16], yet little is known about how anthropogenic landscape features affect these animals. Here, we carried out a landscape genetic study on Przewalski's gazelle (Procapra przewalskii), which inhabits a human-altered area on the Qinghai-Tibet Plateau in central Asia, to investigate the relative genetic effect of different types of landscape features.
The Przewalski's gazelle is one of the most endangered antelope species in the world [17]. It is endemic to China and a flagship species in the eastern part of the Qinghai-Tibet plateau [18]. The species once inhabited Gansu, Inner Mongolia, Ningxia and Qinghai Provinces, China but has experienced a sharp population reduction [19,20]. Now, several hundred of individuals survive in several localities in Qinghai Province near Qinghai Lake ( Figure 1) [21,22]. Accordingly, it is listed as 'Endangered' under the IUCN Red List (2008) [23] and as a Category I National Key Protected Wild Animal Species under the Wild Animal Protection Law of China [20]. Przewalski's gazelle inhabits on a vast grassland without obvious natural barriers between most populations, and the gazelle has high locomotive ability (45-55 kph) [24] and relatively short inter-population distances (8.5-205.8 km, Table 1), however, the total home range sizes of the gazelle are about 800 hm 2 under human disturbance [25], and previous study based on mitochondrial DNA (mtDNA) control region sequences revealed unexpected significant genetic differentiation among four populations [26]. The authors inferred that anthropogenic barriers probably influenced genetic differentiation, but there was no explicit analysis to support the inference [26]. Since new molecular markers and analytic methods are now available, we ask the following questions: what level of genetic differentiation is there in Przewalski's gazelle on the basis of microsatellite markers? Which factors caused genetic differentiation, anthropogenic or other factors? How important is the relative effect of each factor? Here, we considered five factors which possibly influenced genetic differentiation, including human settlement, road, railway, geographical distance and topography (lake and mountain). Given that the distribution and habitat of the gazelle are fragmented by many man-made landscape features (Figure 1), we hypothesize that anthropogenic landscape features influenced population genetic differentiation of Przewalski's gazelle. We expected that significant genetic differentiation would be observed, and certain types of anthropogenic landscape features would significantly correlate with the genetic distances.
In this study, we used 13 nuclear microsatellite loci and landscape genetic methods to test our hypothesis. The aims of this study are: (1) to estimate population genetic differentiation and dispersal pattern of the gazelle; (2) to identify the factors responsible for observed genetic differentiation and to reveal their relative effect. The results provide useful information for conservation of Przewalski's gazelle and, more broadly, yield empirical data for understanding the genetic effect of anthropogenic landscape features on wild ungulates across Eurasia grasslands.

Ethics Statement
Our research on Przewalski's gazelle in the Gonghe Basin and Qinghai Lake Basin was approved by the Chinese Wildlife Management Authority and conducted under the Wildlife Protection Law of the People's Republic of China.

Study area and sample collection
This study was conducted around the Qinghai Lake (98u269-100u549E, 36u139-37u369N). Qinghai Lake is the largest inland saline lake in China, and a perennial lake that freezes in winter. The lake has an area of 4300 km 2 and an average water depth of 16 m (max. 28 m). The study area is 2900-3800 m above sea level, surrounded by several mountains. There are three main types of man-made landscapes fragmenting the distribution and habitat of the gazelle: human settlements (county capital, town and village), roads (national highway, provincial highway and minor road) and railways ( Figure 1). Their characteristics and age were obtained from local chronicles of the four counties in the study area [27][28][29][30], literatures [31][32][33][34], the Fifth National Census Data (2000) and field survey ( Table 2). Samples of Przewalski's gazelle were collected from all known populations (P1-P9, Figure 1 [21,22]. Population P9 is located in the Gonghe Basin and probably isolated from other populations by Qinghainanshan Mountain. The remaining eight populations are located in the Qinghai Lake Basin surrounded by lofty mountains, with Datong Mountain on its north, Riyue Mountain on its east and Qinghainanshan Mountain on its south ( Figure 1). Skin samples were collected from wolf-killed gazelles discovered in the field from 2004-2007 and stored dry at 220uC. Fresh fecal samples were collected in winter (Nov-Dec 2006) to reduce genotyping error rates [35] and were preserved in 100% ethanol. Given our research experience with this species, we have detailed knowledge of the locations known to support the gazelle populations. After a herd of gazelles was located and observed for some time, we approached the gazelles and collected all fresh feces from the ground. Because P4 and P7 have very small population sizes, collecting fresh fecal sample is especially difficult for the two populations. Locations for each sample were recorded using a handheld global positioning system (Garmin Etrex Vista C, Garmin Ltd.).

DNA extraction and microsatellite amplification
Genomic DNA from skin samples was isolated using standard proteinase K digestion and phenol/chloroform extraction procedures [36], followed by a UNIQ-10 column (Sangon) purification. Fecal samples were extracted using QIAamp DNA Stool Mini Kit (QIAGEN) according to manufacturer protocol. Extraction blanks were used as negative controls. We used 13 species-transferred microsatellite markers (AF5, AGLA226, BM1225, CSSM43, IDVGA39, JAB8, RBP3, RT1, T156, TEXAN15, TGLA10, TGLA122 and TGLA378) in this study [37,38] on the basis of consistently amplifying clear and polymorphic products in fecal samples. One homozygote of each locus was sequenced to confirm it was a short tandem repeat (STR) in Przewalski's gazelle. Then, forward primers of the 13 markers were fluorescently labeled with FAM, HEX and TAMRA. We used a reaction volume of 10 ml, containing approximately 10 ng of genomic DNA, 0.2 mM of each primer, 0.1 mg/ml of bovine serum albumin (BSA, Biolabs) and 0.25 U HotStartTaq (QIAGEN). All PCR amplifications were carried out on a Thermo Hybaid MBS 0.2S cycler with an initial denaturation of 15 min at 95uC, followed by 40 cycles at 94uC for Human settlement characteristics were obtained from local chronicles [27][28][29][30], literatures [31][32][33][34], the Fifth National Census Data (2000) and field survey. Road and railway characteristics were acquired from local chronicles [27][28][29][30], local railway station record and field survey. The age of anthropogenic landscape features were obtained from local chronicles [27][28][29][30]

Reliability of genotyping results
We conducted two replicate PCRs for skin samples. Fecal samples were amplified using a modified multiple-tube procedure [39]. In practice, amplification was repeated minimally three times, all heterozygotes were observed in a minimum of two separate reactions and all homozygous showed identical profile of the replicates. Otherwise, we treated the alleles as missing data. Only samples with more than ten loci of allele data were included in statistical analyses. Using GIMLET version 1.3.3 [40], we computed the probability of unrelated individuals or full-sibs bearing an identical multi-locus genotype [P ID or P ID(Sibs) ]. Fecal samples with identity in all alleles or with only a single mismatch were considered to be deposited by the same individual, and duplicates were removed from the data set. The program MICRO-CHECKER version 2.2.3 [41] was used to check for microsatellite null alleles and scoring errors due to large allele drop-out or stuttering. Rates of genotyping error were calculated following the equations of Broquet and Petit [42].

Genetic data analyses
Deviations from Hardy-Weinberg Equilibrium (HWE) for each population and genotypic linkage disequilibrium (LD) across all pairs of loci were estimated in GENEPOP version 4.0 [43]. Bonferroni corrections were applied to tests involving multiple comparisons [44]. We used traditional F-statistics [45] to assess genetic differentiation, with pairwise F ST between populations and overall F ST across the study area. The significance of F ST values were assessed via 10 000 permutation procedure using FSTAT version 2.9.3.2 [46]. We also calculated the standardized measure of genetic differentiation F9 ST (estimated using G0 ST ) to control for the influence of within population heterozygosity [47][48][49][50], using GenoDive version 2.0b19 [51]. Standard errors of F9 ST estimates were obtained through jackknifing, and the significance of F9 ST values were assessed by permutation test (10 000 permutations) [51]. The genetic diversity within each population, measured as H S (equal to the expected heterozygosity), were obtained from GenoDive [51].
The program Barrier version 2.2 [52], which implements a computational geometry method and a Monmonier's Maximumdifference algorithm [53], was used to identify the genetic discontinuities within Przewalski's gazelle. This program provides the locations and robustness of the barriers (namely the genetic discontinuities), and visualizes these on a geographical map [52]. In practice, the nine populations of Przewalski's gazelle were connected by a Delaunay triangulation [54] based on each population's average geographic coordinates of sample locations, followed by the application of Monmonier's algorithm with F ST or F9 ST matrix. To test the robustness of estimated barriers, we used MICROSAT [55] to obtain 100 bootstrap matrices for Barrier significance analysis.
The Bayesian clustering method implemented in software STRUCTURE [56,57] was used to detect population genetic structure based on individual multilocus genotypes. Ten independent runs of K (number of genetic groups) = 1-9 were performed, using correlated allele frequencies [58] and admixture model, with 1 000 000 Markov Chain Monte Carlo (MCMC) repetitions after a 100 000 burn-in period. K was identified using the maximum values of Ln P(D) (the posterior probability of the data for a given K) and DK (the rate of change in the log probability of data between successive values of K) [59]. Bar plots which consist of individual assignment probabilities were applied to cluster the sampling populations into genetic groups.
To assess the genetic variance within and among several possible population groupings, a hierarchical AMOVA (Analysis of Molecular Variance) was performed using Arlequin version 3.11 [60] with 10 000 permutations. We conducted AMOVA analysis on four grouping patterns defined by topography, human settlement, road and railway respectively. These analyses allowed comparing the extent of genetic variance explained by groupings defined by different landscape features. Two additional AMOVA were run on groups identified by inferred genetic discontinuities and genetic structure, respectively.
For migration rate estimation, we used two types of likelihood computation (L home and L home /L max ) implemented in GeneClass version 2.0.h [61] to detect first generation migrants (F 0 ). L home is the likelihood of drawing an individual from the population in which it was sampled, it is appropriate to use when some source populations for immigrants are not sampled [61,62]. And L home / L max is the ratio of L home to the greatest likelihood observed in all sampled populations, it is suitable when all source populations are sampled [61,62]. The likelihood analyses were conducted using a Bayesian method [63], and the probability that an individual is a resident was calculated with a resampling algorithm [62] on 10 000 simulated individuals. To investigate the dispersal pattern of the gazelle, we performed a spatial autocorrelation analysis in GenAlEx 6 [64] which assesses the genetic similarity between pairs of individuals at different geographical distance classes. Because sample sizes were unevenly distributed across distances, variable distance classes were used in the analysis. We ran 1000 random permutations to test the 95% confidence intervals of the null hypothesis, and 1000 bootstraps to estimate 95% confidence intervals for autocorrelation index r of each distance class [65].
The significance of r was tested by comparing the estimated r with the 95% confidence interval about the null hypothesis of a random distribution.

Landscape genetic analyses
To assess which factors influenced genetic differentiation in Przewalski's gazelle, Mantel tests [66] were applied to test the Table 3. Matrix of categorical distances of human settlements (upper diagonal) and roads (lower diagonal) between nine populations.  (Table 1). For landscape features, a categorical matrix was generated describing whether there were ( = 1) or were not ( = 0) certain type of landscape features between populations. Then, the categorical matrix was used in the Mantel tests [68]. In practice, the categorical matrix was applied to the three main types of anthropogenic landscapes (human settlement, road and railway) and topography (including lake and mountain) ( Table 1, 3 and S1). Because geographical distance and landscape features were not independent, partial Mantel tests [69] were used to assess the relative effect of each significant factor inferred in the Mantel tests. All tests were executed using IBD 1.52 [70] with 10 000 permutations to determine statistical significance.
To choose the best model interpreting genetic differentiation based on Mantel results, Akaike's information criterion (AIC)  which represents how far a tentative model is from the true model was used [71,72]. We adjusted AIC values to AIC c values for relatively low number of pairwise comparisons in our data [73]. Lower AIC c scores imply closer approximation to the true model, and generally, models with AIC c values greater than 10 compared to the model with the lowest AIC c value are not supported [72]. Further, we applied Mantel and partial Mantel tests within a causal modeling framework [74] to assess the support for organizational models containing the significant factors inferred in the Mantel tests [using either This involved computing Mantel and partial Mantel correlation coefficients for each organizational model. Then, the observed correlation coefficients (P values) were compared with the expectations of organizational models, and the organizational model with the greatest support was identified. We also used the hierarchical Bayesian method implemented in GESTE version 2.0 [75] to test the effect of geographical distance, human settlement, road, railway and topography on the genetic differentiation of Przewalski's gazelle. This method estimates F ST values for each population and relates them to the factors tested using a generalized linear model, and provides posterior probabilities for each model using a reversible jump Markov Chain Monte Carlo (MCMC) approach. Three independent replicates were executed with default parameter values of the program to ensure consistency of results.

Reliability of genotyping results
We collected 28 skin samples and 182 fresh fecal samples from the field. We amplified 25 skin and 161 fecal samples with more than ten loci. The unbiased P ID of 3.48610 29

HWE and LD
At the 0.05 significant level, all microsatellite loci were in HWE for each population and the entire sample. Thirty seven out of 702 loci pairs for each population (P1, 2; P2, 5; P3, 10; P5, 2; P6, 6; P8, 3; P9, 9) and seven out of 78 loci pairs for the entire sample showed significant LD. After Bonferroni correction, genotypic linkage disequilibrium was not significant.

Population genetic differentiation
The global F ST and F9 ST across the study area was 0.072 (P,0.001) and 0.147 (P,0.001) respectively. Pairwise F ST between sampling populations ranged from 0.0002 to 0.1152, and the range of pairwise F9 ST were from 0.000 to 0.257 (Table 4). From both measures (F ST or F9 ST ), the results indicated significant genetic differentiation between the population pairs except P2-P3 pair, P5-P6 pair, and pairs involving P4 or P7 (Table 4). Withinpopulation genetic diversity (H S ) ranged from 0.462 in P3 to 0.641 in P7 (Table S2). The program Barrier identified five genetic discontinuities (a-e) with high bootstrap support regardless which differentiation measure (F ST or F9 ST ) was used, which separated the nine populations of Przewalski's gazelle into seven groups: P1, P2-P3, P4, P5-P6, P7, P8, P9 (Figure 2). In the Structure analysis, average Ln P(D) maximized at K = 6 genetic groups consisting of: G1 (P1), G2 (P2-P3-P4), G3 (P5-P6), G4 (P7), G5 (P8) and G6 (P9) (Figure 3), whereas the highest value of average DK emerged at K = 3 genetic groups containing: G91 (P1, P5, P6, P7, P8), G92 (P2, P3, P4) and G93 (P9) ( Figure S1). Because the main incongruence was that whether G91 was subdivided, we re-ran STRUCTURE using only the data of G91 to test for further subdivision. The results showed that four genetic groups (the same as G1, G3, G4, G5) did exist within G91 ( Figure S2). Overall, Structure analyses indicated that there were most likely six genetic groups (G1-G6) within Przewalski's gazelle. These genetic groups were the same as AMOVA grouping defined by human settlement. The AMOVA results revealed a significant proportion of genetic variance among groups in all cases, however, the proportion of variance among groups defined by human settlement or genetic structure from Structure (7.68%, P,0.001) and genetic discontinuities from Barrier (7.90%, P,0.001) were much higher ( Table 5). The groupings defined by topography, road and railway also showed significant genetic variance among populations within groups (Table 5).

Dispersal pattern
In the GeneClass analysis, we detected four first generation migrants based on L home /L max . They were Hudong4 (in P2) from P3, Ketu10 (in P3) from P2, Shadao2 (in P4) from P3 and HerG8 (in P6) from P5. Using L home , in addition to the same four F 0 being identified, another first generation migrant (Hudong7) in P2 was detected from P1. Spatial autocorrelation analysis displayed positive and significant r values for the first five distance classes (1, 2, 5, 10 and 20 km), and the highest r value was observed within two kilometers (r was around 0.1, P = 0.001). The correlogram flattened out after 20 km, then r values remained slightly negative up to 80 km, followed by oscillation (Figure 4).

Landscape genetic analyses
In Mantel tests, the results from two genetic distance measures were largely consistent. Significant positive associations were found between genetic differentiation and human settlement, road or topography, but genetic differentiation was not significantly correlated with the presence of railway (Table 6). When using F ST /(12F ST ), Mantel tests also  indicated that 17.04% of the genetic differentiation could be explained by IBD (isolation by distance) ( Figure 5). Regardless of the genetic distance measure, the results of partial Mantel tests showed significant correlation between genetic differentiation and human settlement (Table 6), and AIC c calculations demonstrated that the models concerning human settlement had the lowest AIC c values both in Mantel and partial Mantel tests, thus they were superior to other models (Table 6). From causal modeling analysis with either genetic distance measure, only the isolation by human settlement model was fully supported by all statistical expectations (Table 7 and S3). This model predicts that genetic differentiation is mainly attributable to the effect of human settlement with no significant independent relationships with geographical distance, road or topography. GESTE generated 32 models and gave consistent results across three replicates. The null model which excluded all tested factors was inferred as the best model because it had the highest posterior probability, and the values of other models were very low (Table 8), suggesting little effect of tested factors on observed genetic differentiation. However, the estimate of s 2 for the null model, a measure of model fit, was very high with the upper bound of the HPDI over 1 [0.45, HPDI (0.11; 1.04)], indicating that GESTE failed to identify the true model [75].

Discussion
The F-statistics and F9 ST analyses revealed significant genetic differentiation in the entire sample and between most populations ( Table 4). The results are within our expectation given that the distribution of Przewalski's gazelle has been highly fragmented by anthropogenic landscapes (Figure 1). For population pairs with no significant differentiation, comparisons involving P4 or P7 could be due to its small sample size. P2 and P3 may exchange individuals via barrier-free desert. P5 and P6 are separated by the Qinghai-Tibet railway, yet a bridge with underpasses is present and could allow the migration of gazelles. Although F9 ST presented larger global and pairwise values than F-statistics analysis (Table 4), suggesting a higher level of actual genetic differentiation than indicated by Fstatistics, the differentiation pattern inferred from F9 ST and Fstatistics results were largely consistent. For example, F9 ST and Fstatistics indicated significant genetic differentiation between the same population pairs, and they ranked the pairwise differentiation similarly (Table 4). In the Barrier analysis using either F ST or F9 ST matrices, we found five genetic discontinuities and seven groups ( Figure 2) which are largely in accordance with the observed pattern of genetic differentiation. Genetic structure analysis using Bayesian clustering method detected six genetic groups (Figure 3) nearly identical to the groups inferred by Barrier (except P4), corroborating the results from F-statistics, F9 ST and Barrier analysis. The AMOVA results showed that the grouping patterns according to genetic discontinuities, genetic structure and human settlement explained relatively high proportion of genetic variance among groups and revealed small and nonsignificant among populations within groups (Table 5), thus these groupings seem to be rational for characterizing the distribution pattern of genetic variance. Previous research on Przewalski's gazelle using mtDNA control region sequences revealed a population genetic structure (Yuanzhe -P1 in this study, Hudong-Ketu -P2 and P3 here, Shadao-Gahai -P4 here, and Bird Island -P7 here) consistent with the results of genetic differentiation in this study [26].
GeneClass analysis found five first generation migrants between four pairs of neighbouring populations (P1-P2; P2-P3; P3-P4; P5-P6), indicating low migration rate and restricted dispersal in Przewalski's gazelle. Moreover, the results from GeneClass also corroborated the observed pattern of genetic differentiation, because four F 0 were found between the populations not differentiated. Spatial autocorrelation analysis detected positive autocorrelation at shorter distances as predicted under a restricted dispersal model (Figure 4), suggesting that dispersal events of the gazelle were confined to small geographical scale and most likely to occur between neighbouring populations. The highest autocorrelation index r was found within two kilometers which is approximately equal to the population home range of the gazelle [25], implying that individuals within the same population are closely related. Overall, the results indicated a restricted dispersal pattern of Przewalski's gazelle, congruent with the significant genetic differentiation observed.
Several previous studies have revealed significant genetic effects of anthropogenic landscapes on wild ungulates, such as desert bighorn sheep Ovis canadensis nelsoni [5], European roe deer Capreolus capreolus [76] and Scottish highland red deer Cervus elaphus [77]. In our research, landscape genetic analysis also yielded evidence for significant genetic influences from anthropogenic landscape features on Przewalski's gazelle. Mantel tests found that  (Table 6). Geographical distance was also inferred as a significant factor in Mantel tests when genetic distance was measured as F ST /(12F ST ) ( Figure 5). The results of partial Mantel tests, AIC c calculations and causal modeling analysis with both genetic distance measures further demonstrated that human settlement was the main factor (Table 6, 7 and S3). The AMOVA analyses showed that the grouping defined by human settlement interpreted much higher proportion of genetic variance among groups than those defined by other landscapes (Table 5), thus human settlement appears to be the most important landscape feature among those tested determining the distribution pattern of genetic variance. GESTE analysis did not provide reliable results for identifying the effect of tested factors on observed genetic differentiation, probably because the number of populations in our data was small. Foll and Gaggiotti have suggested that the analytical power of GESTE can be limited by the number of populations analysed (,20 populations), and the method may fail to identify the true model with less than ten populations [75]. The primary role of human settlement could be corroborated by following observations. First, P1, P2 and P3 provide a good example showing the large genetic effect of human settlement (Figure 1). These populations are located on the same side of a road, with a nearly equal interval of 12 km. Several villages separate P1 from P2 and no barrier exists between P2 and P3. In 2this study, we found that P1 and P2 had significant genetic differentiation (Table 4) and the villages between them coincided with the location of inferred genetic discontinuity b (Figure 2), but P2 and P3 presented no differentiation (Table 4). We can see that, even over a short distance of 12 km, human settlements seem to have significantly influenced population genetic differentiation of the species. Second, the people living in this area are nomads. They once hunted gazelles for food in winter, and their livestock compete with the gazelle for food resources [78,79]. Moreover, Przewalski's gazelle are extremely vigilant towards human activity and escape swiftly when people approach them within 500 m [21], and they have shifted their peak foraging time to the mornings and evenings to avoid livestock grazing and human disturbance [20]. Therefore, the physical barrier of human settlement may be amplified by associated human activities such as hunting and livestock grazing. Third, although human population density is very low throughout the study area (four persons per km 2 ), it is much higher within human settlement areas (about 1000 to 1600 persons per km 2 ).
Besides human settlement, attention should also be given to the effect of geographical distance, road and topography. As a restricted dispersal pattern was revealed in Przewalski's gazelle, geographical distance seems to be an important influencing factor for distant populations. The study area is fragmented by roads (Figure 1), and a recent study suggested that roads around Qinghai Lake changed the diurnal activity of the gazelle [80]. Provided that the traffic volume continues to increase along with the development of economy, roads would certainly exert larger influence on the gazelle. Based on our field research experience, the gazelle is unable to traverse the Qinghai Lake or mountain ridge, hence topography appears to play an important role on the distribution of and gene flow between populations of Przewalski's gazelle.
The railway showed no influence on genetic differentiation, probably because of its short history, limited distribution and infrequent traffic (Table 2). Moreover, the bridges and underpasses along the railway are well designed and most likely allow Przewalski's gazelle to safely cross this infrastructure, for example, between P5 and P6. Across the study area, there are other anthropogenic landscape features (e.g. grassland fence) which were not considered here because of very short history (constructed in the last 15 years) and the ability of the gazelle to cross them (e.g. Przewalski's gazelle frequently jump over the grassland fence).
Based on our findings in this study, we propose the following measures to conserve Przewalski's gazelle. First, since human settlement has large effect on Przewalski's gazelle, measures mitigating its negative influence should be carried out. One sensible strategy is collaboration with local community to protect the gazelle, through conservation communication et al. This proposal is probably realistic, given that local people are in favor of protecting the gazelle [81]. Second, according to the identified first generation migrants and the restricted dispersal pattern in the gazelle, suitable habitat for Przewalski's gazelle between the six populations (P1-P6) with short geographical distance should be protected and restored to facilitate the potential migration. Third, besides human settlement, the effect of other anthropogenic landscape features, such as road, railway and grassland fence, should not be neglected and should be minimized through the construction of wildlife underpass and lowering of the fence.
In conclusion, our results indicated significant genetic differentiation, five genetic discontinuities, six genetic groups and restricted dispersal pattern in Przewalski's gazelle, and showed that, under realistic circumstance with compound effect, IBD, topography and some anthropogenic landscape features (human settlements and roads) have been affecting the process of population genetic differentiation. Furthermore, using multiple methods, we demonstrated that human settlement is the main factor shaping current genetic differentiation among those tested. In conservation context, we propose several conservation strategies, such as habitat protection and restoration between close populations, and collaboration with local community. It is possible that other factors, biotic as well as abiotic, may have contributed to the genetic structure of this species, although the analyses presented here seem to point clearly at human settlements being the primary agent responsible for differentiation. To our knowledge, this is the first study assessing the genetic effect of anthropogenic landscape features on wild ungulate in an alpine grassland ecosystem on the Qinghai-Tibet Plateau. Because Przewalski's gazelle is a flagship species [18], our research provides valuable ideas for the management of other wild animals such as larger mammals across the region. Supporting Information