Genetic Diversity and Phylogenetic Evolution of Tibetan Sheep Based on mtDNA D-Loop Sequences

The molecular and population genetic evidence of the phylogenetic status of the Tibetan sheep (Ovis aries) is not well understood, and little is known about this species’ genetic diversity. This knowledge gap is partly due to the difficulty of sample collection. This is the first work to address this question. Here, the genetic diversity and phylogenetic relationship of 636 individual Tibetan sheep from fifteen populations were assessed using 642 complete sequences of the mitochondrial DNA D-loop. Samples were collected from the Qinghai-Tibetan Plateau area in China, and reference data were obtained from the six reference breed sequences available in GenBank. The length of the sequences varied considerably, between 1031 and 1259 bp. The haplotype diversity and nucleotide diversity were 0.992±0.010 and 0.019±0.001, respectively. The average number of nucleotide differences was 19.635. The mean nucleotide composition of the 350 haplotypes was 32.961% A, 29.708% T, 22.892% C, 14.439% G, 62.669% A+T, and 37.331% G+C. Phylogenetic analysis showed that all four previously defined haplogroups (A, B, C, and D) were found in the 636 individuals of the fifteen Tibetan sheep populations but that only the D haplogroup was found in Linzhou sheep. Further, the clustering analysis divided the fifteen Tibetan sheep populations into at least two clusters. The estimation of the demographic parameters from the mismatch analyses showed that haplogroups A, B, and C had at least one demographic expansion in Tibetan sheep. These results contribute to the knowledge of Tibetan sheep populations and will help inform future conservation programs about the Tibetan sheep native to the Qinghai-Tibetan Plateau.


Introduction
Tibetan sheep play agricultural, economic, cultural, and even religious roles in the Qinghai-Tibetan Plateau areas in China and provide meat, wool, and pelts for the local people [1]. The Qinghai-Tibetan Plateau areas are also rich in Tibetan sheep genetic resources, with approximately 17 indigenous sheep populations [2]. Most indigenous Tibetan sheep are not only adapted to their local environment but are also considered important genetic resources and are thus one of the major components of agro-animal husbandry societies. However, most indigenous Tibetan sheep populations are composed of relatively small numbers of individuals, and many populations have been in steady decline over the last 30 years [3]. The climate and landforms of the Qinghai-Tibetan Plateau areas are different from other areas of China. Traffic from other parts of China is blocked; thus, the Tibetan sheep are rarely influenced by external populations. These populations may now be on the verge of extinction and may ultimately be lost, given the rapid destruction of their ecological environment, the continuing introduction of modern commercial Tibetan sheep populations, and the ongoing lack of effective conservation methods [4]. To date, the genetic diversity, phylogenetic relationship, and maternal origin of the Qinghai-Tibetan Plateau populations remain uncertain and controversial.
The study of mitochondrial DNA (mtDNA) polymorphisms has proven to be tremendously useful for elucidating the molecular phylogeny of various species [5][6][7][8] due to the extremely low rate of recombination of mtDNA, its maternal lineage heredity and its relatively faster substitution rate than nuclear DNA [9]. In particular, the control region (CR), also called the displacement-loop region (D-loop) is the main noncoding regulatory region for the transcription and replication of mtDNA. One very useful approach for investigating the history and phylogenic relationships of modern domestic animals is therefore based on mtDNA sequence analysis. The variability and structure of the mtDNA control region makes it possible to describe the genetic polymorphisms and maternal origin of Tibetan sheep, mainly because mtDNA displays a simple maternal inheritance without recombination and with a relatively rapid rate of evolution [10]. The even higher substitution rate in the CR, compared with the heterogeneity rate in the other parts of mtDNA, can be used to optimally characterize intraspecific and interspecific genetic diversity [11][12][13][14][15].
Here, we present an investigation into the mtDNA D-loop variability observed in Tibetan sheep indigenous to the Qinghai-Tibetan Plateau areas. We aimed to increase the number of Tibetan sheep samples by including six available reference genomes from GenBank for our population genetic and phylogenetic analysis of the fifteen Tibetan sheep populations based on the complete mtDNA control region. Our results provide insight into the genetic diversity, phylogenetic evolution, and maternal origin of Tibetan sheep for the conservation and improved management of sheep genetic resources.

Ethic statement
We declare that we have no financial or personal relationships with other people or organizations that can inappropriately influence our work, and there is no professional or other personal interest of any nature or kind in any product, service and/or company that could be construed as influencing the position.

Sample collection
Ten milliliters of blood was collected from the jugular vein of each animal. From the 10 mL samples, 2 mL samples were quickly frozen in liquid nitrogen and stored at -80°C for genomic DNA extraction, as described previously [16]. The total DNA was extracted from the blood using the saturated salt method [17]. The extracted DNA was quantified spectrophotometrically and adjusted to 50 ng/μL. The blood samples were collected from 636 sheep living in the Qinghai-Tibetan Plateau areas in China. The sampled individuals belonged to the fifteen Tibetan sheep populations that are distributed across Qinghai Province (Guide Black Fur sheep, n = 39; Qilian White Tibetan sheep, n = 44; Tianjun White Tibetan sheep, n = 64; Qinghai Oula sheep, n = 44), Gansu Province (Minxian Black Fur sheep, n = 67; Ganjia sheep, n = 58; Qiaoke sheep, n = 71; Gannan Oula sheep, n = 52), and the Tibet Autonomous Region (Langkazi sheep, n = 10; Jiangzi sheep, n = 46; Gangba sheep, n = 85; Huoba sheep, n = 34; Duoma sheep, n = 8; Awang sheep, n = 5; Linzhou sheep, n = 9). The sampling information (population code, sample number, altitude, longitude and latitude, accession number, sampling location, and geographical location) for the fifteen indigenous Tibetan sheep populations is shown in Table 1 and Fig 1. This study did not involve endangered or protected Tibetan sheep populations. All experimental and sampling procedures were approved by the Institutional Animal Care and Use Committee, Lanzhou Institute of Husbandry and Pharmaceutical Sciences, Peoples Republic of China. All samples were collected with the permission of the animal owners.

Data collection
To achieve good coverage of the tested populations, a dataset of six referenced breeds was completed using the six submitted sequences containing the Ovis aries, Ovis vignei, and Ovis ammon mtDNA D-loops for the six individuals in GenBank (Table A in S1 File). These six breeds were from six international geographic regions and included Omusimon, Ovignei, Oammon, OasiaA, OeuropeB, and Omexic. The GenBank accession numbers for these reference sequences are AY091487, AY091490, AJ238300, AF039578 (haplogroup A), AF039577 (haplogroup B), and AY582801, respectively [10,18,19].

Polymerase chain reaction and nucleotide sequencing
One pair of polymerase chain reaction (PCR) primers and sequencing primers was designed based on the 5' and 3' conserved flanking sequences of the complete mtDNA D-loop using the Primer Premier 5.0 software [20] and synthesized by BGI Shenzhen Technology Co., Ltd. (Shenzhen, China). The nucleotide sequences of forward primer CsumF was 5'-GGCTGGGAC CAAACCTAT-3', and the nucleotide sequence of reverse primer CsumR was 5'-GAACAACC AACCTCCCTAAG-3'. PCR was performed in a thermal cycler (Mastercycler gradient, Eppendorf, Germany) with a total reaction volume of approximately 30 μL, containing 2 μL genomic DNA (10 ng/μL), 3 μL (3 pM) each primer, 3 μL 10×Ex Taq reaction buffer, 2 μL dNTP (2.5 mM), 0.2 μL Taq DNA polymerase (5 μL/U) (TaKaRa, China), and 16.8 μL ddH 2 O. The PCR conditions were as follows: initial denaturation for 5 min at 94°C, 36 cycles of denaturation at 94°C for 30 s, annealing at 56°C for 30 s, and extension at 72°C for 1.5 min. The final extension step was followed by a 10 min extension at 72°C. The PCR amplification products were subsequently stored at 12°C until use.
The amplified D-loop fragment was purified using a PCR gel extraction kit from Sangon Biotech Co., Ltd. (Shenzhen, China) and sequenced directly using a BigDye Terminator v3.1 cycle sequencing ready reaction kit (Applied Biosystems, Darmstadt, Germany) in an automatic sequencer (ABI-PRISM 3730 genetic analyzer, Applied Biosystems, CA, USA). PCR for the sequencing was performed in an automatic sequencer with a total reaction volume of approximately 5 μL containing 3 μL genomic DNA (10 ng/μL), 1 μL (3 pM) of each sequencing primer, 0.5 μL BigDye, and 0.5 μL ddH 2 O. The sequencing conditions were as follows: initial denaturation for 2 min at 95°C, 25 cycles of denaturation at 95°C for 10 s, and annealing at 51°C for 10 s. The final extension step was followed by a 190 s extension at 60°C. The PCR sequencing products were subsequently stored at 12°C until use.

Data analysis
The sequences were arranged for multiple comparisons using Clustal Omega [21] and were aligned using ClustalW and BLAST [22]. These results were compared with other sequences obtained from GenBank. The reference sequences for tree construction were taken from the  [24]. To identify differences between the geographic regions using the AMOVA program, four groups were established. The phylogenetic and molecular evolutionary relationships, average number of nucleotide substitutions per site between populations (D xy ), net nucleotide substitutions per site between populations (D a ), ME phylogenetic haplotype and clustering tree, and genetic distance were assessed using Molecular Evolutionary Genetics Analysis (MEGA) version 6.0 [25]. We also sketched network and mismatch distribution graphs using the median-joining method implemented in the NETWORK version 4.6.1.2 software to assess the haplotype relationships [26].

Polymorphic site and sequencing analysis of the complete control region
Based on the reference sequences from GenBank accession numbers (AY091487, AY091490, AJ238300, AF039578, AF039577, AY582801), all of the sequences were aligned with 1274 comparative sites (707 had gaps or missing data, and 567 had no gaps or missing data), and 350  Table 1.
haplotypes were obtained from the 642 sequenced individuals (636 Tibetan sheep and 6 reference sequences). The length of the sequences from the fifteen Tibetan sheep populations of 636 individuals varied considerably, between 1031 and 1259 bp, although the majority were between 1180 and 1183 bp ( Table B in S1 File). A total of 196 variable sites were obtained from the sequences, including 63 singleton variable sites (62 double variants and 1 triple variant) and 133 parsimony-informative variable sites (124 double variants, 7 triple variants, and 2 quadruple variants). There were 158 transitions and 38 transversions within the 196 variable sites, of which 15 sites were found to have both transitions and transversions. The most commonly observed substitution caused a transition mutation. With the exception of the insertion or deletion of several nucleotide sites, the observed variations in the length of the mtDNA D-loop sequences of the Tibetan sheep mainly resulted from variability in the number of 75 bp tandem repeat motifs (between three and five repeats).
The nucleotide composition of all the haplotypes was 32.961% A, 29.708% T, 22.892% C, 14.439% G, 62.669% A+T, and 37.331% G+C. The A+T haplotype was substantially more common than the G+C haplotype, showing an AT bias (Table C in S1 File). The largest haplotype group (haplogroup A) consisted of 490 individuals and 259 haplotypes; the next largest haplotype groups (haplogroup B and haplogroup C) consisted of 145 individuals and 43 haplotypes (64 individuals and 43 haplotypes and 81 individuals and 47 haplotypes, respectively). The smallest haplotype group (haplogroup D) consisted of 1 individual and 1 haplotype. The number of haplotypes, individuals, and frequency detected in each Tibetan sheep population of haplotype group varied from 1 to 49, from 0 to 62, and from 0 to 0.875, respectively ( Table 2). The haplotype diversity and nucleotide diversity were calculated separately for each Tibetan sheep population (Table 2)  Genetic distance and average number of nucleotide differences Table 3

Genetic differentiation
To examine the genetic differentiation between the fifteen Tibetan sheep populations, we calculated Wright's F-statistics of subpopulation within total (F ST ) and genetic differentiation coefficient (G ST ) ( Table 4). We also calculated the gene flow (N m ) ( Table D in

Phylogenetic relationship
To extend our knowledge of the phylogenetic relationship of the fifteen Tibetan sheep populations, a phylogenetic tree was constructed using minimum evolution (ME), neighbor joining using the Maximum Composite Likelihood method (Fig 2) and an unweighted pair-group method with arithmetic means (UPGMA) dendrogram based on the complete mtDNA D-loop sequences of 642 individuals (Fig A, Fig B, and Fig C in S1 File) and 350 haplotypes (Fig 3) from fifteen Tibetan sheep populations and six reference breeds. The six methods produced nearly consistent topological structures and similar support levels; therefore, only the ME tree is presented (Figs 2 and 3). According to the ME tree, NJ tree, UPGMA tree, and median-joining network dendrogram ( (Fig A in S1 File), followed by haplogroups B (Fig B in S1 File) (64) and C (Fig C in S1 File) (81); however, only one animal from the Linzhou sheep (LZ03) belonged to haplogroup D. The Duoma sheep were composed of two haplogroups, the Awang sheep were composed of one haplogroup, and the remaining 13 Tibetan sheep populations were composed of three haplogroups ( Table 2). The four references breeds-OasiaA, OeuropeB, Omusimon, and Omexic-belonged to haplogroups A and B. The other two reference breeds-Omusimon and Ovignei-clustered within a group (Fig 3 and Table F in S1 File. The AMOVA revealed a variation of 4.46% among the populations and of 95.54% within the populations; this finding was significant at P<0.05. The F ST was 0.045, which indicated that 4.5% of the total genetic variation was due to population differences, and the remaining 95.5% came from differences among individuals in each population.

Population expansions
Because the sample size for most of the populations was more than 30 individuals, the detection of population expansion was performed at the individual population level (data not shown) and in all haplotype sequences. The mismatch distribution analysis of the complete dataset (lineages A, B, C, D, and fifteen Tibetan sheep populations of mtDNA D-loop) is shown in Fig 4 and (Table G in S1 File). The charts of the mismatch distribution for the samples of the fifteen Tibetan sheep populations and the total samples were multimodal. However, the mismatch distribution for Linzhou sheep was a unimodal function. The mismatch distribution of the complete dataset showed that there were two major peaks, with maximum values at 4 and 27 pairwise differences and two smaller peaks at 45 and 51 pairwise differences (Fig E in S1 File). These results suggest that at least two expansion events occurred during the population demographic history of the Tibetan sheep population. The mismatch distribution analysis revealed a unimodal bell-shaped distribution of pairwise sequence differences in lineages A, B and C, but that of the lineage D was a sampling function.

High mtDNA D-loop diversity of Tibetan sheep populations
The haplotype diversity and nucleotide diversity of the total individuals were 0.992±0.010 and 0.019±0.001, respectively. The fifteen Tibetan sheep populations in our study showed a high level of haplotype and nucleotide diversity. This finding is consistent with archeological data and other genetic diversity studies [15,[27][28][29], but the haplotype diversity found here was higher than that found in a previous study [30], and the nucleotide diversity found here was lower than that found in a previous study [4]. These results indicate a relatively higher level of genetic diversity in the fifteen Tibetan sheep populations compared with other sheep populations [1,4,31]. For example, the haplotype diversity and nucleotide diversity values of Turkish sheep breeds distributed in a Turkish population were 0.950±0.011 and 0.014±0.001 [31]. However, according to Walsh's work on the required sample size for the diagnosis of conservation units [32], a sample of 59 individuals is necessary to reject the hypothesis that individuals with unstamped ("hidden") character states exist in the population size. Thus, the sample size necessary to reject a hidden state frequency of 0.05 is 56 when sampling from a finite population of 500 individuals. Our genetic diversity estimation is therefore a precise reflection of Tibetan sheep due to the large sample size used in this study. For the Linzhou sheep, Langkazi sheep, Huoba sheep, Qinghai Oula sheep, Guide Black Fur sheep, Tianjun White Tibetan sheep, Ganjia sheep, Qiaoke sheep, Gangba sheep, and Gannan Oula sheep with broad distribution, a high genetic diversity could only be observed with such a large sample size and wide collection area. However, an even higher diversity may be found if even more samples were used, and a further investigation of the genetic diversity of these fifteen Tibetan sheep populations is still worth further research. These Tibetan sheep populations experienced a genetic bottleneck during the 20th century and are classified as the most rare populations of sheep [33]. In addition, the positive Ewens-Watterson and Chakraborty's values were significantly different among the fifteen Tibetan sheep populations, suggesting a previous decline in the population size of the mtDNA D-loop diversity. This finding was consistent with the results of a previous study [33]. Such genetic diversity may be caused by an increased mutation rate in the mtDNA D-loop, the maternal effects of multiple wild ancestors, overlapping generations, the mixing of populations from different geographical locations, natural selection favoring heterozygosis or subdivision accompanied by genetic drift [30].

Maternal origins of the Tibetan sheep populations
The sequence motifs from the 1180 bp to the 1183 bp region of the mtDNA D-loop form the basis for the four major haplogroups (A-D) in the Tibetan sheep mtDNA haplotypes. Of these groups, haplogroup D is quite rare. The Tibetan sheep haplotypes were found to belong to all four major haplogroups, although only 0.157% belonged to haplogroup D, and these sheep were exclusively from the Linzhou sheep haplotype. This finding demonstrated that Tibetan sheep populations possess abundant mtDNA diversity and therefore a widespread origin of their maternal lineages. This study revealed a significant biogeographical association of the Asian Ovis mtDNA haplotypes with haplogroup A. Furthermore, the thoroughbred Tibetan sheep has been proposed to be shared in the haplogroup A, and the contribution of Asian sheep breeds to this population has also been reported. In this study, the overall sequences of all fifteen Tibetan sheep populations, including the fourteen Tibetan sheep populations respectively other than Duoma sheep and Awang sheep, were found in the common haplogroups B and C. It is generally believed that domestic sheep have two maternal lineages (haplogroup A and haplogroup B) based on earlier mtDNA analysis [4,10,18,34]. Recently, a new maternal lineage (haplogroup C) was found in Chinese domestic sheep [27,30]. The ME phylogenetic tree and median-joining analyses in our study revealed the presence of four mtDNA haplogroups in the Tibetan sheep populations. Of these groups, the haplogroup of lineage A was predominant, and the haplogroups of lineage B and lineage C were the second most common. In this paper, the proportion of haplotypes of lineage D was 0.157%, further demonstrating that lineage D is the most rare of the mtDNA lineages. Our findings were consistent with the results of previous studies on domestic sheep breeds in China [35,36]. Previous studies identified three mtDNA haplogroups in both China [29,36,37] and other countries [38,39,40]. The four mtDNA haplogroups of lineages A, B, C, and D found in the Tibetan sheep populations in the Qinghai-Tibetan plateau areas further supported the hypothesis of multiple maternal origins in Chinese domestic sheep.

Genetic differentiation of Tibetan sheep populations
The F ST value represents the level of genetic differentiation within a given population. Thus, there is "little differentiation" at a value of 0.05, "moderate differentiation" at values of 0.05-0.25, and "great differentiation" at values >0. 25. In this study, the AMOVA analysis also revealed the distinct population of Qinghai-Tibetan Plateau areas among other Tibetan sheep populations with a significant positive variance. Gene flow (N m ), also known as gene migration, refers to the transfer of alleles from one population to another. N m haplotype values >1 and N m sequences <1 indicate a poor gene exchange, such that genetic drift will result in substantial local differentiation [41,42]. The low G ST value, combined with the low N m of sequences used in this study, indicate that the great differentiation mainly resulted from the independent evolution of each isolated population and substantial local differentiation caused by the genetic drift [43]. An important factor leading to this result is likely the lower effective population sizes, as the Gannan Oula sheep, Qiaoke sheep, Ganjia sheep, and Qianlian White Tibetan sheep live in canyons and valleys and therefore have a limited ability to migrate and correspondingly lower population sizes relative to the other Tibetan sheep populations. As the effective population size declines, the nucleotide substitutions have a greater probability of reaching fixation [44,45]. In addition, the estimated divergence time (data not shown) among the fifteen Tibetan sheep populations was consistent with the Pleistocene climate fluctuations and the uplift of the Qinghai-Tibetan Plateau, indicating that known paleogeographic factors might have played important roles in the speciation of Tibetan sheep.

Genetic relationships among the Tibetan sheep populations
Our study showed that the fifteen Tibetan sheep populations native to the Qinghai-Tibetan Plateau are clustered into four groups: 490 Tibetan sheep represent the maternal origin of the haplogroup of lineage A, 64 Tibetan sheep represent the maternal origin of the haplogroup of lineage B, 81 Tibetan sheep represent the maternal origin of the haplogroup of lineage C, and 1 Tibetan sheep represents the maternal origin of the haplogroup of lineage D. This genetic relationship displayed a high consistency with traditional classification schemes and the results of previous studies [27,[46][47][48][49]. All fifteen Tibetan sheep populations belong to four maternal origins. The genetic differentiation of the fifteen Tibetan sheep populations was mainly the result of geographic isolation, natural selection, different living conditions, and breeding history. Because Tibetan sheep are a portable food and wool resource, the commercial trade and extensive transport of sheep along human migratory paths might help account for the observed pattern by promoting genetic exchange. Other study methods, such as genetic approaches, including the degree method and the phylogenetic relationship clustering method, also indicated that indigenous sheep were the maternal origin of haplogroups A, B, C, and D [46,48].

Population expansion of Tibetan sheep populations
Because the sample sizes of most of the populations were less than 34 individuals, the detection of population expansion was performed at the level of the individual populations (data not shown). The mismatch distribution analysis of the complete dataset, haplogroups A, B, C, D, and fifteen Tibetan sheep populations of the mtDNA D-loop, is presented in Fig 4 and (Table G in S1 File). The complete dataset of all Tibetan sheep populations had a significantly large negative Tajima's D value and F S value (Tajima's D = -0.466, p = 0.020; F S = -7.484, p = 0.001). This result was consistent with a demographic model showing two large and sudden expansions, as inferred from the mismatch distribution. The mismatch distribution of the complete dataset suggested that there were two major peaks with maximum values at 4 and 27 pairwise differences and two smaller peaks at 45 and 51 differences. These results suggest that at least two expansion events occurred in the population demographic history of the Tibetan sheep living on the Qinghai-Tibetan Plateau. The mismatch distribution analysis revealed a unimodal bell-shaped distribution of the pairwise sequence differences in haplogroups A, B, and C. However, the distribution of lineage D was a sambong function. Mismatch analysis of haplogroups A, B, and C suggested that single population expansion events occurred in the demographic history of the Tibetan sheep populations. This finding was similar to the previously reported results [29].

Conclusion
China holds abundant populations of Tibetan sheep, with significant mtDNA haplotype diversity observed in the sheep of the Qinghai-Tibetan Plateau areas. Here, the large-scale mtDNA D-loop sequences analysis of fifteen Tibetan sheep populations has provided evidence for four maternal haplogroups with high diversity. Phylogenetic analysis showed that all four previously defined haplogroups (A, B, C, and D) could be identified in the 636 tested individuals of the fifteen Tibetan sheep populations, although the D haplogroup was only found in the Linzhou sheep. The estimation of demographic parameters from the mismatch analyses shows that haplogroups A, B, and C had at least one demographic expansion in the Tibetan sheep of the Qinghai-Tibetan Plateau areas.
Supporting Information S1 File. The UPGMA phylogenetic tree show that the 490 sequences of 15 Tibetan sheep populations ( Fig A). The UPGMA phylogenetic tree show that the 64 sequences of 12 Tibetan sheep populations (Fig B). The UPGMA phylogenetic tree show that the 81 sequences of 14 Tibetan sheep populations (Fig C). Median-joining networks for the mtDNA D-loop in the control region show that 636 sequences of 15 Tibetan sheep populations and six reference breeds fall into five distinct cluster haplogroup A, haplogroup B, haplogroup C, haplogroup D and haplogroup E, respectively. The majorities of the 490 individuals were grouped in haplogroup A, followed by haplogroup B (64) and C (81); however, only one animal from LZ03 belonged to haplogroup D. The AW population was composed of one haplogroup, the DM population was composed of two haplogroups, and the remaining 13 Tibetan sheep populations were composed of three haplogroups (Fig D). The mismatch distribution of the complete dataset of the mtDNA types of Tibetan sheep of the four lineages on the Qinghai-Tibetan Plateau areas showed that there were two major peaks, with maximum values at 4 and 27 pairwise differences and two smaller peaks at 45 and 51 pairwise differences (Fig E). Mitochondrial genomes of the 6 reference breeds included in the phylogenetic analyses in this study ( Table A). The length of the complete mtDNA D-loop sequence in fifteen Tibetan sheep populations (Table B). Base pair composition of mtDNA D-loop of fifteen Tibetan sheep populations (Table C). Gene flow (Nm) of the sequence (above the diagonals) and Nm of the haplotype (below the diagonals) between fifteen Tibetan sheep populations (Table D). Dxy (the average number of nuc. subs. per site between populations) (above the diagonals) and Da (the number of net nuc. subs. per site between populations) (below the diagonals) of the difference in the number of nucleotides per site differences between fifteen Tibetan sheep populations (Table E). Hierarchical analysis of the molecular variance (AMOVA) of the D-loop region of mtDNA for fifteen Tibetan sheep populations (Table F). Neutrality tests for fifteen Tibetan sheep populations (Table G). (XLSX)