Non-invasive genetic analysis indicates low population connectivity in vulnerable Chinese gorals:concerns for segregated population management

Jing Yang,Guo-Fen Zhu,Jian Jiang,Chang-Lin Xiang,Fu-Li Gao,Wei-Dong Bao,＊

1 College of Biological Sciences and Technology,Beijing Forestry University,Beijing 100083 China

2 College of Biology and Food Science,Hebei Normal University for Nationalities,Chengde Hebei 067000,China

3 Beijing Songshan National Nature Reserve Administration,Beijing 102115,China

4 Inner Mongolia Saihanwula National Nature Reserve Administration,Daban Inner Mongolia 025000,China

ABSTRACT Detailed information on the size and genetic structure of wildlife populations is critical for developing effective conservation strategies,especially for those species that have suffered population decline and fragmentation due to anthropogenic activities. In the present study, we used a non-invasive approach combining fecal pellet sampling with mitochondrial DNA and nuclear DNA microsatellite marker analysis to monitor and compare the population structure of the Chinese goral (Naemorhedus griseus) in Beijing and northeast Inner Mongolia in China. Of the 307 fecal samples confirmed to be from N. griseus, 15 individuals (nine females and six males) were found in the Beijing population and 61 individuals (37 females and 24 males) were found in the Inner Mongolian population. Among these 76 individuals,we identified eight haplotypes and 13 nucleotide polymorphic sites from mtDNA and 45 alleles from 10 microsatellite loci. Spatially structured genetic variation and a significant level of genetic differentiation were observed between the two populations. In both populations, the sex ratios were skewed toward females, indicating high reproductive potential, which is crucial for population recovery and conservation of this patchily distributed vulnerable species. We suggest that managing the two populations as evolutionarily significant units with diverse genetic backgrounds could be an effective solution for present population recovery,with the possible relocation of individuals among different groups to help ensure future goral species prosperity.

Keywords: Genetic diversity; Mitochondrial DNA;Microsatellite marker; Naemorhedus griseus; Noninvasive sampling;Sex ratio

INTRODUCTION

Determining the size, genetic diversity, and connectivity among populations is important for understanding the ability of target species to adapt to different habitat conditions (Banks &Lindenmayer, 2014; Dharmarajan et al., 2013) and for monitoring and developing conservation strategies for endangered species (Lande, 1988). Among the many factors affecting wildlife population dynamics, such as life history traits, spatial behaviors, predation pressure, human disturbance, and food resource abundance, the sex ratio is a fundamental biological index that can indicate population growth (Székely et al., 2014) and status (Lindstr?m et al.,2002). For example, a greater number of females than males can imply that a population is increasing rather than decreasing (Cameron, 2004). However, determining the sex of individuals is often difficult in the field because heterogeneous habitats may block an observer's line of sight and target animals may hide or lack obvious external sexspecific characteristics (Lewis et al., 2002). Advances in molecular biology have allowed the use of animal feces, hair,eggshells, and food residue as non-invasive sources of DNA,which can be used in conservation genetics research(DeSalle & Amato, 2004). For example, fecal samples have been used in studies on rare and endangered animals involving species identification, sex determination, individual identification, effective population size estimation, and population genetic structure analysis (Haag et al., 2009;Mondol et al.,2009).

The Chinese goral (Naemorhedus griseus), a goat-like wild ungulate of similar body size, is mainly distributed in central and northern China (Hrabina, 2015; Smith & Xie, 2008; Mead,1989). Excessive medicinal exploitation and habitat fragmentation have both resulted in rapid declines and isolation of goral populations throughout most of its range(State Forestry Administration, 2009). Currently, the Chinese goral is listed as a Class II state key protected wildlife species in China and has also been listed as Vulnerable (VU) on the International Union for Conservation of Nature Red List(Duckworth et al., 2008) and under Appendix I protection by the Convention on International Trade in Endangered Species of Wild Fauna and Flora. However, no accurate estimates currently exist for Chinese goral populations aside from a previous field observation related to population size and habitat selection (Chen et al., 1999). A systematic national survey on this species indicated that habitat fragmentation has forced the remaining populations to retreat in small numbers to patchy areas in nature reserves, resulting in severe isolation and increased extinction risk of local populations (State Forestry Administration, 2009).Furthermore, little is known regarding the genetic structure or molecular biology of these small-group goral populations living under the threat of habitat fragmentation, despite the importance of such knowledge in planning effective conservation strategies.

We previously conducted field surveys on wild animal diversity in the mountain ranges around Beijing and at a national nature reserve in Inner Mongolia, which revealed that goral populations were experiencing slow recovery. This recovery highlighted the need to investigate the biological traits of the populations for efficient conservation of this endangered species. Furthermore, the Inner Mongolian population was separated from all other geographical groups and may exhibit a different genetic background. Thus, we studied the structure and genetic connectivity of the goral populations to detect whether biological differences exist between the two groups, which may have valuable implications for isolated population management. In this study,mitochondrial DNA(mtDNA) gene fragments and nuclear DNA microsatellite markers combined with non-invasive fecal sampling were used to examine the minimum population size,sex ratio, and genetic connectivity of two geographically separated populations of Chinese goral in the Beijing Songshan National Nature Reserve (BS) and Inner Mongolian Saihanwula National Nature Reserve (IMS). The BS population, which inhabits the mountain ranges around Beijing and is connected with nearby populations in Hebei Province,experienced severe population decline as a consequence of historical hunting for food and medicine before the promulgation of the Wildlife Protection Law in 1989. The present population is slowly recovering due to effective local government management of illegal hunting and livestock grazing (Chen et al., 2002). Several recent camera trapping surveys confirmed its re-appearance (Lan & Jin, 2016; Liu et al., 2012) and daily activity patterns in the mountain ranges of Beijing (Zhang et al., 2017). The IMS goral population inhabits the southern segments of Daxing'anling (Greater Khingan Mountain) of Inner Mongolia and is isolated in a national nature reserve surrounded by human settlements. Several ecological studies have been conducted by our research team on this population, including on food selection, sexual hormone variation, and habitat suitability evaluation (Tang et al.,2018,2019).

The aims of the current study were to: (1) assess monitoring protocols, especially regarding sex structure and individual identification, by non-invasive genetic sampling to predict future population dynamics and guide development of conservation strategies; (2) evaluate the baseline genetic status of the two geographically separated populations for efficient management; and (3) acquire in-depth information in terms of restricted distribution area and population isolation for the IMS population to understand its genetic development trends under limited habitat conditions.

MATERIALS AND METHODS

Study areas and goral populations

This study was undertaken in the Beijing Songshan (BS)National Nature Reserve and Inner Mongolian Saihanwula(IMS) National Nature Reserve in China (Figure 1). The BS reserve is located in the Yanqing district, northwest of Beijing City (E115°43'-115°50', N40°29'-40°33'), and covers an area of 6 212 ha, mainly to preserve natural pine (Pinus tabulaeformis) forest, secondary deciduous broad-leaved forest, and various types of wild animals and plants, including the Chinese goral. The goral populations in Beijing are very patchy, with individuals found almost everywhere in mountain areas (Chen et al., 2002), but in small-sized groups (Liu et al., 2012; Zhang et al., 2017). The regional population is connected with southern conspecific populations in Hebei and Shanxi provinces through the Taihang Mountains (State Forestry Administration, 2009). The IMS reserve is located north of Balinyouqi Banner in Chifeng City (E118°18'-118°55',N43°59'-44°27') and covers an area of 100 400 ha. As a forest nature reserve, IMS protects a rich ecosystem and numerous endangered wild animals, including red deer(Cervus elaphus), Eurasian lynx (Lynx lynx), gray wolf (Canis lupus), great bustard (Otis tarda), and black stork (Ciconia nigra). The Chinese goral population is the only remaining group in the Greater Khingan Mountains, and is in a state of severe isolation, with no viable populations nearby. Historical corridors likely included the southern Greater Khingan Mountains connected with the Yanshan Mountains in the south (including Beijing) (Hrabina, 2015; Mead, 1989).However, these dispersal pathways became disconnected by human settlements on the plains between the mountains over centuries of anthropogenic movement (State Forestry Administration,2009).

Figure 1 Locations of Beijing Songshan and Inner Mongolian Saihanwula National Nature Reserves in China

Sampling procedure

In November 2013, Chinese goral habitat survey and fecal sample collection were undertaken by reserve staff with the permission of the Departments of Natural Resources Management of the BS and IMS administrations. The sampling teams searched for goral pellets randomly in major habitats based on previous field observations.All fresh pellets were collected without considering outer shapes or distances between the pellet stacks to avoid sex bias. Only fresh (based upon color and moisture) fecal samples were taken. Newly purchased Ziploc bags and tools were used for sample collection to avoid possible DNA cross-contamination and all samples were stored in anhydrous ethanol at -20°C within 6 h of collection. The collection period lasted 10 d from 5 and 10 November in BS and IMS, respectively. GPS positions and collection times were recorded for each fecal sample. In addition, a muscle tissue sample from a dead adult female in IMS and a blood sample from a rescued male in Beijing were used as positive controls.

DNA extraction

We extracted DNA from the muscle and blood samples using a TIANamp Genomic DNA Kit (Tiangen Biotech Co. Ltd.,Beijing, China) and from fecal samples using an EZNA Stool DNA Kit (Omega, Dorivalle, GA, USA) following the manufacturer's protocols, with some modifications (Zhang et al.,2009)to improve the amount and quality of DNA extracted.

Species identification based on mtDNA

Two mtDNA fragments were used to identify species and analyze population genetic structure, respectively. Species were determined from the fecal samples via PCR using a partial mitochondrial cytochrome-b gene with goral speciesspecific primers, as per Kim et al. (2010) (i.e., BJGL-F1: 5'-CTGCCTAATTCTACAAATCCTA-3' and BJGL-R: 5'-TAGGAG GATTACCCCAATA-3'). Each experiment included a positive and negative control (DNase/Rnase-free deionized water template rather than DNA). Specific primers (i.e., DF249 5'-ACCCCATGCATATAAGCATG-3' and DR496 5'-GGATACGC ATGTTGACAAGG-3') were used to amplify a 227 bp segment of the mtDNA control region and were designed using Primer Premier 5 (Lalitha, 2000) in reference to the published mtDNA control region from the South Korean long-tailed goral(Naemorhedus caudatus) (GenBank accession Nos.:EU259152, EU259154, EU259155, EU259165, EU259169,EU259171, and EU259176) (Kim et al., 2008, 2010). These segments were used to analyze the basic genetic structure of the two goral populations. The PCR conditions are shown in Supplementary Materials 1A. Both PCR products were isolated by electrophoresis in a 2% (w/v) agarose gel, and the control region PCR products were then further purified using a Wizard SV Gel and PCR Clean-Up System (Promega,Madison, WI, USA) before final sequencing by SinoGenoMax Co.Ltd.(Beijing,China).

Individual identification based on microsatellite loci

Ten microsatellite loci (SY12A, SY12B, SY48, SY58, SY71,SY76, SY84, SY84B, SY112, and SY129) were used for individual identification (An et al., 2005, 2010). Fluorescence labels (HEX or FAM) were added to the 5' terminus of the forward primers of all microsatellite loci. As fecal DNA in this analysis was very dilute, a two-step PCR procedure was used to increase the amplification rate and reduce the genotyping error rate (Piggott et al., 2004; Taberlet et al., 1999) (see Supplementary Materials 1B for PCR conditions). To guard against contamination, 75% alcohol was sprayed over the experimental table before and after setting up each PCR run.Each DNA sample underwent a minimum of five independent PCR runs for each microsatellite locus to reduce the effects of possible allelic dropout or false alleles (Broquet et al., 2007;Mondol et al.,2009;Taberlet et al.,1996).

The PCR products were assayed by electrophoresis on 2.5% agarose gels, and the amplified target fragments were delivered to SinoGenoMax Co., Ltd. (Beijing, China) for genotyping analysis using an ABI-PRISM 310 DNA Analyzer(Applied Biosystems, USA). Similarly, fluorescence of different DNA sizes of microsatellite alleles was determined using GeneMarker v2.4.0 (SoftGenetics LLC., State College, PA,USA)to confirm the genotyping results.

Sex discrimination

Nested PCR based on amplification of a Zn-finger gene on the X and Y chromosomes was used for sex identification, as successfully applied for Korean long-tailed goral (Kim et al.,2008, 2010). Detailed PCR procedures are provided in the Supplementary Materials 1C.

The PCR products were assayed by electrophoresis on 2.5% agarose gels to confirm sex. If the X and Y alleles were both detected or only the Y allele was detected in a sample,sex was determined to be male; if only the X allele was detected, sex was determined to be female. Samples without amplification products were considered failures. Each sample underwent a minimum of five independent PCR runs, and at least three consistent results were required for successful determination to reduce possible error from allelic loss(Scandura,2005).

Genetic data analysis

The obtained mtDNA control region sequences were analyzed with Clustal W v2.0 for sequence alignment (Larkin et al.,2007) and were proofread by eye. DnaSP v5.0 (Librado &Rozas, 2009) was used to identify haplotype frequency and number of polymorphic sites, haplotype diversity (h),nucleotide diversity (π), average number of nucleotide differences (k), coefficient of genetic differentiation between populations (FST), and significance level of differences. Gene flow (Nm) was calculated by MEGA5 software via FSTand was used to analyze nucleotide composition and calculate genetic distances between haplotypes based on the Kimura double parametric method (Tamura et al., 2011).The haplotypes were compared with the Korean long-tailed goral homologous sequences (GenBank accession Nos.: EU259152-EU259176)to evaluate phylogenetic relationships. Arlequin v3.5(Excoffier & Lischer, 2010)was used to determine Fu's Fs (Fu,1997) and Tajima's D (Tajima, 1989) to assess whether nucleotide polymorphisms deviated from expectations under neutral theory.

The neighbor-joining method in MEGA5 was employed to reconstruct phylogenetic relationships among partial mtDNA control region haplotypes. To acquire the phylogenetic tree root, homologous sequences of sheep (Ovis aries) KJ954145 were used as an outgroup. The confidence coefficient of the systematic tree was tested via bootstrap resampling with 1 000 simulated samples. In addition, the maximum parsimony median-joining network, which can effectively illustrate the relationships among sequences with low diversity, was used to further visualize genealogical relatedness among haplotypes with NETWORK v5.0.0.1(Bandelt et al.,1999).

Micro-Checker software (Van Oosterhout et al., 2004)was applied to assess the null alleles of each microsatellite locus and to identify reliable genotyping for each sample. The frequencies of the null alleles were obtained using the Brookfield method (Van Oosterhout et al., 2004) employed in Micro-Checker, and the adjusted allele frequencies were then used to recompute heterozygosity (HEnull) and observed heterozygosity (HOnull) values in GenAlEx v6.5 (Peakall &Smous, 2010). Samples amplified with less than five loci were judged as failures and were excluded from subsequent analyses (Mondol et al., 2009). Arlequin v3.5 was used to examine the Hardy-Weinberg equilibrium and linkage disequilibrium (LD) and to calculate LD for all microsatellite loci, intra-population genetic differentiation coefficients (FST),and levels of significance. If a significant presence of null alleles was detected, we recalculated the FSTwith corrected genotype frequencies in GenAlEx v6.5. GenAlEx v6.5 and Cervus v3.0 (Kalinowski et al., 2007) were used to calculate allele number (Na), number of effective alleles (Ne), observed heterozygosity (Ho), expected heterozygosity (He), and polymorphism information content (PIC) of each microsatellite locus.

RESULTS

Species identification, individual identification, and sex determination

In total, 40 and 332 fresh fecal samples were collected from BS and IMS, respectively. Among these, the partial mitochondrial cytochrome b gene was successfully amplified in 33 samples from BS and 274 samples from IMS, which were thus identified as Chinese goral based on adequate DNA, with a species identification rate of 82.50% and 82.53%,respectively.

Genotyping results derived from the 10 pairs of microsatellite primers recognized 76 individuals from the 307 tested samples after removal of failed amplifications and duplicates, with 15 and 61 individuals recognized from the BS and IMS populations, respectively. The cumulative probabilities of two unrelated individuals sharing the same multi-locus (PID) for all loci estimated in Gimlet software were 0.000 04 and 0.000 01 for the BS and IMS populations,respectively. The corresponding probabilities of two full siblings (PIDsib) were 0.000 6 and 0.000 3, proving that the 10 loci used here were sufficient for accurately identifying individuals.

In total, 15 and 151 samples were identified as female and 14 and 92 were identified as male in BS and IMS,respectively. The sex identification success rates were 87.88% and 88.69% for the two groups of identified goral samples. After correcting for repeated individuals, we identified nine females and six males in the BS population,with a female to male ratio of 1.50: 1, which did not significantly deviate from 1: 1 (χ2=2.77, P＞0.05), and 37 females and 24 males in the IMS population, with a female to male ratio of 1.54:1, which also did not significantly deviate from 1:1(χ2=0.6,P＞0.4).

Genetic diversity of mtDNA control region

After individual identification, we obtained 76 partial mitochondrial control region sequences from the two populations. Thirteen polymorphic sites were detected from the 227 bp gene fragments by sequence alignment, with a rate of 5.73% for the mtDNA sequences. Within the 13 polymorphic sites, 13 transitions, three transversions, and no insertion/deletions were detected. All mutation sites were analyzed with Tajima's neutrality test, with an obtained Tajima's D value of 1.48(P＞0.10).

Eight control region haplotypes were identified, including three (BJ-A to BJ-C) from BS and five (NM-A to NM-E) from IMS. We found no shared haplotypes between the two geographic populations. The nucleotide sequences of BJ-A to BJ-C and NM-A to NM-E for the eight haplotypes were submitted to GenBank under accession Nos. KP330565 -KP330572.

Genetic parameters such as haplotype diversity (h),nucleotide diversity (π), and average number of nucleotide differences (k) indicated moderate genetic diversity for the two goral populations(Table 1).

Table 1 Genetic diversity parameters and demographic analyses of partial mtDNA control region from two goral populations in China

Genetic structure of microsatellite loci

We screened 10 pairs of polymorphic microsatellite primers to determine the genotypes of the populations. Among the identified Chinese goral feces, 26 samples from BS and 212 from IMS provided reliable and complete genotypes at the 10 loci, with genotyping success rates for BS and IMS of 77.79%and 77.37%, respectively. Forty-five alleles were detected in 10 polymorphic microsatellite loci, 25 were shared between the two populations, eight were exclusive to BS, and 12 were exclusive to IMS. The distribution frequencies of the 45 alleles in the two populations were not identical. The genetic diversities calculated by the microsatellite dataset are shown in Table 2. In the IMS population, the Micro-Checker results,with a 95% confidence interval, indicated that null alleles existed at three of the 10 loci (i. e., SY12B, SY71, and SY84B). The frequency of null alleles differed for the three microsatellite loci (i.e., 7.68% for SY12B, 7.44% for SY71, and 10.16% for SY84B). Thus, the values of observed heterozygosity for the IMS population slightly increased when the dataset was adjusted based on the null allele frequencies,but there was no change in the expected heterozygosity(Supplementary Materials 2). No evidence for genotyping misjudgment caused by null alleles, strip discontinuity, or missing alleles was found within the 10 microsatellite loci in the BS population.

Table 2 Genetic diversity among 10 microsatellite loci observed in the two goral populations

Results further indicated that the SY12B (P=0.004) and SY84B (P=0.024) loci from the IMS goral population deviated clearly from the Hardy-Weinberg equilibrium, suggesting independent evolution of this population, whereas the other loci showed no deviation from the two populations. After 45 linkage disequilibrium tests with the 10 loci, SY84B and SY112 from the BS population (χ2=8.42, P＜0.05) and SY58 and SY84 from the IMS population showed linkage disequilibrium (χ2=7.118, P＜0.05), with no linkage disequilibrium found in other loci. After sequential Bonferroni correction for pairwise comparisons, the differences in linkage disequilibrium loci were not significant, indicating that all loci were inherited independently and could be an accurate reflection of the genetic diversity status in the studied goral populations.

Genetic differentiation of haplotypes for two populations

The genetic distances between the eight mtDNA control region haplotypes ranged from 0.004 4 to 0.060 3, with an average of 0.031 8. The average genetic distances between haplotypes within the BS and IMS populations were 0.014 9 and 0.012 5, respectively. Maximum and minimum genetic distances were 0.022 4 and 0.008 9 for the BS population and 0.022 5 and 0.004 4 for the IMS population,respectively.

The inter-population genetic differentiation coefficient was estimated at FST=0.774 (P＜0.001), suggesting a significant level of genetic differentiation between the two populations.The gene flow value for the two populations was Nm=0.07,showing that the inter-population genetic distance was large and that the number of migrating individuals was only 0.07 for each generation.These results thus revealed a lack of genetic exchange between the two populations. The neighbor-joining phylogenetic tree showed that the haplotypes were clustered into two large clades, the BS and IMS goral populations were clustered into one clade and the South Korea populations were clustered into the other clade (Figure 2A). The Chinese goral populations were further separated into two small independent branches, showing distinct genetic differentiation and a lack of gene flow. The median-joining network also grouped the haplotypes into a same topology, with several mutation steps existing between the two Chinese populations(Figure 2B).

Figure 2 Integrative genetic analyses of mtDNA control region haplotypes in Chinese goral

Phylogenetic discrepancy between populations

The FSTgenetic differentiation coefficient and Nei's genetic distance for the BS and IMS populations were 0.104 and 0.336, respectively. Considering the existence of null alleles,the recomputed FSTand Nei's genetic distances were 0.065 and 0.334, respectively, with adjusted genotypes. These results indicated that the two populations had a moderate level of genetic differentiation. Cluster analysis based on STRUCTURE software showed that the BS and IMS populations were two separate clades with no interaction,indicating that genetic differentiation was significant, and the phylogenetic relationship was distant for the two populations(Figure 3).

DISCUSSION

Estimating the size of a wild animal population is central for managing game species and conserving rare and endangered animals.The Chinese goral is typically very alert and sensitive to disturbance. Furthermore, its fur is often indistinguishable from the surrounding rocks on the mountaintops it inhabits,which has caused great difficulties in traditional transect surveys. However, appropriate techniques to extract and analyze DNA from feces can help identify individuals in the field and avoid bias from repeated counting of a single individual (Haag et al., 2009). In the current study, we identified at least 61 animals in the 100 km2sampling area in IMS and 15 animals in a 30 km2area in BS. These results,combined with the female-biased sex ratios, suggest a high recovery potential for these separate goral populations;however, further research into the genetic variations of their populations is needed to reveal the underlying mechanism.

Figure 3 STRUCTURE clustering based on microsatellite genotypes of Chinese goral populations

Sex ratio analysis can be used to evaluate the structure and future development of an animal population (Kruuk et al.,1999); it can also provide basic data for planning and implementing effective management strategies to conserve rare and endangered wildlife (Bradshaw et al., 2003;Majumder et al., 2017). In this study, a sex-linked Zn-finger gene was used for sex structure analysis. Results showed a similar female-biased sex ratio in both populations, i.e., 1.54:1(female to male) in IMS and 1.50:1 in BS.Although there was no statistical difference from the theoretical 1:1 ratio, the high female bias in the sex structure of the wild goral populations is expected to promote population recovery. A number of hypotheses have been proposed to explain adaptive sex-ratio variation in mammals (Cockburn et al., 2002). One commonly applied model for ungulates is the local resource competition(LRC) hypothesis, which states that competition for environmental resources may shape the primary sex ratio,and predicts that mothers with poor body condition should reproduce the sex that is more likely to disperse and consequently reduce local competition for resources between mothers and their offspring within the natal area (Clark, 1978).In our study areas, food resources are abundant (Li, 2005; Liu et al., 2012) and the females are in good condition (as confirmed by field camera trapping photos). However, the bare rock microhabitats to which gorals retreat when facing mammalian predators are limited. Thus, the high proportion of females within the populations may indicate a demographic variation due to the limited microhabitat resource pressure.Our ongoing GPS tracking captured six gorals at the IMS reserve, which were all females, thus confirming their high proportion. Furthermore, the preliminary home range results showed that sub-adult females (n=3) shared their small home ranges with adult females. Thus, the high female to male ratios in our study fit the LRC model prediction, in which healthy mothers reproduce more philopatric daughters than dispersal sons in environments with abundant food resources and limited bare rock microhabitats. Although our results found a greater number of females in the two Chinese goral populations, thus indicating high reproductive potential and population growth, further in-depth field monitoring on the social structures and breeding behaviors is needed to reveal actual population growth trends, especially for the isolated population at IMS.

Average genetic distance and nucleotide diversity between haplotypes serve as important indicators of mtDNA genetic diversity. In mammals, a population is considered genetically diverse if the mean genetic distance between haplotypes is greater than 0.01, with higher values indicating greater genetic diversity (Nei, 1987). The overall genetic distances among haplotypes in the BS and IMS populations were 0.014 9 and 0.012 5, respectively, indicating rich genetic diversity in the two goral ranges. Compared to the average genetic distance among haplotypes, nucleotide diversity π can better reflect mtDNA diversity in a population (Nei, 1987). Nucleotide diversity in mammals can range from 0.003 to 0.04 with a maximum value of 0.1, with a higher π value indicating greater genetic diversity (Neigel & Avise, 1993). In this study,nucleotide diversities were 0.009 7 and 0.010 7 for the BS and IMS populations, respectively, indicating moderate genetic diversity. Compared to other Caprinae subfamily members,such as the Helan Mountain bharal (Pseudois nayaur) with π=0.003 92 (Wang et al., 2006), red goral (Naemorhedus baileyi)with π=0.000 88 (Xiong et al., 2013), and serow (Capricornis sumatraensis) with π=0.024 9 (Liu et al., 2013), the goral populations in our study showed an intermediate level of nucleotide diversity.

Haplotype diversity can reflect recent changes in a population. In this study, the mtDNA control region showed a high level of haplotype diversity and moderate level of nucleotide diversity in the two geographically separated populations. Thus, it is reasonable to speculate that the two populations may have arisen from rapid growth of groups with small effective population sizes, similar to the reported population dynamics of the takin (Budorcas taxicolor) with high haplotype diversity and low nucleotide diversity (Li et al.,2003). Fossils of Naemorhedus species have been found in the middle and lower Pleistocene in China, upper Paleolithic deposits in Far East Russia, and the third and fourth interglacial periods in Korea (Mead, 1989), areas which are partially covered by the present-day distribution ranges,indicating that Naemorhedus experienced a historical distribution range recovery after the quaternary glacial.

The existence of null alleles is one of the biggest drawbacks for microsatellite markers as they can generally increase FSTand genetic distance values with decreasing within-population genetic diversity (Paetkau et al., 1997; Slatkin, 1995). Herein,among the 10 microsatellite markers in the two studied populations, three loci (i.e., SY12B, SY71, and SY84B) (see Supplementary Table S1) from the IMS population showed the possible existence of null alleles, which could be attributed to the low quality and quantity of fecal DNA (Waits et al., 2001).After adjustment of allele frequencies in these loci, the observed heterozygosity (HO) for the three single markers increased remarkably, and the FSTdecreased from 0.104 to 0.065 between the two populations. This reduction may be slightly overestimated due to the evaluation method based on prior probabilities; however, the correction results still need to be carefully considered.

The estimated genetic differentiation between the two studied goral populations was FST=0.774 (in mitochondrial genes), which may reflect different demographic histories. The IMS population is the only remaining reproductively viable population in the Greater Khingan Mountains (State Forestry Administration, 2009), with severe isolation and habitat fragmentation from anthropogenic activities accelerating its independent evolution. Therefore, strengthening future monitoring would be an effective way to surveil the dynamics of genetic diversity on this isolated goral population for the purpose of maintaining long-term population viability.

The population size and sex ratio results provide baseline data for implementing effective population conservation strategies. The techniques used in this study, which include non-invasive fecal sampling followed by DNA extraction and analysis of population structure, are simple and accurate, with a relatively high success rate compared to other studies on herbivores (Broquet et al., 2007). Our research was performed after the first winter snowfall with an average air temperature below 0 °C and the fecal samples were fixed in anhydrous alcohol within several hours of collection to delay DNA degradation. Moreover, a two-step PCR procedure was used to increase the PCR amplification products and reduce the genotyping error rate. The relatively high identification rate confirmed the validity and reliability of the genotyping procedure, which should provide important guidance for future non-invasive genetic monitoring on goral population structures.

The mtDNA and microsatellite diversity analyses both showed significant levels of genetic differentiation between the BS and IMS goral populations. This implies that the distribution range of historically intact goral populations was segregated by expansion of human settlements and disappearance of connecting groups, leading to population isolation (State Forestry Administration, 2009). Habitat fragmentation and population segregation are severe factors threatening the effective conservation of many endangered wild animals throughout the world (Keil et al., 2015) and increasing the risk of extinction of isolated populations.Although it is difficult to assess whether the independent evolution observed in the current study will have positive or negative effects on the long-term existence of populations ,we would recommend that the two Chinese goral populations be preserved and monitored separately as evolutionarily significant units to maintain species genetic diversity within different management areas. If needed, based on further genetic monitoring, the translocation of individuals among geographically separated populations may improve the genetic diversity of the small segregated groups.

SUPPLEMENTARY DATA

Supplementary data to this article can be found online.

COMPETING INTERESTS

The authors declare that they have no competing interests.

AUTHORS'CONTRIBUTIONS

J.J., C.L.X., and W.D.B. drafted the experiments; J.Y. and G.F.Z. performed the experiments and wrote the paper; F.L.G. analyzed the data.All authors read and approved the final version of the manuscript.

ACKNOWLEDGEMENTS

We thank Ji-Gui Wu, Baart, and the staff of the two nature reserves for collecting the goral fecal samples. Grateful thanks are given to Alice Laguardia from BJFU for valuable comments and editing of the original manuscript. We thank NPG Language Editing for proof-reading the manuscript.