Insights into the genetic diversity of indigenous goats and their conservation priorities

Objective An experiment was conducted to evaluate genetic diversity of 26 Chinese indigenous goats by 30 microsatellite markers, and then to define conservation priorities to set up the protection programs according to the weight given to within- and between-breed genetic diversity. Methods Twenty-six representative populations of Chinese indigenous goats, 1,351 total, were sampled from different geographic regions of China. Within-breed genetic diversity and marker polymorphism were estimated calculating the mean number of alleles, observed heterozygosities, expected heterozygosities, fixation index, effective number of alleles and allelic richness. Conservation priorities were analyzed by statistical methods. Results A relatively high level of genetic diversity was found in twenty-four population; the exceptions were in the Daiyun and Fuqing goat populations. Within-breed kinship coefficient matrices identified seven highly inbred breeds which should be of concern. Of these, six breeds receive a negative contribution to heterozygosity when the method was based on proportional contribution to heterozygosity. Based on Weitzman or Piyasatian and Kinghorn methods, the breeds distant from others i.e. Inner Mongolia Cashmere goat, Chengdu Brown goat and Leizhou goat obtain a high ranking. Evidence from Caballero and Toro and Fabuel et al method prioritized Jining Gray goat, Liaoning Cashmere goat, and Inner Mongolia Cashmere goat, which agree with results from Kinship-based methods. Conclusion Conservation priorities were determined according to multiple methods. Our results suggest Inner Mongolia Cashmere goat (most methods), Jining Gray goat and Liaoning Cashmere goat (high contribution to heterozygosity and total diversity) should be prioritized based on most methods. Furthermore, Daiyun goat and Shannan White goat also should be prioritized based on consideration of effective population size. However, if one breed can continually survive under changing conditions, the straightforward approach would be to increase its utilization and attraction for production via mining breed germplasm characteristics.


INTRODUCTION
Goats, one of the most ancient livestock species, were domesticated ~10,000 years ago in Southwestern Asia [1].It is widely accepted that Chinese indigenous goat breeds originated from the plateau of Southwest China and nearby regions in central Asia.About 39% of docu mented mammalian livestock breeds are goats, but 7% of known goat breeds are extinct and 20% are threatened.Longterm natural and artificial selection, imposed by environmental changes and animal husbandry, have resulted in 58 indigenous and 8 cultivated goat breeds.A considerable number of desirable traits occur in Chinese indigenous and domesticated goats, such as extensive adaptability to stressful environments, outstanding disease and cold resistance, strong coarse fodder resistance, and high prolificacy.These traits are most likely multigenic and therefore represent a diverse natural gene pool.Interesting examples include the Jining gray goat (JNQ, associated with genes which determine the color of lambskin) and the Zhongwei goat (ZWS), the only furbearing goat in the world (associated with genes respon sible for the unique fur phenotype).This genetic resource is invaluable for future goat breeding and utilization.
Recently, Chinese indigenous goat breeds have been threat ened by the introduction of exotic goat breeds, which have typically been selected for optimal production of meat, wool, or other products.Because shortterm profit maximization in centivizes the replacement or crossbreeding of Chinese indigenous goat breeds with more productive breeds, the indigenous gene pool is at risk.Both population size reduction and genetic variation dilution are major threats.There has been an in creasing tendency to neglect the unique qualities and variety offered by indigenous breeds in favor of exotic breeds and their crosses.However, it is now recognized that it is im portant to establish conservation priorities and strategies to conserve genetic diversity within and betweenbreeds, primarily to avoid further losses of genetic resources.Gene tic diversity is one the most important factors that determines whether a breed survives and flourishes, or ultimately faces depression and extinction.The genetic diversity of domestic animals is crucial to meet current production needs in various environments, prospective and changing breeding objectives, and sustained genetic improvement.In response, various conservation strategies have been implemented.However, in order to make the best use of limited conservation funds, it is necessary to prioritize specific indigenous breeds for conservation [2].
Several conservation approaches have focused on breeds with maximum conservation value (i.e., breeds with high levels of genetic variation) [3,4].The Weitzman approach, which has been widely employed to establish conservation priorities, uses betweenbreed genetic diversity to prioritize breeds that are highly distant from others based on genetic distance [4,5].However, the Weitzman method has been criticized because it does not take withinbreed diversity into consideration [5,6].Other factors, including rare alleles that occur at anomalously high frequency due to inbreeding, strict genetic isolation, or founder effects, can compromise the effectiveness of the Weitzman method [7].Another approach relies primarily on withinbreed genetic diversity, which is usually calculated in terms of average expected heterozygosity (H E ).The strategy prioritizes a breed if its removal from a population results in a maximum loss of global average H E .However, the draw back to this method is that it is insensitive to betweenbreed genetic diversity.Alternative methods that consider between and withinbreed diversity, based on coancestry or kinship, have been developed to prioritize breeds for conservation [8,9].Although this approach ideally includes molecular and pedigree data, pedigree information is not available for many breeds, and thus in practice the method mainly depends on molecular data.To combine between and withinbreed di versity, various weights can be assigned to balance the two components [5,1012] F ST (Wright's fixation index in total population) can be used to adjust the betweenbreed diversity derived from the Weitzman method, and 1F ST to adjust the withinbreed diversity derived from the H E change in a meta population [5].For example, Piyasatian and Kinghorn assigned betweenbreed components five times more weight than withinbreed components [11].Ginja et al [13] and Cañón et al [14] conducted conservation priority analyses for cattle using different methods.No one method has yet emerged as the best, but together they provide a comprehensive view of genetic diversity that can inform a conservation program.
Conservation priorities have been studied using the methods described above for domestic animals, including cattle [14], pig [10], sheep [15], chicken [16], horse [17], and goat.How ever, priorities for Chinese indigenous goat breeds have not been examined, except for one study in which the Weitzman approach was used to rank 12 mutton goat breeds [18].In the work presented here, we conducted a comprehensive analy sis of conservation priorities for 26 Chinese indigenous goat breeds, based on 30 microsatellite markers.

Animal care
Protocols for blood and tissue sampling were approved by the Biological Studies Animal Care and Use Committee of the National Animal Husbandry Services, Beijing, People's Republic of China.All experimental procedures followed guidelines established under the Law of Animal Husbandry in the People's Republic of China (No. IASCAASAE03,12 Dec 2016).The experimental procedure was approved by the Institutional Animal Care and Committee at National Animal Husbandry Service.

Sampling
In total, 1,351 samples from 26 representative Chinese indige nous goat breeds (Supplementary Figure S1) were used in this study.Goats from 19 provinces or autonomous regions in China were sampled, covering high altitude regions and plains.The blood samples were collected at major production centers meeting the characteristics and features of the breed, from healthy animals that were unrelated within three generations.
Thirty microsatellite markers recommended by the In ternational Society for Animal Genetics (ISAG)/Food and Agriculture Organization of the United Nations (FAO) were used to genotype the 1,351 samples.The markers are CSRD247, ETH10, MAF209, OARAE54, SRCRSP15, SRCRSP3, SRCRSP5,  MAF70.DNA was isolated from cryopreserved blood samples using the DNeasy Blood and Tissue Kit (Qia gen, Valencia, CA, USA).Polymerase chain reaction (PCR) amplifications were conducted using fluorescencelabelled primers.PCR products were analyzed by capillary gel elec trophoresis using an ABIPRISM 3130xl Genetic Analyzer (Applied Biosystems, Foster City, CA, USA) and LIZ 500 internal size standards (Applied Biosystems, USA).Data were collected and analyzed using GENEMAPPER 3.5 (Applied Biosystems, USA).

Genetic diversity analysis
Withinbreed genetic diversity and marker polymorphism were estimated using GENALEX [19] and the Excel Micro satellites Toolkit to calculate the mean number of alleles (MNA), observed heterozygosities (H O ) and H E .Allelic richness (R t ) over all loci for each breed was computed with FSTAT v 2.9.3.2 [20].In this study, estimation of N E(0.05) were performed with linkage disequilibrium method and random mating model (minor allele frequency 0.05) using LDNe program [21].
The methods of conservation analyses used in this study have been demonstrated previously [13,14].To better under stand the statistical results, three approaches were employed: i) a method designed to minimize the overall kinship coeffi cient (Core Set method); ii) a method that considered between breed diversity alone (Weitzman approach); iii) a method aimed at combining the within and betweenbreed compo nents of overall genetic diversity.These are described in detail below.

Core Set methods
The Core Set method [22] is based on measures of coancestry or kinship, and aims to eliminate the genetic overlap between populations in a core set.The coefficient of kinship was used as a measure of genetic similarity [22].Negative contribution estimates were avoided through an iterative process that assigns the lowest value as zero, and recalculates contributions after a removal of the population.A coefficient matrix of kinships was estimated using three methods to correct for alleles iden tical: i) markerestimated kinships (MEK), in which variation is based on weighted equal drift similarity (WEDS) [23]; ii) MEK estimation based on WEDS, using a bootstrap proce dure (bootstrap 250); iii) a variation of the MEK method based on log linear regression obtained with a weighted loglinear model (WLM) [24].Analyses of conservation priorities based on similarity matrices were performed using the FORTRNA application, which was developed and kindly shared by Eding and Meuwissen [22,24].
Pairwise kinship distances were computed from the MEK matrix following Eding et al [22].Neighbornet phylogenies for all breeds were constructed with kinship genetic distance using SPLITS TREE soft (version 4.12.6)[25].Average mole cular coancestries (fm), based on allele frequency, were analyzed using MOLKIN3 [26].Genetic relationships were used to classify the goat breeds and generate contour plots of kinship coefficients (MEK based on WEDS and f m ) using the colorRamps package in R.
Withinbreed contributions to diversity were computed using the average withinbreed H E .The partial contribution made by individual breeds was computed as the proportional variation in average internal heterozygosity of the complete data set after removal of each breed (PC He ).

Weitzman approach
The Weitzman method was also used to calculate the partial contribution made by each breed to total genetic diversity (PC Weitz ), based on Reynolds genetic distance (on diversity, Weitzman).This method estimates betweenbreed diversity and ignores withinbreed diversity.Analyses conducted using the Weitzman method were computed using a FORTRAN program developed by Ollivier and Foulley [5].Reynolds ge netic distances were computed using MOLKIN3 software [26].

Combined methods
When making decisions about conservation priority, both within and betweenbreed genetic diversity should be taken into account [5].Three approaches were used to calculate the two components of global diversity for a metapopulation: i) Ollivier and Foulley [5] proposed that the aggregate diversity (PC Fst ) should weight the betweenpopulation diversity using the factor F ST , and weight the withinpopulation diversity using the factor 1F ST ; ii) Piyasatian and Kinghorn proposed that the betweenpopulation diversity component should receive five times more weight than the withinpopulation diversity (PC 5:1 ) [11]; iii) Caballero and Toro [9] and Fabuel et al [10] assigned equal weights to the withinpopulation coancestries and genetic distances.In this case, calculations were performed using MOLIKN3 [26].Pairwise correlation analysis was per formed using the PICANTE package of R with the Pearson method.The level of allelic diversity in this study was similar to that reported for 9 Chinese cashmere goats [27], but lower than in the Chinese Cashmere [28] and Swiss [29] goat populations.The heterozygosity of most of the breeds examined in this study was also lower than the values reported in these studies, particular the H O .Loss of genetic diversity among populations can be caused by genetic introgression [30], environmental change [31] and inbreeding.In our study, F IS had a positive value in all breeds, demonstrating a deficit of heterozygosity (Table 1).F IS is usually used to obtain a deeper understand ing of the degree of endangerment and inbreeding, and is also considered as a criterion for conservation priority.In present study, 5 breeds had F IS values less than 0.05; 16 breeds fell into the range 0.05 to 0.15, and 5 breeds in the range 0.15 to 0.25.These levels are higher than previously reported in Chinese cashmere goats [27].Effective population size (N E ) is also a valuable indictor for evaluating conservation priorities.

RESULTS AND DISCUSSION
Effective population size (N E0.05 ) were computed based on linkage disequilibrium with minor allele frequency 0.05 (Table 1).N E0.05 ranged from 32.7 DYS up to 569.7 ZWS.According to Franklin's 50/500 rule of thumb [32], effective population size under 50 suggested that population is facing to serious genetic threaten.Small population are easy to inbreeding depression.The lowest population size of DYS means that it is necessary to implement conservation program.According to allele statistics, lower allele richness, H O , H E , N EA , MNA, and F IS of Daiyun indicated a relatively high inbreeding.The reasons behind is human intervention of mating process which influencing by introduction of breeds with high perfor mance and intensive breeding.So, DYS and Shannan White goat (SNB) (N E0.05 = 47.5) should be paid more attention and considered in conservation program.
The neighbornet phylogeny of kinship distances shows the relationships among the 26 Chinese indigenous popula tions (Figure 1a).The phylogeny of kinship distances was similar to that generated using the Reynolds distance (Sup plementary Figure S2a).Goats from northern China, such as the Inner Mongolia Cashmere goat (MGR), Tibetan goat, ZWS, XJS, Chaidamu goat (CDS), and SNB clustered together, and were close to the Huanghuai goat (HWS) and JNQ.Ex cept for the Shannan white goat, these breeds belong to the North China cluster [32].The DYS, Fuqing goat (FQS), and Xiangdong Black goat (XDH) clustered together with a long branch length, like Longlin goat (LLY) and Chengdu Brown goat (CDM) on one side, with Longling yellow goat (LLS), Zhaotong goat (ZTS), Yuling goat (YLS), and Maguan poll goat (MGS) on the other side.The Matou goat (MTS) and Yichang White goat (YCB) clustered together, partly consis tent with a previous study by Wang et al [33] that placed MTS, YCB, and SNB into the Central China cluster.Although Wang et al [33] also proposed that LLS, MGS, YLS, ZTS, CDM, and LLY belonged to the Southwest China cluster, our results clus ter LLY and CDM within a separate branch.The discrepancy may reflect the fact that more loci were used in the present study, supporting a more accurate analysis.The majority of goat populations cluster in agreement with their geographi cal origin.
Withinbreed diversity could also be estimated using kin ship coefficients with either the MEK derived from genotypes or the average coancestries (f m ) obtained from allele frequencies.
To aid in the analysis, contour plots were plotted to visualize both within and betweenbreed kinships.Each population was sorted based on its genetic proximity defined in the phy logenetic neighbornet graph (Figure 1b; Supplementary Figure S2b).light blue and light yellow areas also reveal goat breeds with relationships.The two groupings consist of the LLS, MGS, YLS, and ZTS (0.136<MEK<0.179, 0.348<f m <0.380), and the DYS, FQS, XDH, YCB, and MTS (0.092<MEK<0.160, 0.328< f m <0.407).

Analyses of conservation priorities for Chinese indigenous goat
Based on different approaches, analyses of conservation pri orities for Chinese indigenous goat populations is presented in Tables 2 and 3. Fourteen (out of 26) populations were es timated to make a null contribution to total genetic diversity, based on the Bootstrap, WEDS, and WLM kinship methods (Table 2).The three highly prioritized breeds were the JNQ, Liaoning Cashmere goat (LNS), and MGR (0.1361<WEDS< 0.1622, 0.1324<Bootstrap<0.1613and 0.1071<WLM<0.2621).
A similar ranking of conservation priority was obtained using all three methods.This may indicate that high withinbreed genetic diversity exists in these populations.
A considerable number of populations with negative value for PC He (14 breeds) were observed, based on the proportional contribution of each breed to the average H E for the whole population (Table 2).For these populations, there is a gain in total diversity after they are excluded.As expected, the major ity of these breeds have highly inbred status (high withinbreed kinship coefficients) and relatively low heterozygosity (low H E ) (Figure 1b; Table 1).The XDH, DYS, FQS, and MTS have the most negative PC He values (-0.218<PCHe <-0.195) (Table 2).The most important conservation priorities, ranked high to low, are the JNQ, MGR and LNS (0.258<PCHe <0.363) based on PC He value (Table 2).Breeds with intermediate PC He values include the Tibet goat, HWS, Lubei White goat (LBB) and CDS (0.173<PCHe <0.227) (Table 2).
The conservation priorities assigned using the Weitzman  2).Ollivier and Foulley [5] proposed that both within and betweenbreed components should be taken into consideration (PC Fst ) [5].In this approach, the F ST value of the whole popu lations is used to weight PC Weitz (F ST = 0.138 in this study) and 1F ST is used to weight PC He .The PC Fst method prioritizes breeds by kinship (i.e., the HWS with PC Fst = 0.969 and WEDS = 0.2017), and also prioritizes breeds that contribute highly to PC Weitz (i.e., the CDM with PC Fst = 1.034 and PC Weitz = 6.71).The breed with the highest priority is the MGR (PC Fst = 1.235,WEDS = 0.1622, WLM = 0.1071, PC He = 0.307 and PC Weitz = 7).Although the Tibet goat makes an intermediate contribu tion to total diversity as estimated by various methods (WEDS, WLM, PC He , and PC Weitz ), it also has a relatively high PC Fst value (PC Fst = 0.951).However, the ranking derived from the PC 5:1 method (as proposed by Piyasatian and Kinghorn) is similar to the priority estimated using PC Weitz , because it pri oritizes the CDM, MGR, and LZS, due to the increased weight given to the betweenbreed component.
The results obtained from the combined approach proposed by Caballero and Toro [9] and Fabuel et al [10] are shown in Table 3.The breeds making the highest contributions to global coancestry were the FQS, MTS, and LZS this is the result of high withinbreed coancestry (0.423<f ii <0.499) and the rela tively low distance from all the other populations (0.095<D nei < 0.120).Although the DYS has a relatively high f ii value (0.499), its mean genetic distance was also larger (D nei = 0.147) and thus its contribution to f was relatively low (0.0119).The dif ference between f ii and D nei is responsible for calculating the contribution to global coancestry.The absolute contribution to total diversity prioritized the JNQ (0.0356), CD (0.0343) and Xijiang goat (0.0335).Ranking by absolute contribution yield ed results similar to assigning priority based on proportional contributions to genetic diversity.Using absolute contribu tions, priority was assigned to the JNQ (4.808), CDS (4.632), and XJS (4.524).The ZWS has the lowest contribution (1.931), perhaps due to its relatively low sample size (n = 25).When the proportional contribution to genetic diversity was estimated without reference to sample size, only the JNQ (4.240) and CDS (4.078) maintain their higher priorities.The Tibet goat also obtained a high priority (4.051).In contrast, the DYS made the smallest contribution (PC1 = 3.052 and PC2 = 3.362).Finally, total genetic diversity was analyzed after removing one subpopulation at a time from the total population.Re moval of the JNQ resulted in the largest reduction in total genetic diversity (GDT|i = 0.7378, loss/gain = -0.4%),fol lowed by the MGR (GDT|i = 0.7382, loss/gain = -0.3%)and LNS (GDT|i = 0.7386, loss/gain = -0.3%).
Linear correlation coefficients calculated for pairwise con tributions derived from different methods are shown in Table 4.No negative correlation occurs between any pairwise con tributions, even between Weitzman and H E contributions (0.321).The correlations between PC 5:1 and other methods are similar to the correlation between PC Weitz and others be cause of the excess weight given to betweenbreed components in PC 5:1 .Compared to PC weighted , higher correlations were ob served between PC unweighted and other methods.As there are too many null contributions in the three core set methods, it is not useful to calculate the correlation between these and the other methods, despite the higher correlation between core set methods and PC He .

Marker polymorphisms, within-breed diversity and breed relationships
Summary statistics describing microsatellite marker polymor phisms and genetic diversity per breed are presented in Table 1 and Supplementary Table S1, respectively.A total of 430 alleles were observed in 26 Chinese indigenous goat popula tions, with locus BM6444 showing the largest observed (N A , 12.385±0.60)andeffective(NE,6.202±0.299)allelenumber(SupplementaryTableS1).In contrast, the smallest values for N A and N E were 2.308±0.092and1.393±0.06atlocusMAF209.The mean H O for all loci is 0.577±0.008,rangingfrom0.269±0.041 at locus MAF209 to 0.785±0.031atlocusOARFCB48.The mean H E is 0.638±0.006,rangingfrom0.25±0.03atlocusMAF209 to 0.827±0.01atlocusBM6444.The F ST value across the 30 loci has a relatively high mean (0.141±0.010), indicat ing that there is genetic differentiation among the 26 Chinese indigenous goat populations.There are 23 loci exhibit posi tive F IS values, and the positive mean F IS is 0.082±0.029.This may indicate nonrandom mating, and these loci may be mor phological or productive traits under selection.The MNA was relatively high in the 26 Chinese indigenous goat populations, of which the largest was in the Xinjiang goat (XJS, 8.23±2.93),and the lowest in the Daiyun goat (DYS, 4.40±2.25)(Table1).However, the highest and lowest values for the effective number of alleles were 4.369±0.299 in the JNQ and 2.424±0.217 in the DYS.The mean Rt was 5.61± 2.04, ranging from 7.03±2.27 in the XJS to 3.89±1.91 in the DYS.Additionally, higher estimates of expected and H O were observed in the JNQ (H E , 0.737±0.024;H O , 0.646±0.012),whereas the values of H E and H O were 0.507±0.039and 0.434± 0.013 in the DYS, representing the smallest heterozygosity (Table 1).All populations in which H E exceeded H O yield positive values for F IS .