
Citation: | Pengcheng Wang, Ping Hu, Jinping Zhang, Lixia Zhang, Jing Zhang, Zhengwang Zhang. 2023: Using non-destructive sampling to evaluate the population genomic status of captive Brown Eared Pheasants. Avian Research, 14(1): 100078. DOI: 10.1016/j.avrs.2023.100078 |
Evaluating the genetic status of threatened species is an essential task in conservation genetics. However, the genetic status of threatened species has been mostly evaluated through techniques that fail to estimate genetic diversity at the whole genomic level. Next generation sequencing can meet this demand, but high quality samples such as blood or muscle tissues are required. However, it is difficult to collect such samples from threatened species because sampling work may impact their health. Therefore, it is essential to design a workflow to evaluate the whole genomic status of threatened species using non-destructive sampling. Even though non-destructive sampling has been used in traditional barcoding technique, the barcoding technique cannot evaluate the whole genomic status. Brown Eared Pheasant (Crossoptilon mantchuricum) is an endangered species, with captive populations maintained in Taiyuan Zoo, China, and Europe. However, the genetic diversity, inbreeding pattern, and mutation load of these two populations are unclear. To uncover the genetic status of these two captive populations, we applied 2b-RAD technology to evaluate the genomic status of these populations using feathers as samples. The feathers could be collected by non-destructive sampling. The results indicate that the Taiyuan Zoo population has a lower genetic diversity and higher inbreeding coefficient than the European population. The Taiyuan Zoo population has lethal mutations when homozygous. The current project uses a non-destructive sampling technique to evaluate the whole genomic status of the two captive populations, providing a paradigm for conservation genetics, which will facilitate the development of conservation biology.
To protect threatened species from global climate change and disturbances resulting from human activities, it is necessary to carry out conservation biology studies (Barnosky et al., 2011; Hohenlohe et al., 2021). Evaluating the quality of gene pools is an important task in conservation biology because it helps us to identify priority conservation units (Vandersteen Tymchuk et al., 2010; Nadyeina et al., 2014). The gene pool always refers to the total of the all genetic variations in a population. Generally, the gene pool quality of populations with high genetic diversity is relatively high. In addition, analyzing the gene pool quality could facilitate conservation workers to identify and evaluate the source population and fitness, respectively, of the reintroduced population for the purpose of their conservation (Jeremie et al., 2013; Mueller et al., 2022). However, it is difficult to evaluate the genetic quality of the populations accurately and comprehensively with limited genetic markers because the populations generally harbor numerous variations.
Allozymes (Li et al., 2010), microsatellites (Gu et al., 2013), mitochondrial and nuclear loci (Wang et al., 2017) have been used in conservation genetic researches. For example, mitochondrial DNA has been used to evaluate the global genetic diversity pattern of terrestrial amphibians and mammals (Miraldo et al., 2016). Microsatellites have been used to identify the conservation units of the Blue Eared Pheasant (Crossoptilon auritum; Gu et al., 2013). However, conservation units sometimes cannot be accurately identified using only a limited set of genetic markers (Zhao et al., 2013). For example, it is suggested that Giant Panda (Ailuropoda melanoleuca) has two management units based on mitochondrial and microsatellite loci (Zhang et al., 2007), but the whole genome data suggest there are three distinct lineages (Zhao et al., 2013). The whole genome data could not only tell us the genetic diversity across whole genome but also provide the detailed information of demographic history (Zhao et al., 2013; Xue et al., 2015). What is more, it could help researchers to investigate the inbreeding coefficient (Xue et al., 2015; Wang et al., 2020), mutation load (Xue et al., 2015; Wang et al., 2020), and predict the trend of the population (Gu et al., 2021; Chen et al., 2022). It is necessary to use whole-genome data in conservation genetics research.
With the development of genomic sequencing technology, researchers could use population genomic data to evaluate the genetic diversity of each population and to identify the conservation units (Zhao et al., 2013; Wang et al., 2020). For example, the population genomic data of the Brown Eared Pheasant (C. mantchuricum) indicated that the species has a relatively low genetic diversity and has three conservation units (Wang et al., 2020). The samples used for these kinds of studies are usually blood or animal remains from the field. However, collecting blood from threatened species not only threatens their survival, but is also legally restrictions in some countries. For example, capturing and collecting blood will cause stress reaction of animals. Moreover, it is difficult to find animal remains in the wild. Hence, it is necessary to design a framework that employs whole-genome sequencing to evaluate the genetic diversity of threatened species using non-destructive sampling. Non-destructive sampling has been widely used in previous studies, but they did not evaluate the genetic status at the whole genomic level (Wasko et al., 2003; Schmaltz et al., 2006).
The 2b-RAD technology could facilitate the development of conservation genetics based on non-destructive sampling. It is a restriction site-associated DNA (RAD) genotyping method based on type IIB restriction enzymes and is suitable for high-throughput genotyping (Wang et al., 2012). This method has the advantages of high technical repeatability and cost-effectiveness (Wang et al., 2012). Moreover, trace DNA (100–200 ng) could be used for constructing 2b-RAD libraries and identifying single-nucleotide polymorphisms (SNPs; Cui et al., 2015). The whole genome sequencing DNA libraries always need complete DNA sequences, such as the DNA length is larger than 10 kb. But, the samples that cannot provide complete DNA sequence, such as feathers and feces, may be used to construct 2b-RAD DNA libraries. Furthermore, these types of samples can be collected in a non-destructive manner. Although the 2b-RAD technology has been widely used in population genomics (Gao et al., 2022; Liu et al., 2022; Ruocco et al., 2022), its application in conservation genetic studies is scarce. The RAD sequencing technology has been widely used in conservation genetic researches (Lasalle et al., 2022; Xie et al., 2022), but few studies applied 2b-RAD sequencing technology, which may promote the development of conservation genetics.
The Brown Eared Pheasant (C. mantchuricum) is a threatened species that is classified as vulnerable (VU) in the IUCN red list (BirdLife International, 2016). It is a native species in China, with narrow and fragmented habitats (Zheng, 2015). It is usually found in mountain forests in the temperate zone belt at 1500–1800 m above sea level. Its current distribution is limited to parts of Shaanxi, Shanxi, and Hebei provinces in China, and Beijing, the capital of China. It has been reported that the Brown Eared Pheasant has a low genetic diversity (Wang et al., 2020). Organizing the reintroduction effort may aid in the recovery of this endangered species. In this regard, it is necessary to identify the population source and evaluate the genetic diversity of such species to successfully implement their reintroduction programs. Currently, there are several captive populations of Brown Eared Pheasant. In China, the Taiyuan Zoo has the biggest captive populations of Brown Eared Pheasant, which are raised in one feedlot. The Brown Eared Pheasants were raised in different farms in Europe. These populations may be potential gene pools of this species because they have considerable population size. However, the genetic diversity, inbreeding pattern, and mutation load of these populations are unclear.
The current study employed the 2b-RAD technique to sequence the genomic data from feather samples and used it to: 1) calculate the genetic diversity; 2) estimate the kinship pattern; and 3) evaluate the mutation load of this species. All these analyses will help us understand the genetic status of the captive Brown Eared Pheasant. Moreover, the current project will provide a framework to evaluate the whole-genomic genetic status of threatened species using non-destructive sampling. The framework would improve the development of conservation genetics.
To compare the genomic status of Taiyuan Zoo and European populations of captive Brown Eared Pheasants, we collected feathers (primaries and tail feathers) from 66 and 76 individuals of the populations, respectively (Appendix Table S1). There were three feathers pulled from living individuals and stored at room temperature. The root of feather was used to extract DNA using the TIANamp Micro DNA Kit (DP316) of TIANGEN BIOTECH (Beijing) Co., LTD. The DNA extraction process followed the protocol of the DNA Kit. The extraction of DNA, construction and sequencing of 2b-RAD DNA library was performed by Shanghai OE Biotech Co., Ltd (Shanghai, China). The concentration and amount of DNA was measured after extracting DNA. The samples that have more than 100 ng DNA were used to construct the 2b-RAD DNA library. A type IIB restriction enzyme (BsaXI) was used to construct the libraries. An average of 14,674,497 (95% confidence interval (CI) = 13,908,242–15,440,752) clean reads per sample was aligned against the autosomes of the Brown Eared Pheasant (accession no. GWHAOPW00000000 in CNCB-NGDC) using Bowtie v. 1.2.3 (Langmead and Salzberg, 2012). Samtools v. 1.14 was used to calculate the mapping rate and coverage of each non-zero depth position (Li, 2011). The gstacks function in Stacks v. 2.41 was used to identify the SNPs (Rochette et al., 2019). The SNPs were then filtered using the following criteria: the minimum percentage of individuals that process the locus is 95%, minimum minor allele frequency is 0.1, maximum observed heterozygosity is 0.6, and the maximum proportion of missing data in each site is 40%. To evaluate whether the DNA quality will influence the number of no-missing SNPs in each sample, the vcftools v. 0.1.13 was used to extract the no-missing SNPs of each sample (Danecek et al., 2011). The DNA were grouped into four groups, which have [8000, 10,000), [10,000, 11,000), [11,000, 12,000), [12,000, 13,000) no-missing SNPs in each sample, respectively. The Welch two sample t-test was used to test whether there is significant difference in DNA amount between any paired groups.
To compare the genetic diversity of the two populations of Brown Eared Pheasants, we first calculated the number of heterozygous SNPs in each sample. Here, we define heterozygosity as the number of heterozygous SNPs in each sample divided by the total number of SNPs. We then attempted to estimate the inbreeding levels in the two populations. The 'het' function in vcftools v. 0.1.13 was used to calculate the inbreeding coefficient (Danecek et al., 2011).
To compare the genetic diversity of Brown Eared Pheasant between captive populations and the wild populations. We downloaded re-sequence data of the wild populations from National Genomics Data Center (NGDC; Wang et al., 2020). The re-sequence data is from 19 samples, which include the Shanxi and Shaanxi populations. The accession ID of the data could be found in Appendix Table S2. The pipeline that used to identify SNPs in Wang et al. (2020) was used to analyze the re-sequence data. We also set that the minimum percentage of individuals that process the locus is 95%, minimum minor allele frequency is 0.1, and maximum observed heterozygosity is 0.6 in this SNP set. After we got the clean SNP set, we employed above method to calculate heterozygosity and inbreeding coefficient of each sample.
The independent sample t-test was used to test whether the captive populations and wild populations have similar genetic diversity levels and inbreeding coefficients. To test whether heterozygosity has a negative relationship with the inbreeding coefficient, the RMA method was used to compute the linear regression between the heterozygosity and inbreeding coefficient in the captive populations. The 'lmodel2' R package was used to estimate the linear regression.
We then attempted to evaluate the pedigree of the Brown Eared Pheasants. The KING software v. 2.2.8 was used to estimate the kinship coefficient and proportion of SNPs with zero identical-by-state (IBS0; Manichaikul et al., 2010). Because the input file for the KING software should be in the PLINK binary format, the vcf file was first converted to plink format using vcftools v. 0.1.13 (Danecek et al., 2011). Thereafter, the "make bed" function in Plink v. 1.9 (Purcell et al., 2007) was used to generate the input file for the KING software. The "kinship" function in KING was then used to perform the analysis. Pairs having an estimated kinship coefficient of 0.177–0.354, 0.0884–0.177, 0.0442–0.0884, and < 0.0442 were classified as having first degree, second degree, third degree, and no close relationship, respectively, according to the manual of the KING software.
The mutation load of the two captive populations was evaluated using snpEff v. 5.1d, which was used to annotate and predict the function effect of each variation in the coding region (Cingolani et al., 2012). Because the sequencing depth and coverage of the wild populations is different with that of the captive populations, we did not compare mutation load between captive populations and wild populations. The "build" function in snpEff was used to build the SNP database of the Brown Eared Pheasant. After constructing the database, the vcf file was used as the input file for snpEff. After annotation, snpSift v. 5.1d (Ruden et al., 2012) was used to extract the missense, synonymous, and loss of function variants. To estimate the effect of missense and loss of function variants on the fitness of the species, we calculated the ratio of the number of homozygous missense mutations to the number of homozygous synonymous mutations, and the ratio of the number of homozygous losses of function mutations to the number of homozygous synonymous mutations, respectively, in each sample. The ratios for the heterozygous sites were also calculated. Missense and loss of function variants that have a negative effect on species fitness are usually heterozygous because homozygous variants reduce the survival rate. The independent sample t-test was used to test whether the two captive populations have similar ratios. The workflow of the above analysis is summarized in Fig. 1.
The mean number of no-missing SNPs is 11571 (95% CI = 11380–11762). The DNA amount of the samples that have less than 10,000 SNPs is significantly less than that of the samples that have more than 11,000 SNPs (Fig. 2). However, the more DNA content does not necessarily result in more SNPs in the corresponding samples (Fig. 2; Appendix Table S1). The average mapping rate and coverage was 95.24% (95% CI = 94.81%–95.66%) and 1.92 (95% CI = 1.90–1.93), respectively. The 2b-RAD DNA libraries were successfully constructed for 60 and 39 samples from Taiyuan Zoo, China, and Europe, respectively (Appendix Table S1). The probability of successful DNA-library construction is 90.91% and 51.32% for Taiyuan Zoo and Europe populations, respectively. The current study obtained 43,534 high-quality SNPs for the wild populations. Because the re-sequence data of the captive populations and the wild population is from 2b-RAD and whole genome DNA library, respectively, we did not combine the two SNP sets together.
The genetic diversity of the Taiyuan Zoo population was significantly lower than that of the European population (p < 0.01; Fig. 3A). Both populations have lower genetic diversity than the wild population (Fig. 3A). The mean heterozygosity of the Taiyuan Zoo and European populations was 0.20 (95% CI = 0.19–0.21) and 0.25 (95% CI = 0.22–0.27), respectively. The mean heterozygosity in the wild population is 0.30 (95% CI = 0.25–0.34). The mean inbreeding coefficient of the sampled Brown Eared Pheasants in Taiyuan Zoo was 0.40 (95% CI = 0.37–0.43). Compared to the Taiyuan Zoo population, the European population had a significantly lower inbreeding coefficient (0.17; 95% CI = 0.11–0.23; Fig. 3B). Both captive populations have significant higher inbreeding coefficient than the wild population (0.02; 95% CI = −0.09–0.12; Fig. 3D). The inbreeding coefficients of both captive populations had a significantly negative relationship with their heterozygosity (Fig. 3C and D).
Based on the kinship coefficient, we found that 61.58% of the sampled pairs in Taiyuan Zoo had no close relationship (Fig. 4A, C). In Taiyuan Zoo, the percentage of pairs having first-, second-, and third-degree relationship was 13.90%, 16.38%, and 8.14%, respectively (Fig. 4A, C). The degree of relationship between the individuals in Europe was lower than that between the individuals in Taiyuan Zoo, with 78.00% of the sampled pairs showing no close relationship (Fig. 4B, D). The percentage of pairs having first-, second-, and third-degree relationship was 3.64%, 10.80%, and 7.56%, respectively, in the European population.
The captive population of Taiyuan Zoo had an average of 13.58 (95% CI = 12.14–15.03) and 20.40 (95% CI = 18.87–21.93) homozygous and heterozygous missense variants, respectively. None of the samples from Taiyuan Zoo had the homozygous loss of function mutation, while the average number of heterozygous losses of function variants in this population was 0.60 (95% CI = 0.35–0.85). The average number of homozygous and heterozygous synonymous variants was 15.08 (95% CI = 14.21–15.96) and 37.93 (95% CI = 35.25–40.62), respectively. In the European population, the average number of homozygous and heterozygous missense variants was 44.82 (95% CI = 42.07–47.57) and 38.26 (95% CI = 34.91–41.60), respectively, while the average number of homozygous and heterozygous synonymous variants was 55.74 (95% CI = 51.53–59.96) and 59.49 (95% CI = 53.24–65.74), respectively. The average number of homozygous and heterozygous variants in this population was 1.08 (95% CI = 0.80–1.35) and 0.54 (95% CI = 0.28–0.79), respectively. The ratio of the number of homozygous missense mutations to the number of homozygous synonymous mutations was similar between the Taiyuan Zoo and European populations (Fig. 5A). However, the ratio of the number of heterozygous missense mutations to the number of heterozygous synonymous mutations in the Taiyuan Zoo population was significantly lower than that in the European population (Fig. 5A). The ratio of the number of homozygous losses of function mutations to that of homozygous synonymous mutations in the Taiyuan Zoo population was also significantly lower than that in the European population (Fig. 5B). Meanwhile, the ratio of the number of heterozygous losses of function mutations to the number of heterozygous synonymous mutations was similar between the Taiyuan Zoo and European populations (Fig. 5B).
Understanding the genetic status of threatened species could help us identify the gene pool for conservation efforts. But it is unclear the whole genomic status of captive populations of Brown Eared Pheasant. Here, we combined non-destructive sampling and 2b-RAD technique to evaluate their genetic status. The captive population of Brown Eared Pheasants in Taiyuan Zoo showed lower genetic diversity and higher inbreeding levels than the European population (Fig. 3). Brown Eared Pheasants are raised on several farms in Europe. Because the farmers have opportunities to exchange the birds, inbreeding can be avoided in this population. However, in Taiyuan Zoo, the Brown Eared Pheasants are bred in one feedlot, which increases the opportunities for inbreeding. Frequent inbreeding may be the reason for the low genetic diversity in the Taiyuan Zoo population. It is possible that the sampling process may also contribute to the current results as the possibility of samples being collected from closely related individuals in Taiyuan Zoo cannot be ruled out. Nevertheless, we need to select non-related individuals with high genetic diversity for breeding programs.
The Brown Eared Pheasant population in Taiyuan Zoo harbors recessive lethal mutations. This population harbors heterozygous loss of function mutations (Fig. 5), which suggests that these mutations are recessive lethal mutations. The absence of homozygous loss of function mutations in this population suggests that this mutation type reduces the survival rate. The Europe population harbored homozygous loss of function mutations (Fig. 5), suggesting that this mutation type did not have a deleterious effect on the survival of this population. The European population also harbored more heterozygous missense mutations than the Taiyuan Zoo population (Fig. 5), suggesting that part of the missense mutations in the European population are recessive deleterious mutations. Accurately identifying lethal mutations in threatened species is essential for conservation efforts as it could facilitate the expansion of their populations.
The current study provides insight into conservation management for breeding and future reintroduction process of Brown Eared Pheasants. First, the Taiyuan Zoo should improve genetic diversity of breeding population because the population has low genetic diversity (Fig. 3A). Because mutations in different populations are generally different, introducing Brown Eared Pheasants from other captive populations and wild populations into Taiyuan Zoo could improve genetic diversity of the captive population. Second, it is necessary to avoid inbreeding in captive populations, because the current project found both captive populations have relatively serious inbreeding (Fig. 3B). Inbreeding will lead to accumulation of deleterious mutations in the populations (Xue et al., 2015). Last but not least, evaluating the genetic diversity, inbreeding coefficient and the missense mutations is an important task before selecting individuals for reproduction and reintroduction. Establishing genetic information archive cards of captive individuals may help managers implement the above conservation suggestions.
The workflow of the current study (Fig. 1) could be used in conservation genetics. Since trace DNA (100–200 ng) can be used for constructing 2b-RAD DNA libraries (Wang et al., 2012; Cui et al., 2015). The degraded samples, such as feathers from birds, excrement and hair from mammals, and the skin shed by reptiles, may be used for constructing 2b-RAD DNA libraries (Barbanti et al., 2020). Collecting these kinds of samples will have minimal effects on the survival of threatened species. Researchers could pick up such samples without catching animals, which is benefit for the researchers and animals. Certainly, there will be duplicate samples if we pick up scattered samples. However, the King software could identify the duplicate samples through the relationship inference (Manichaikul et al., 2010). Moreover, the programs used for analyzing 2b-RAD data are commonly used software, such as vcftools. Conservation workers do not need specialized bioinformatics skills to use these programs to analyze the genetic status of endangered species. These advantages will facilitate the widespread use of our workflow in conservation biology. Researchers will obtain population genetic parameters such as genetic diversity following the workflow in the current project (Fig. 1). The whole genomic status of captive populations of threatened species could be evaluated for conservation management using the workflow developed in current study. The current study provides a paradigm for future conservation genetic studies on endangered species based on non-destructive sampling.
The raw sequence data reported in this paper have been deposited in the Genome Sequence Archive in National Genomics Data Center, China National Center for Bioinformation (NGDC, CNCB, https://ngdc.cncb.ac.cn/) under the BioProject PRJCA014106.
PW conceived the experiments, analyzed the data, wrote the manuscript, and reviewed the drafts; PH and JZ analyzed the data and revised the drafts; LZ collected the samples; JZ and ZZ conceived and designed the experiments, and reviewed drafts of the paper. All authors read and approved the final manuscript.
This study was examined and approved by the National Forestry and Grassland Administration of China and the Ethics and Animal Welfare Committee at the College of Life Science, Beijing Normal University, China (Approval Number CLS-EAW-2015-012).
The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Zhengwang Zhang reports financial support was provided by National Natural Science Foundation of China. Pengcheng Wang reports financial support was provided by Priority Academic Program Development of Jiangsu Higher Education Institutions, China. Jing Zhang reports financial support was provided by Beijing Zoo Management Office, China.
We give our best thanks to John Corder, World Pheasant Association (WPA), UK, and Ivan Roels, WPA, Benelux, for their help in collecting the samples in Europe. We also appreciate the WPA members who provided the samples for the project.
Supplementary data to this article can be found online at https://doi.org/10.1016/j.avrs.2023.100078.
Ruden, D., Cingolani, P., Patel, V., Coon, M., Nguyen, T., Land, S., et al., 2012. Using Drosophila melanogaster as a model for genotoxic chemical mutational studies with a new program, SnpSift. Front. Genet. 3, 35.
|
Wang, P., Burley, J.T., Liu, Y., Chang, J., Chen, D., Lu, Q., et al., 2020. Genomic consequences of long-term population decline in Brown Eared Pheasant. Mol. Biol. Evol. 38, 263–273.
|
Zheng, G., 2015. Pheasants in China. Higher Education Press, Beijing.
|