Genetic Diversity and Structure Analysis of Donkey Population Clusters in Different Indian Agro-climatic Regions

This study conducted with a panel of 24 polymorphic microsatellites revealed high number of alleles and heterozygosity in all the Indian donkey clusters available in different agro-climatic regions. All the markers are highly polymorphic as reflected from high allele number and heterozygosity, thus suitable for genotyping studies. Genetic diversity within each donkey population as well as between populations of different regions was also quite high indicating no extinction threat to population. The genetic differentiation suggested that 89.59% genetic variation corresponds to difference among individuals and 10.41% is attributed to differences among population clusters. However conservation and preservations of donkeys is important as their overall population is decreasing rapidly in country. Even the 24 microsatellites utilized for individual assignment did not reveal 100% assignment of donkeys to their clustered population. Breed relationship analysis revealed closeness among Leh and Bihar donkeys which seems to be unique as geographically these populations are at distant places and mixing of these donkeys cluster is not feasible. Closeness of donkey population from Haryana, Rajasthan, Gujarat and Baramati regions may be due to sharing of common breeding tract and movement of donkeys in these agro-climatic regions with their owners during harsh and odd environmental conditions for their survival. Population structure analysis revealed that donkey population from seven regions can be divided in two main clusters: first cluster having donkey population from Bihar and leh while second cluster included donkeys from Gujarat, Rajasthan, Spiti. Population from Baramati and Haryana had more than 50% individual population as admixed while Spiti donkeys had no admixed individuals. The current study aimed to provide insight into the genetic relationships and diversity between various indigenous donkey populations will offer a valuable reference for rational strategies in breed assignment to local non-descript donkeys, donkey conservation and breeding programs.


Introduction
India possesses a vast genetic resource of equines (1.13 millions) including horses, donkeys, mules and wild asses, distributed over different geographical locations. As per the 19 th livestock census (2012), the total donkey population in the country has decreased by 27.22% over the previous census (2007) and this decrease in population may be due to their decreased demand and lesser utility. Donkeys (0.32 million) are owned mainly by poorest of the poor for earning their livelihood by using them as pack or transport animal in rural and difficult terrains including hot and cold deserts. Inspite of their large and segmented population in different agro-climatic regions in India, these donkeys are known as local, non-descript donkeys without any proper breed characterization.
This most unprivileged animal, donkey "Beast of Burden", contributes maximum for the livelihood of their owners [1]. Recently, one breed of donkey, "Spiti" inhabiting the Himachal Pradesh area of India has been registered based on physical, phenotypic and reproductive traits (INDIA_DONKEY_0600_SPITI_05001).
Beside this, recently variations in phenotypic characters of different donkeys populations have been recorded [2,3].
Since donkey population has suffered a substantial decline in its size, it is expected that levels of inbreeding and risk of extinction may be high which necessitates a genetic study to assess variability among these donkey populations for their conservation along with their breed assignment, better utilization, health, production and management practices etc.
In the present study, we carried out genetic characterization of indigenous donkeys populations available in seven different geographic and agro-climatic regions in India (Gujarat, Haryana, Baramati (Maharashtra), Bihar, Rajasthan, Spiti (HP), Leh (Ladakh,

Journal of Biodiversity & Endangered Species
J&K) to assess genetic diversity within and between populations using molecular markers.
The present study deals with the genetic clustering, extent of genetic diversity and structure between and within seven Indian donkey populations (phenotypic clusters).
Importance of this equine power cannot be overlooked as presently it is difficult to exactly assess and predict its future requirement and availability of fossil fuel.

Animal selection
Out of fifteen agro-climatic regions in India, 286 local and apparently healthy donkeys were selected from seven locations representing six different regions. Efforts were made to give an overall representation to that particular population (cluster) by including a few donkeys from different locations in that particular region. Animals selected for this study were adult, male, non-descript and apparently healthy from each cluster (geographic region).

Blood samples and molecular techniques
Blood samples (5-8 ml) from the jugular vein of all the donkeys were collected aseptically into vaccutainers coated with EDTA (0.5 mM, pH 8.0). DNA was extracted from whole blood samples as described by Sambrook [12]. DNA was isolated and aliquoted in small lots for further utilization after quality check.
Multiplex PCRs were carried out in 5 sets having 4 to 5 microsatellite primers in each set in such a way that label of the same product size differed for their easy recording.
At the end of the reaction, 5.0 μl of stop dye (95% formamide, 0.25% bromophenol blue and 0.25% xylene cyanol) was added and 2 μl of PCR products was loaded on 2% agarose gel, electrophoresed and visualized over ultraviolet (UV) light after ethidium bromide staining to detect the amplification.
PCR multiplex (four to five primers) was combined with 0.3 μl of Liz 500 as an internal lane standard (Applied Biosystems, U.S.A.) and 9.20 μl of hi-di formamide per sample.
The resulting mixture was denatured by incubation for 5 min at 95°C. These denatured samples were run on automated DNA sequencer of Applied Biosystems (ABI 3100 Avant). The electropherograms drawn through Gene Scan were used to extract DNA fragment sizing details using Gene Mapper software, version 3.0 (Applied Biosystems).

Statistical analysis
For all the microsatellites loci, the observed and effective numbers of alleles [13] were calculated using POPGENE software.
The observed and expected heterozygosity estimates were calculated according to Levene and Nei [14,15] as implemented in POPGENE software [14].
A closely related diversity measure, the polymorphism information content (PIC), was estimated [15] which provides a value for utilizing a microsatellite marker for diversity analysis.
The Hardy Weinberg exact test was performed using GENEPOP software [16]. The Fis values were estimated for each of the loci for F statistics [17,18] while overall Fis, Fst, Fit and gene flow (nm) were estimated using Popgene software.
The correspondence analysis was carried out using GENETIX software [19]. The Neighbor-joining tree based on the inter-individual genetic distances in allele sharing was constructed using algorithm of Saitou and Nei [20] as implemented in Phylip software [21].
Genetic distances between different population clusters were obtained using POPGENE software.
Population or breed differentiation was further investigated by using a Bayesian clustering approach implemented in Structure version 2.3.4 [22].
This program generates clusters of individuals based on their multilocus genotypes.
We used an admixture model with a burn-in of 50,000 iterations and 150,000 Markov chain Monte Carlo repetitions to estimate the probable number of genetic clusters (K) without giving any prior information.

Genetic diversity within and between different populations
Twenty four microsatellite markers were successfully amplified with 286 DNA samples of all the seven indigenous donkey population clusters.
Allele count, allele size, shanon index, heterozygosity and Fis values across different loci were recorded (Table 1)  In total, 264 alleles were obtained across all the donkey populations with observed number of alleles (na) per locus ranging from 6 (LEX68; LEX 73) to 22 (AHT05) with a mean value of 11.008 ± 0.160 and effective number of alleles from 1.5618 (COR022) to 9.005 (COR069) with a mean value of 3.8685 ± 0.0795.
For all the loci, value of effective number of alleles was less than that of observed alleles. All the loci (24) had more than 5 number of observed alleles while effective number of alleles were less than 2 alleles at 3 loci, between 2-5 at 6 alleles and rest of the loci had more than 5 effective number of alleles.
All the studied microsatellites loci were polymorphic, which indicates that the microsatellites used were highly suitable for genetic diversity analysis. The observed value of heterozygosity was lowest Except three loci (LEX33, LEX034 and HTG 07), the expected heterozygosity at all other loci was higher than observed heterozygosity indicating that most of the loci showed deviation from Hardy-Weinberg equilibrium.
Inbreeding coefficient values (Fis) ranged from -0.0169 (HTG 07) to 0.7774 (ASB 17) with a mean values as 0.2365 ± 0.011. Three loci (LEX33, LEX 34 and HTG 07) had negative values indicating heterozygosity deficiency while values for nine loci was less than 0.2 indicating low inbreeding among the indigenous populations.

Genetic differentiation and relationship between different donkey populations
The overall means of Fit, Fst and Fis obtained from jack-knifing over loci were significantly different from zero ( Table 2).
An overall significant deficit of heterozygosity (Fis) equal to 13.37% had occurred in surveyed loci within samples. The Fst values ranged from 0.0244 (VHEQ054) to 0.4461 (ASB17) with an overall genetic differentiation of 10.41% among donkey populations.
The overall global deficit of heterozygote across populations (Fit) amounted to 22.39%.
Almost all the loci significantly contributed to the genetic differentiation, deficit of heterozygote within samples and overall global deficit of heterozygote of all the seven indigenous donkey populations.
The gene flow analysis revealed that there was a large effective number of migrants among donkey populations namely between Rajasthan and Spiti, Bihar and Leh donkeys which reveal large effective number of migrants (Nm >1.0) for most of the microsatellite loci analysed in the present study.
The estimates of Fst between each pair of breeds (    Significant deviations from Hardy-Weinberg Equilibrium [17] were also detected in all the populations for most of the loci. The data was analysed for testing for Hardy Weinberg equilibrium using Genepop v 4.0. The two tests conducted were Weir and Cockerham and Robertson and Hill [18]. Using these tests, 9, 10, 7, 10, 7, 6 and 6 loci out of 24 loci for Gujarat, Haryana, Maharashtra, Rajasthan, Spiti, Bihar and leh donkey populations showed deviation from Hardy Weinberg Equilibrium, respectively. Nine loci (HMS06, HTG06, LEX54, LEX33, NVHEQ054, HTG07, ASB02, COR018 and LEX73) in Gujarat donkey population showed deviation from Hardy Weinberg Equilibrium as their P value under heterozygote deficit was less than 0.01. Except locus ASB02, none of the other microsatellite showed deviation in all the donkey populations while rest deviated from 2 to 3 populations only.
The results point to non-random mating of gametes leading to deviation of loci from HW equilibrium and presence of population structure due to inbreeding.

Population relationship and structure analysis
Topology of various donkey populations was prepared on the basis of genetic distances estimated on allele sharing basis (Table 3). On comparing all the donkey populations, donkey from Spiti and Rajasthan were observed to be very close to each other as genetic distance among them was the least (0.0378) whereas donkeys from Leh and Gujarat areas were far apart with a genetic distance of 0.5259. The neighbor-joining algorithm was used for the construction of Phylogenetic tree (Figure 1). Phylogenetic tree indicated that Rajasthan donkey population was very close to Spiti population while donkey population from Bihar was in close proximity with Leh donkey populations. Donkeys populations can be divided in 2 groups: One is having Bihar and Leh donkeys as a separate group whiles other group having donkeys from Gujarat, Maharashtra Rajasthan, Spiti and Haryana. The correspondence analysis (multivariate analysis) was carried out using each individual as a component (Figure 2).  Axis 1 contributes 36.86% while Axis 3 contributes 20.89% to the total variations. It was revealed that donkey populations from Haryana and Gujarat formed a separate group while donkey from Bihar and Leh were quite close to each other. Similarly donkey populations from Spiti and Rajasthan were quite close to each other along with donkeys from Baramati (Maharashtra). The neighbor-joining dendrogram of allele sharing distance formed clearly defined cluster for all the populations (Figure 3). Leh and Bihar populations shared their genetic pools. For further measuring the population structure and extent of admixture, we used structure to find out unbiased structure without using prior knowledge regarding the number of donkey populations. Using ΔK based on the rate of change in the log probability of data between successive K values, the most likely number of clusters (K) present in data was detect. ΔK detected the uppermost hierarchical level (K=2) of the structure (Figure 4). At this (K =2), all the seven donkey populations, clustered into three groups ( Figure 5). First group consisting of donkey populations from Leh and Bihar, second cluster comprised of donkey populations from Gujarat, Rajasthan and Spiti while third cluster comprising donkey populations form Haryana and Maharashtra had admixed population of above two clusters.

Genetic diversity within and between different donkey populations
In any declining animal population or breed, assessment of genetic variability is very important for any breed conservation programme for avoiding the extinction of that breed or animal population [8]. Further for any population genetic study, the precision of the estimated genetic diversity is a function of the number of loci analyzed, the heterozygosity of these loci and the number of animals analyzed from each population [25]. This is the first report on genetic diversity among various donkey populations clustered in different agro-climatic regions in India. Average number of alleles and expected heterozygosity varied among different clusters were quite high, which indicated the presence of high genetic variability among them and suitability of these microsatellites as effective markers to analyze the genetic diversity and phylogenetic relationships between different population clusters. Further among different population, donkeys from Gujarat, Haryana and Baramati regions had high number of alleles than other populations which indicated high genetic variability in these clusters than others. Several studies in exotic donkey breeds namely Spanish, Baudet du Poitu and Catalonian breeds revealed no appreciable difference in genetic variability [5,26,27]. Presence of such a high genetic variability in local donkey clusters might be due to inclusion of highly polymorphic microsatellite loci and random sampling from different areas within the same agro-climatic zones (to capture larger variability) and use of automated DNA sequencer [27]. Lower observed heterozygosity than expected indicated the existence of population structure and local inbreeding effects. Among all the donkey clusters, overall difference between expected and observed heterozygosity was least in donkey cluster from Bihar followed by Leh, Spiti and Gujarat populations which indicated limited inbreeding or random union of gametes. Locus ASB02 revealed Hardy Weinberg deviation in contrast to earlier observations in Spanish donkeys in which it failed to amplify [5].
Out of twenty four microsatellites used, overall inbreeding coefficient values (Fis) in three loci (LEX33, LEX 34 and HTG 07) had negative values indicating out-breeding while values for nine loci was less than 0.2. The population wise analysis revealed that Bihar and Spiti populations had significantly higher inbreeding compared to other populations. Overall estimates of F-statistics was significantly (P < 0.05) different from zero. The average Fis value of all the populations was positive and significant (P < 0.05). The overall population had 13.37% deficit of heterozygotes. The genetic differentiation as revealed by Fst suggested that 89.59% genetic variation corresponds to difference among individuals and 10.41% is attributed to differences among population clusters.

Population structure and structure analysis
Though donkeys available in different agro-climatic zones have never been categorized into breeds, therefore clustered population of donkeys in each region was considered as a separate identity for defining them as a breed. Phenotypic studies recently carried out indicated that donkeys from Leh, Spiti and Bihar areas were of similar size but different from rest of donkeys belonging to Rajasthan, Gujarat, Baramati and Haryana areas [3]. Breed relationship analysis (Phylogenic tree) has clustered donkey populations (Leh and Bihar) separately from other populations (Spiti, Baramati, Gujarat, Rajasthan and Haryana). Closeness among Leh (Cold desert-Zone 1) and Bihar (fertile Middle Gangatic plains-Zone 4) donkeys seems to be unique as these populations are geographically at distant places and mixing of these donkeys cluster is not feasible except for the possible shifting of donkeys in these places long ago. Similarly clustering of Spiti donkey population (cold desert) with donkeys from Rajasthan area (Hot desert) and other areas is not explainable due to similar reasons. However, closeness of donkey population from Haryana, Rajasthan, Gujarat and Baramati (Maharashtra) may be due to sharing of common breeding tract and movement of donkeys in these agroclimatic regions with their owners during harsh and odd environmental conditions for their survival.
The 24 microsatellites were utilized for individual assignment but did not reveal 100% assignment of donkeys to any clustered population. The accuracy of individual specific assignment procedures depends on several factors such as degree of reproductive isolation and genetic differentiation of breeds in question [28] as well as mixing of donkeys during animal sale and purchase fairs, their movement during lean period with their owners in Gujarat, Maharashtra, Haryana and Rajasthan. These seem to be the possible reasons for wrongly assigned individuals. The results revealed that microsatellites may not be used for breed assignment studies in donkeys due to large gene flow, little differentiation and absence of strong population structure.
Population structure analysis enables estimation of a hidden structure that is the number of various clusters (K partitions) without using apriori information about individuals or populations or breeds. We carried out the estimation of K value as suggested by Evanno [29] using ΔK based on the rate of change in the log probability of data between successive K values. ΔK detects the uppermost hierarchical level of structure. The height of the modal value is as an indicator of the strength of signal detected by Structure. At K=2, donkeys from Bihar and Leh formed a single cluster with a few individuals having admixture from rest of populations. Although due care was taken while choosing the representative animals based on their phenotypic characters but at genetic level it is evident that admixed individuals were invariably present in the populations. Maintaining the donkey population in isolation and following systematic breeding plans may avoid mixing of breed or populations can enhance its purity. Second cluster included donkeys from Gujarat, Rajasthan and Spiti. Spiti donkeys had no admixed individuals while the Gujarat and Rajasthan population revealed a few admixed individual. Like Leh and Spiti ponies, donkeys from these two areas also had little or no admixing [30]. This may be attributed to their area specific isolation. Donkeys populations from Haryana and Baramati (Maharashtra) had more than 50% admixed individuals which could be due to mating of these donkeys in their neighborhood areas. Since Population structure is purely based on statistical analysis and does not consider any preinformation regarding breed or individual, it can be accepted to a limited extent as far as donkey population from Bihar and Leh donkeys are concerned. Phenotypically, donkeys from Spiti region are very similar to Leh and Bihar donkeys [31][32][33] but Bayesian clustering approach has clustered this population with Gujarat and Rajasthan population. This genetic variability and structure analysis study will help in identifying the donkey population as a separate breed including their phenotypic and other reproductive traits also.

Conclusion
This study conducted with a panel of 24 polymorphic microsatellites revealed high number of alleles and heterozygosity in all the Indian donkey clusters available in different agro-climatic regions. All the markers are highly polymorphic as reflected from high allele number and heterozygosity, thus suitable for genotyping studies. Genetic diversity within each donkey population as well as between populations of different regions was also quite high indicating no extinction threat to population. The genetic differentiation suggested that 89.59% genetic variation corresponds to difference among individuals and 10.41% is attributed to differences among population clusters. However conservation and preservations of donkeys is important as their overall population is decreasing rapidly in country. Even the 24 microsatellites utilized for individual assignment did not reveal 100% assignment of donkeys to their clustered population. Breed relationship analysis revealed closeness among Leh and Bihar donkeys which seems to be unique as geographically these populations are at distant places and mixing of these donkeys cluster is not feasible. Closeness of donkey population from Haryana, Rajasthan, Gujarat and Baramati regions may be due to sharing of common breeding tract and movement of donkeys in these agro-climatic regions with their owners during harsh and odd environmental conditions for their survival. Population structure analysis revealed that donkey population from seven regions can be divided in two main clusters: first cluster having donkey population from Bihar and leh while second cluster included donkeys from Gujarat, Rajasthan, Spiti. Population from Baramati and Haryana had more than 50% individual population as admixed while Spiti donkeys had no admixed individuals. The current study aimed to provide insight into the genetic relationships and diversity between various indigenous donkey populations will offer a valuable reference for rational strategies in breed assignment to local non-descript donkeys, donkey conservation and breeding programs.