Mitochondrial DNA diversity and maternal origins of Pakistani donkey

Domestic donkey plays a key role as a draft animal in rural economy of Pakistan where its population is increasing every year. The complete mtDNA control region of forty randomly sampled donkeys was PCR- amplified and sequenced bi-directionally using specific primers. Distinct mtDNA haplotypes obtained in the current study (KY446001−KY446011) were subjected to haplotype (h) and nucleotide diversity (π) measures using DnaS as well as to phylogenetic, Network, and AMOVA analyses. There were a total 27 polymorphic sites present within 11 unique mtDNA haplotypes from the studied 40 animals from different regions. Neighbor-joining network and median-joining network both illustrated the splitting of all these haplotypes into two well-defined Nubian and Somali lineages, confirming African maternal origin of Pakistani domestic donkey. Diversity parameters h (0.967± 0.037) and π (0.02917± 0.00307) were found to reveal high levels of genetic diversity in Pakistani donkeys. AMOVA demonstrated only 1% of genetic differences between two mtDNA maternal lineages, pointing to lack of population substructure in Pakistani donkeys as is the case with worldwide domestic donkey population. Pakistani donkeys have African maternal origin and high levels of mtDNA diversity. High genetic diversity may be due to non-selective breeding and heteroplasmy. We herein provide the first report on mtDNA diversity of control region in Pakistani domestic donkey.


Introduction
Donkey has been domesticated about 5000 to 6000 years ago when human stared the use of wild animals for transportation, food and domestic uses (Rossel et al., 2008). The domestication of the Equidae indicates a major cultural shift away from agrarian life styles towards more migration and trade (Kimura et al., 2011). The domestication of donkey has been reported from two different linage origins; Equus africinus africinus and Equus africinus Abstract Domestic donkey plays a key role as a draft animal in rural economy of Pakistan where its population is increasing every year. The complete mtDNA control region of forty randomly sampled donkeys was PCR-amplified and sequenced bi-directionally using specific primers. Distinct mtDNA haplotypes obtained in the current study (KY446001−KY446011) were subjected to haplotype (h) and nucleotide diversity (π) measures using DnaS as well as to phylogenetic, Network, and AMOVA analyses. There were a total 27 polymorphic sites present within 11 unique mtDNA haplotypes from the studied 40 animals from different regions. Neighbor-joining network and median-joining network both illustrated the splitting of all these haplotypes into two well-defined Nubian and Somali lineages, confirming African maternal origin of Pakistani domestic donkey. Diversity parameters h (0.967± 0.037) and π (0.02917± 0.00307) were found to reveal high levels of genetic diversity in Pakistani donkeys. AMOVA demonstrated only 1% of genetic differences between two mtDNA maternal lineages, pointing to lack of population substructure in Pakistani donkeys as is the case with worldwide domestic donkey population. Pakistani donkeys have African maternal origin and high levels of mtDNA diversity. High genetic diversity may be due to non-selective breeding and heteroplasmy. We herein provide the first report on mtDNA diversity of control region in Pakistani domestic donkey.

Mitochondrial DNA diversity and maternal origins of Pakistani donkey
Diversidade de DNA mitocondrial e origens maternas do burro do Paquistão S. Earnist a , S. Nawaz a , I. Ullah a , M. A. Bhinder b , M. Imran a , M. A. Rasheed c , W. Shehzad a and M. Y. Zahoor a *  Sequences were visualized with Chromas Lite 2.1 software and aligned by using MEGA V6 software (Tamura et al., 2013). Eleven haplotypes were mapped from new sequences (Accession Nos. KY446001-KY446011). These haplotypes and 150 previously reported sequences of the domestic and ancient donkeys mtDNA D-loop (Equus kiang, Equus hemionus kulan, Equus hemionus onager, Equus africanus somalicus, Equus africanus africanus) from different region of the world and horse sequence (Equus caballus, AF354440) as an out group was used in phylogenetic and median joining tree construction. We have estimated nucleotide diversities (π) and haplotype (h) within the donkey populations by using DnaSP (Librado and Rozas, 2009) and MEGA V6 to estimate genetic distance between different clades. Neighbor joining trees was constructed based on individual sequences and haplotypes data were constructed by MEGA V6 using Tamura-Nei distances and Kimura 2-parameter distances for D-loop. We developed Median-joining network diagram of all reported sequence along with new sequences by using Network (5.0.0.); to find out the genetic relationships between haplotypes of the D-loop regions (Bandelt et al., 1999). We also used AMOVA for nucleotide variance if they have effect on population (Messina et al., 2018).

Results
We mapped 11 haplotypes resulting from 27 polymorphic sites after analyzing 399 bp of control region readable in all sequences after alignment (Table 1). We found substitution and transversion mutations that show high rate of mutation in donkey population. The transition to transversion ratio of 23:4 demonstrated a strong transition bias. Such bias has also been reported in other domesticated species (Zahoor et al., 2016).
The nucleotide and haplotype diversity of 11 new haplotypes were analyzed as; haplotype (gene) diversity, Hd: 0.967± 0.037 and nucleotide diversity (per site), Pi: 0.02917± 0.00307. Whereas NJ tree with E. asinus somalicus and other E. asinus africanus species reveal different values of genetic diversity. These haplotypes of control region sequence and 27 polymorphic sites demonstrates high polymorphism and genetic diversion in population.
The bootstrap values show that clade1 and clade 2a have common ancestral convergence and whereas clade 2b (E. asinus somalicus) has more distinct divergence. Genetic tree describes E. asinus africanus (Nubian wild ass) and E. asinus somalicus (Somalia wild ass) have common ancestor somehow more close to Nubian wild ass. The other wild ass Asiatic species separate from both these species clade that shows a different sub-species. Inclusive genetic distance and within subpopulations have been described in Table 2. AMOVA value was only 1% for genetic differences between two mtDNA maternal lineages, which shows the Pakistani population have no proper substructure in as is the case with worldwide domestic donkey population.
Genetic distance between clades show that clade 2b have high rate of diversity between clade 1, and clade 2a. It indicates the ancestral species differentiation in local animals. Combinations of all reported breeds sequences somalicus and has been taken place in different sites of Africa (Fesseha, 1991).
The use of donkey as a domestic animal in Indian subcontinent evident in Harappa, ancient city contains the ruins of a Bronze Age fortified city (2600-1900 BC), which was part of the Indus Valley Civilization, presently in Punjab and Sindh provinces of Pakistan (Porter et al., 2016). Donkey is transport source as well as its bypruducts as medicinal and cosmetic uses in Pakistan and meat source in neighbor countries i.e. China. There is no breed description data available of donkey in Pakistan (FAO, 2003). In west border area two breeds are known named as Sperki and Shingari in local language.
Mitochondrial DNA (mtDNA) have advantage of maternal inheritance, extremely low paternal leakage, present in large amount in cell, haploid form of DNA (Cummins et al., 1997). These features make mtDNA a useful and one of the most frequently used markers in molecular systemics (Anderson et al., 1981). It is preferred for genetic diversity, population structure, geography and animal evolution studies (Ankel-Simons and Cummins, 1996;Merwad et al., 2014). There are some related studies in our region for different species (Mustafa et al., 2018;Firyal et al., 2013). We have analyzed genetic diversity of donkey from Pakistani on the base of D-loop region of mtDNA.

Material and Methods
Donkey is used a source of transportation throughout the Pakistan. Animals were selected from the Lahore Multan, some Sindh regions and Jhelum cities donkey markets where animals came from overall the country. Blood samples from 40 animals were collected from the jugular vein of each animal and shifted into 0.5M EDTA coated tubes and was stored in the freezer at 4C 0 . Pictures and area information of donkeys were also taken representing all varieties of country.
DNA extraction was carried out by inorganic method from blood samples (Sambrook and Russell, 2001). Extracted DNA was confirmed on agarose gel electrophoresis and concentration was measured by nano-drop spectrometry. The sequence of mitochondrial complete d-loop of Equus asinus was accessed from gene bank database (www.ncbi. nlm.gov). Primers were designed by Primer3 software from GenBank: X97337.1 for amplification of specific sequence of mitochondrial d-loop region.
Following primers pair was used for amplification of product of 1100 bp ([Forward primer: 5' TAGCTCCACCATCAACACCC 3'] [Reverse: 5' GGCATTTTCAGTGCCTTGCT 3']). The optimized PCR of 25 µL recipe for this reaction was as follow DNA 2µL (30-50ng), dNTPs (10mM) 2.5 µL, Mgcl 2 2 µL, buffer 2 µL, primer forward 1 µL (10 pM), primer reverse 1 µL (10 pM), Polymerase (5U) 3µL and water 14.2 µL. Thermo cycler setting adjusted as; 95 c o for 5 minutes 94 c o , for 30 second, 57 c o for 30 second, 72 c o for 1 minute and final extension at 72 c o was 10 minutes for 25 cycles. The amplicons were sequenced through di-deoxy chain termination method using both orientation of forward and reverse primers and fluorescently labeled products were analyzed on ABI genetic analyzer 3700. Domestication estimation by median joining network tree ( Figure 2) along with other reported sequences indicates multiple domestication events beside of number of sequences used for median joining tree whereas Pakistani donkey don't come in star like haplotype which show that Pakistani donkeys have not been re-domesticated in this region while relates with ancestors domesticated in Africa. If we see overall domestic pattern it shows multiple domestication events in different regions in Africa. This multiple domestication may have increased its genetic diversity and expansion. Whereas E. asinus somalicus adjusted in different clade which is separately domesticated and show highly divergent species.

Discussion
This is the first report on equus asinus mitochondrial DNA in Pakistan as Indus region remained major focus area for trade and defense line from many centuries. There is no breed characterization in the country and we have found 11 haplotypes falling into two clades. Genetic diversity in these samples is; Pi: 0.02917± 0.00307 as compared to other related populations i.e. Ethiopian 0.02070 ± .003, Spanish 0.0006 to 0.0169, Chinese 0.0228 ± 0.0117 ancient Chinese donkeys 0.01815 Turkey Anatolian π=0.00176 ± 0.001 (Lei et al., 2007;Han et al., 2014;Kul et al., 2016).
There are two clades divergence, one specifying Equus asinus africinus, the other Equus asinus somalicus. Bootstrap values 90 of divergent group classify E. asinus somalicus in separate group. While first clad has two sub-clades one relates to E. asinus africinus and second sub-clade have ancestral group E. asinus africinus and E. asinus somalicus. Multiple sub-clades of this group are describing founder affect and expansion from a set of founder haplotypes of mitochondrial DNA. A higher value of Tajima D statistic  donkey which is also seen in our study but heteroplasmy have to addressed along with sequence based genetic analysis (Xu et al., 1996). The most common route of migration was through trade from Africa, Middle & Near East, and then to subcontinent. So, Pakistani donkeys have been originated from E. asinus africinus linage. In conclusion it shows intrusion of donkey to Indus valley from different regions of world through past invasions. The novel haplotypes mapped in Pakistani donkey shows distinctness of local animal to describe indigenous breeds. It demonstrates heteroplasmy of mtDNA that can be helpful in conservation and survival of donkey's population. SNPs mapped can be used as marker for indigenous donkey breeds demarcation. The study will also help in stabilization of donkey population and resolving issue of breed description.
(D: 0.37967) are quite evident for population expansion from founder haplo-groups. Amusingly, all worldwide haplotypes were intermingled among one another, supporting a multiple domestications and re-domestication of all present day donkeys. The haplotypes of this region are present in both sub-clades along with haplotypes from other regions of the world.
Along with high genetic diversity Network Median Joining Tree suggest that Pakistan species fall in African domesticated haplotypes which don't show star like structure (Beja-Pereira et al., 2004;Lei et al., 2007;Kefena et al., 2014;Han et al., 2014). It indicates that donkeys are not domesticated in Indus valley on the other hand there is no re-domestication, no selective pressure and no genetic drift. But they have high genetic diversity which is may be due to non-selective breeding and heteroplasmy in