Sequencing Reveals Population Structure and Selection Signatures for Reproductive Traits in Yunnan Semi-Fine Wool Sheep (Ovis aries)

Yunnan semi-fine wool sheep are among the most important cultivated sheep breeds in China. However, their population structure, genetic characteristics and traits of interest are poorly studied. In this study, we systematically studied the population characteristics and selection signatures of 40 Yunnan semi-fine wool sheep using SNPs obtained from whole-genome resequencing data. A total of 1393 Gb of clean data were acquired. The mapping rate against the reference genome was 91.23% on average (86.01%–92.26%), and the average sequence depth was 9.51X. After filtering, 28,593,198 SNPs and 4,725,259 indels with high quality were obtained. The heterozygosity rate, inbreeding coefficient and effective population size of the sheep were calculated to preliminarily explore their genetic characteristics. The average heterozygosity rate was 0.264, the average inbreeding coefficient was 0.0099, and the effective population size estimated from the heterozygote excess (HE) was 242.9. Based on the Tajima’s D and integrated haplotype score (iHS) approaches, 562 windows and 11,356 core SNPs showed selection signatures in the Yunnan semi-fine wool sheep population. After genome annotation and gene enrichment analysis, we found traces of early domestication in sensory organs, behavioural activity and the nervous system as well as adaptive changes in reproductive and wool traits under selection in this population. Some selected genes related to litter size, including FSHR, BMPR1B and OXT, were identified as being under selection. Specific missense mutations of the FSHR gene that differed from the reference genome were also identified in the population, and we found some SNP variations that may affect litter size. Our findings provide a theoretical basis for the conservation and utilization of Yunnan semi-fine wool sheep. Furthermore, our results reveal some changes common to sheep after domestication and provide a new opportunity to investigate the genetic variation influencing fecundity within a population evolving under artificial selection.

Introduction

Yunnan semi-fine wool sheep are among the most important cultivated breeds in Yunnan Province, China. This breed has excellent wool quality, high adaptability and strong robustness. Currently, it is widely used to produce wool and meat products (Wang, 2009). From 1954 to 1971, to improve the wool fineness of ZhaoTong sheep, they were successively crossbred with Rambouillet hybrid sheep, Caucasian sheep, XinJiang fine-wool sheep, and New Zealand Romney sheep, among others. In 1977, Lincoln sheep were introduced to improve wool length, and in 1979, the ideal cross combination was obtained (National Commission On Livestock and Poultry Genetic Resources, 2011). In 2000, this breed was approved by the Chinese Livestock and Poultry Breed Approval Committee, becoming the first approved breed of coarse semi-fine wool sheep in China. Due to the long-term cultivation and growth of these sheep in alpine regions, they have developed high adaptability to cold climates, high altitudes and hypoxic conditions. In general, Yunnan semi-fine wool sheep are medium in body size, and their wool quality is comparable to that of Romney and Lincoln sheep. Due to its competitive advantages and rapid adaptation to the local environment, the population of Yunnan semi-fine wool sheep has been continuously expanded, with the population size exceeding 100,000 individuals in 2005 (Yuan and Sun, 2014). However, its conception rate and lamb survival rate are not high, and more than 90% of ewes give birth to one lamb per parity. Therefore, the improvement of reproductive traits is important for Yunnan semi-fine wool sheep.

Selection signatures are the imprints left on the genome of a species under long-term natural and artificial selection during the process of evolution. Genome selection signatures are an effective tool for studying the adaptability of domestic animals in different natural environments and exploring the genetic mechanisms underlying phenotypic differences. Various methods are currently used to assess genetic diversity and detect selection signatures. These methods can be divided into three categories depending on the different types of examined genomic information: methods based on the allele frequency spectrum, e.g., Tajima’s D (Tajima, 1989) and the composite likelihood ratio (CLR) (Nielsen et al., 2005); methods based on linkage disequilibrium (LD), e.g., extended haplotype homozygosity (EHH) (Sabeti et al., 2002) and the integrated haplotype score (iHS) (Voight et al., 2006); and methods based on population differentiation, e.g., the fixation index (FST) (Weir and Cockerham, 1984). In this study, we employed the iHS and Tajima’s D tests to detect selection signatures within the population.

With the development of sequencing technologies and the improvement of biological analysis methods, selection signature detection has been used to mine the imprints left by domestication and artificial selection and the genes associated with certain important traits (Porto-Neto et al., 2014; Liu et al., 2016; E et al., 2019; Abied et al., 2020; Zhao et al., 2020; Zhong et al., 2020; Zhu et al., 2020). However, due to certain limitations, it is difficult to obtain a complete picture of the germplasm and selection characteristics of indigenous sheep breeds in China. Therefore, the objectives of this study were to explore the population characteristics and genetic structure of Yunnan semi-fine wool sheep, to screen selection regions related to important traits and to identify variants by using genome-wide selection signature detection. Furthermore, the present study will provide a theoretical basis for the breeding and conservation of Yunnan semi-fine wool sheep.

Materials and Methods

Ethics Statement

The ethical committee of the Yunnan Animal Science and Veterinary Institute (Kunming city, Yunnan Province, China) approved all experiments in this study (201909006). In addition, during the study, all authors strictly complied with the Regulations on the Administration of Laboratory Animals (Order No. 2 of the State Science and Technology Commission of the People’s Republic of China, 1988) and the Regulations on the Administration of Experimental Animals of Yunnan Province (the Standing Committee of Yunnan Provincial People’s Congress 2007.10). There was no use of human participants, data or tissues.

Animals and Phenotypes

A total of 40 female Yunnan semi-fine wool sheep with similar body weights (approximately 50 kg) and ages (3 years old) were sampled from the Laishishan sheep farm in Qiaojia City, Yunnan Province (N26°54′55.29″, E102°55′40.02″), which is one representative Yunnan semi-fine wool sheep farm in Yunnan Province, and genetic exchange (rams or sperm) was carried out with other Yunnan semi-fine wool sheep farms. The samples were unrelated according to the pedigree. Based on their litter sizes in two successive parities, these ewes were classified into two groups: ewes with a litter size of 1 (20 individuals) and those with a litter size of 2 (20 individuals). In addition, the sequencing data of 4 Zhaotong sheep and 3 Romney sheep were downloaded from the National Center for Biotechnology Information (NCBI) Sequence Read Archive (SRA) database to further explore the population structure of Yunnan semi-fine wool sheep.

Sample Collection, DNA Extraction, and Whole-Genome Sequencing

Ear tissues of the selected sheep were collected, transferred to sterile centrifuge tubes, and stored in a freezer at −80°C. DNA was extracted from the sample of each individual and analysed following standard experimental procedures to guarantee quality. Sequencing libraries of each individual were constructed by Beijing BerryGenomics Biotechnology Co., Ltd., using the Illumina NovaSeq6000 PE150™ platform.

Genomic Data Processing and SNP Calling

Trimmomatic software (version 0.38) (Bolger et al., 2014) was applied to remove adaptors and low-quality reads. After quality control, the reads of each sample were mapped to the GCF_002742125.1_Oar_rambouillet_v1.0 genome with the ‘mem’ algorithm of BWA software (version 0.7.17) (Li and Durbin, 2009). The bam and sorted bam files were generated using SAMtools (version 1.9) (Li et al., 2009). Then, the Genome Analysis Toolkit (GATK, version 3.7.0) (McKenna et al., 2010) pipeline was used to call SNPs and indels.

To ensure high quality of the variations for the following analysis, BCFtools (version 1.10.2) (Danecek and McCarthy, 2017) was used to filter the minimum or maximum number of alleles (bcftools view -m 2 -M 2) from autosomes and chromosome X. In addition, SNPs with minor allele frequencies lower than 0.01 or missing rates greater than 0.2 were excluded by VCFtools-0.1.16 (Danecek et al., 2011). To identify the functions of the variants, SnpEff 4.3 software (Cingolani et al., 2012) was utilized to annotate the filtered SNPs.

Population Structure Analysis

Based on the filtered SNPs of Yunnan semi-fine wool sheep and their parental sheep breeds, principal component analysis (PCA), admixture analysis and LD decay analysis were performed to explore the population structure of Yunnan semi-fine wool sheep. For PCA, PLINK (v1.90b4) software (Purcell et al., 2007) was used to generate the input file (.bed.bim.fam). gcta64 software (v1.26.0) (Yang et al., 2011) was applied to construct a genetic relationship matrix (–make-grm) and calculate principal components. For admixture analysis, ADMIXTURE (version 1.3.0) (Alexander and Lange, 2011) was employed to carry out unsupervised hierarchical clustering. PopLDdecay software (version 3.41) was used to calculate LD decay, LD was evaluated on the basis of the squared coefficient of correlation (r²) between loci (Zhang et al., 2019), and the parameters were set as follows: -MaxDist 300 -OutType 1. All of the results were presented in graphs created using R software.

Heterozygosity Rate, Inbreeding Coefficients and Effective Population Size

As important population parameters, the heterozygosity rate, inbreeding coefficient and effective population size were calculated in this study. PLINK (v1.90b4) software (Purcell et al., 2007) was used to estimate the heterozygosity rate, which was calculated with the formula $\frac{N\left(NM\right)-O\left(HOM\right)}{N\left(NM\right)}$, where $O\left(HOM\right)$, $E\left(HOM\right)$ and $N\left(NM\right)$ represent the observed number of homozygotes, expected number of homozygotes and number of nonmissing genotypes, respectively, which were obtained by using PLINK. The inbreeding coefficients (F_ROH) of individuals were estimated based on runs of homozygosity (ROHs). PLINK software was used to detect ROHs with the following parameters: 1) minimum density of 1 SNP per 100 kb, 2) minimum ROH length of 1000 kb, 3) 1 heterozygote allowed per window, 4) 5 missing calls allowed per window, 5) minimum number of SNPs in a window of 50, and 6) 0.05 missing calls allowed per window. The following formula was used to calculate F_ROH:

$F_{ROH}=\frac{{{\displaystyle \sum }}_{i}RO{H}_i}{L_{auto}} \left(\mathrm{M}\mathrm{c}\mathrm{Q}\mathrm{u}\mathrm{i}\mathrm{l}\mathrm{l}\mathrm{a}\mathrm{n}\mathrm{e}\mathrm{t}\mathrm{a}\mathrm{l}.,2008\right)$

where ${\displaystyle \sum }_{i}RO{H}_i$ is the total ROH length on autosomes and $L_{auto}$ is the total length of autosomes. The effective population size (Ne) can indicate dynamic changes in population size, and Ne was estimated according to the random mating model based on the heterozygote excess (HE) method using the default parameters in NeEstimator (version 2.1) software (Do et al., 2014). NeEstimator calculates Ne as follows:

$Ne=\frac{1}{2D}+\frac{1}{2\left(D+1\right)} \left(\mathrm{P}\mathrm{u}\mathrm{d}\mathrm{o}\mathrm{v}\mathrm{k}\mathrm{i}\mathrm{n}\mathrm{e}\mathrm{t}\mathrm{a}\mathrm{l}.,1996\right)$

where Ne is the effective population size and D is Selander’s index, calculated as follows: $D_j\left(i\right)=\frac{H_j^{obs}\left(i\right)-H_j^{\exp }\left(i\right)}{H_j^{\exp }\left(i\right)}$. $H_j^{\exp }\left(i\right)$ is the expected heterozygote frequency of allele i at locus j. The actual D value is the weighted mean value, calculated as follows: $D=\frac{{\displaystyle \sum w_j{D}_j}}{{\displaystyle \sum _{j=1}^k\left(n_j-1\right)\sqrt{N_j}}}$, where w_j is the weight of locus j and is calculated as $w_j=\sqrt{N_j}\frac{n_j-1}{n_j},$ where $N_j$ is the sample size, and $n_j$ is the number of alleles at locus j.

Detection of Selection Signatures

Selection signature detection within populations can reveal selection dynamics and the history of evolution. In this study, two metrics, Tajima’s D (Tajima, 1989) and the iHS (Voight et al., 2006), were utilized to identify selection signatures in the whole genome of Yunnan semi-fine wool sheep. Tajima’s D is based on neutral mutation theory. When Tajima’s D is equal to 0, it indicates that the population is in a neutral evolutionary state, whereas Tajima’s D values lower or greater than 0 indicate that the population has experienced purifying selection events/population expansion or balancing selection, respectively. Tajima’s D values were calculated in nonoverlapping 50-kb sliding windows. Only the windows with the 1% highest and the 1% lowest Tajima’s D values were identified as subject to selection. Different from the single locus-based Tajima’s D metric, the iHS is mainly based on haplotypes. A negative iHS value indicates that the mutated allele may be affected by positive selection. In this study, haplotypes were first constructed by using SHAPEIT (Delaneau and Marchini, 2014), and iHS statistics were then calculated using Selscan software (version 1.2.0) (Szpiech et al., 2014). Thereafter, the iHS values were standardized to follow a standard normal distribution, and the statistics exceeding the 0.1% quantile ($\left|Z\left(iHS\right) value\right|=3$) indicated that the core SNP experienced selection.

Functional Annotation of Selected Regions

Based on the core SNPs or windows under selection detected based on the iHS and Tajima’s D analyses, selected regions were defined for the bioinformatics analysis. From each core SNP detected based on iHS analysis, the region was extended 20 kb upstream and downstream to define the selected region. Under the Tajima’s D approach, each selected 50-kb sliding window was extended 150 kb upstream and downstream to define the selected region. Gene information in the selected regions was obtained using the Biomart database (http://asia.ensembl.org/biomart/martview/). Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway and Gene Ontology (GO) enrichment analyses were performed for further gene function analysis using DAVID Bioinformatics Resources 6.8 (Huang et al., 2009) (https://david.ncifcrf.gov/). The GO terms and KEGG pathways with P values less than 0.05 were considered significant.

Mutation Analysis of Selected Genes Relevant to Reproductive Traits

In this study, 40 Yunnan semi-fine wool sheep were classified into two groups on the basis of litter size: 1 lamb (20 individuals) or 2 lambs (20 individuals) in two successive parities. To investigate whether SNP mutations or frequency differences existed between the two groups, three significantly selected genes (FSHR, BMPR1B, and OXT) related to reproductive traits were studied. PLINK (v1.90b4) software was employed to perform a standard case (individuals with two lambs)/control (individuals with one lamb) association analysis using Fisher’s exact test with the option plink --file gene --fisher.

Results

Genetic Variants in Yunnan Semi-Fine Wool Sheep

In the present study, a total of 40 female Yunnan semi-fine wool sheep were sequenced. After filtering out contaminated reads, including adaptors, low-quality reads and reads with an N ratio greater than 10%, 1393 Gb of clean data were acquired, and 35,178,304 raw SNPs and 5,256,453 raw indels were obtained. The mapping rate against the reference genome was 91.23% on average (86.01%–92.26%), and the average sequence depth was 9.51X (8.09X∼12.10X) (see Supplementary Material S1 for details). After filtering, 28,593,198 SNPs and 4,725,259 indels with high quality were retained.

Following annotation, 27,623,037 SNPs on autosomes and 970,161 SNPs on chromosome X (Supplementary Material S2) were obtained. These SNPs were partitioned according to their locations: intergenic (15,897,566; 56.03%), intronic (9,491,317; 33.45%), exonic (177516; 0.63%) and other gene regulatory regions (Figure 1A). Among the SNPs in exonic regions, as shown in Figure 1B, synonymous mutations constituted the overwhelming majority (60.07%), followed by missense mutations (39.14%), and the other types of observed mutations included stop-gain, stop-loss, and start-loss mutations (less than 1%). The transition/transversion (Ts/Tv) ratio was 2.4721.

FIGURE 1

FIGURE 1. SNP annotation of Yunnan semi-fine wool sheep. (A) Genomic annotation of SNPs according to SnpEff. “Others” contained a few other “splice site acceptor,” “splice site donor” and “splice site region” sites. (B) Pie chart of SNPs annotated in exonic regions.

Genetic Diversity of the Yunnan Semi-Fine Wool Sheep Population

To describe the genetic diversity of the Yunnan semi-fine wool sheep population, two metrics were used. For the whole population, the average heterozygosity rate and inbreeding coefficient were 0.264 (0.231–0.374) and 0.0099 (0.0004–0.0239), respectively (see Supplementary Material S4 for details), indicating that the heterozygosity rates and inbreeding coefficients varied greatly within the population. The effective population size estimated from the HE of the Yunnan semi-fine wool sheep population using the lowest allele frequency of 0.05 as the threshold was 242.9 (95% CI: 232.6–254.2).

Population Structure Analysis

The population structure is illustrated in Figure 2. ADMIXTURE analysis indicated that the ancestors of Yunnan semi-fine wool sheep included Romney and Zhaotong sheep, with Romney sheep contributing more (K = 2); Figure 2A also demonstrates other ancestors of Yunnan semi-fine wool sheep (K = 3). Moreover, PCA (Figure 2B) further showed that the population of Yunnan semi-fine wool sheep was aggregated and separated from Romney and Zhaotong sheep, indicating that it is a new breed. LD decay analysis of two groups with litter sizes of 1 or 2 lambs (Figure 2C) confirmed the absence of obvious population stratification within Yunnan semi-fine wool sheep, and the difference between litter sizes might be due to ongoing selection.

FIGURE 2

FIGURE 2. Genetic differentiation of 40 Yunnan semi-fine wool sheep. (A) ADMIXTURE analysis of Yunnan semi-fine wool sheep and their parental sheep breeds (Romney and Zhaotong sheep). (B) PCA of 40 Yunnan semi-fine wool sheep and their parental sheep breeds using 710,412 common SNPs located on autosomes. (C) LD decay curves (based on r²) of two groups (single or double lambs).

Detection of Genome-Wide Signatures of Selection

As shown in Figure 3, regardless of whether the iHS or Tajima’s D approach was used, the observed selection signatures were distributed across all chromosomes. Under the iHS approach, the 0.1% quantile ($\left|Z\left(iHS\right) value\right|=3$) was used as the threshold to extract selected core SNPs. A total of 11,356 core SNPs showing the strongest selection were identified, accounting for 1% of all selected SNPs. The iHS values ranged from −7.69 to 8.26 (Figures 3A,C). Considering the LD among SNPs, the selected SNPs within each 50-kb interval were combined into one selected region. In total, 3512 selected regions were identified. Chromosomes 1, 2, 3 and 4 harboured the most selected regions, i.e., 422, 323, 297 and 206 selected regions, respectively (Figure 3F).

FIGURE 3

FIGURE 3. Selective sweep analysis of 40 Yunnan semi-fine wool sheep. (A) Distribution of iHS values for selected regions. (B) Distribution of Tajima’s D values for selected regions. (C) Manhattan plot of $Z\left(iHS\right)$ values of sheep. The dashed line denotes a threshold of $\left|Z\left(iHS\right) value\right|$ = 3 (p = 0.001). (D) Manhattan plot of Tajima’s D values. The black line is the “suggestive line”: the top 1% of values are >2.97345, and the bottom 1% of values are < −1.63182. (E) The 108 common selected genes detected by both the iHS and Tajima’s D approaches. (F) The distribution of the selected SNPs on each chromosome obtained by the iHS and Tajima’s D approaches. In the figure, “27” on the X-axis represents chromosome X.

A total of 56,115 nonoverlapping 50-kb windows were assessed, and only the top and bottom 1% of the windows with high Tajima’s D values were determined to be regions of selective sweeps (Figures 3B,D). A total of 562 windows showing the strongest selection and 449 equivalent selected regions were identified by combining the selected windows. Selected regions were identified on all chromosomes (Figure 3F). Among all chromosomes, chromosomes 1, 2, 3 and X harboured the most selected windows, i.e., 42, 41, 41 and 41 windows, respectively. The most significant regions were located on chromosome 20 from 16.70–16.85 Mb, with an average Tajima’s D value of −2.77432, and this region harboured 495 SNPs. Five known genes were identified within this putative selective sweep region, namely, CUL7, KLC4, MEA1, MRPL2 and PTK7, which may regulate cell activity and microtubule motor activity.

Functional Enrichment of Candidate Genes

According to the two methods of selection signature detection, 1856 genes were found in the 20-kb regions upstream and downstream of the core SNPs according to the iHS test, and 1753 selected genes were obtained in the 350-kb intervals upstream and downstream of the selected windows based on the Tajima’s D test. We combined the genes obtained via the two methods described above, and a total of 108 genes were identified (Figure 3E). To evaluate the functions of the selected genes obtained by the two methods, gene annotation and enrichment analyses (GO and KEGG) were conducted.

For the iHS, the enrichment analysis identified a total of 53 significant GO terms (Supplementary Material S3; Figure 1), namely, 13 cellular component terms, 17 molecular function terms, and 23 biological process terms. Notably, most of the selected genes were enriched in the categories of immune system processes, intracellular signal transduction and neural synapses. In the KEGG analysis results, 19 pathways were annotated, as shown in Figure 4A. Some known reproduction-related, domestication-related and wool-related pathways were found to be significantly enriched, including the vascular smooth muscle contraction, calcium signalling, salivary secretion, gastric acid secretion, long-term depression and long-term potentiation pathways.

FIGURE 4

FIGURE 4. KEGG enrichment analysis of selected genes. (A) KEGG analysis of the selected genes based on the (A) iHS and (B) Tajima’s D approaches.

Among the results based on Tajima’s D method, 18 significant GO terms were identified (Supplementary Material S3; Figure 2), including olfactory receptor activity, immune process and dendrite morphogenesis terms. A total of 12 KEGG pathways were obtained (Figure 4B), including folate biosynthesis related to reproductive traits and the GABAergic synapse pathway related to inhibition of the nervous system.

Overall, our results revealed traces of early domestication in Yunnan semi-fine wool sheep, such as modifications in sensory organs, behavioural activity and the nervous system, as well as adaptive changes in economic traits, including reproductive and wool traits, under artificial selection (Table 1).

TABLE 1

TABLE 1. GO terms or KEGG pathways associated with domestication, reproduction and immune system processes.

Mutation and Association Analysis of Genes Related to Reproduction

Genes related to reproductive traits, such as follicle-stimulating hormone receptor (FSHR), bone morphogenetic protein receptor type 1B (BMPR1B) and oxytocin/neurophysin I prepropeptide (OXT), were identified as being under selection in this study. OXT was identified by the Tajima’s D test, and BMPR1B and FSHR were identified by the iHS test. Based on the whole-genome sequencing results obtained for Yunnan semi-fine wool sheep, the polymorphic loci in the exons of the BMPR1B, FSHR and OXT genes were analysed by comparison with the reference genome (Table 2).

TABLE 2

TABLE 2. Mutations in the BMPR1B and FSHR genes in exonic regions of the Yunnan semi-fine wool sheep genome.

In the FSHR gene, a total of 2437 SNPs were identified. Six mutations were detected in its exonic region (Table 2), four of which were missense mutations (c.28G < A, c.692G > A, c.783T > G, and c.854C > T, known mutation loci). At the same time, we calculated the distribution of minor alleles in the two groups based on litter size. Twenty-three SNPs with significant differences as detected by Fisher’s exact test were finally obtained (Figure 5A). Although these SNPs were all located in intronic regions, they are probably related to reproductive performance in the two groups.

FIGURE 5

FIGURE 5. SNPs of the FSHR and BMPR1B genes differing between the single- and double-lambing groups identified using Fisher’s exact test. Manhattan plot of (A) FSHR and (B) BMPR1B SNPs differing between the single- and double-lambing groups. The red line indicates p = 0.05. The location of the top 10 significant loci on the (C) FSHR and (D) BMPR1B gene structures.

In the BMPR1B gene, among the 2645 identified SNPs, 7 SNP mutations were located in the exonic region, only one of which was a missense variant (c.52G > A, known mutation locus, rs605658565), and 44 SNPs showed significant differences between the single- and double-lambing groups (Figure 5B). In the OXT gene, there were only 3 identified SNPs (Chr13:54216488, 54216953, and 54216953). Although none of the SNPs showed significant differences, we found that in the two-lamb group, the frequency of minor alleles increased de novo.

Discussion

Population Situation of Yunnan Semi-Fine Wool Sheep

Yunnan semi-fine wool sheep were the first coarse-grade semi-fine wool sheep breed cultivated in China. In the last 2 decades, the breed has become an indispensable component of sheep germplasm resources in China. This study is the first to explore the population characteristics, genetic structure and selection signatures based on sequencing data of this breed. We used the heterozygosity rate and inbreeding coefficient to evaluate genetic diversity within the population. ROHs are continuous homozygous segments, and they can reflect the inbreeding level of a population at the genome level. Using ROHs to calculate the inbreeding coefficient (F_ROH) is more accurate than estimating the inbreeding coefficient from pedigree data (Kim et al., 2015). In this study, the average inbreeding coefficient (0.0099: 0.0004–0.0239) was lower than that of other sheep breeds, e.g., Chaka sheep (F_ROH = 0.032) (Cheng et al., 2020), Italian Alpagota sheep (F_ROH = 0.053) (Mastrangelo et al., 2018)) and other species, e.g., white leghorn (F_ROH = 0.28) (Zhang et al., 2020) and Laiwu pigs (F_ROH = 0.133) (Fang et al., 2021), indicating that the population structure of Yunnan semi-fine wool sheep is ideal. Combining the average inbreeding coefficient and heterozygosity rate, we could preliminarily infer that although these two parameters varied greatly within the Yunnan semi-fine wool sheep population, sufficient genetic diversity was maintained overall during cultivation and artificial selection. This conclusion is well supported by the Ne results. In the long-term selection process, the greater the selection intensity is, the smaller the effective population size and the lower the genetic diversity of the population. According to previous studies, when the effective population size is below 50 individuals, the population is threatened (Meuwissen and Woolliams, 1994). According to our results, the Ne of Yunnan semi-fine wool sheep is large enough to indicate a low risk of excessive inbreeding at present.

The breeding of Yunnan semi-fine wool sheep is carried out under a strict system: new genetic materials cannot be introduced into the core group, and sheep in propagation groups must be obtained from foundation seed farms, which maintains the good breeding status of these sheep. However, the herd size in breeding farms is small, e.g., only 1010 and 270 sheep in the Lashishan breeding farm and Xiaohai breeding farm in Qiaojia County, respectively, both of which are representative breeding farms of Yunnan semi-fine wool sheep in Yunnan Province. Such a small population size will lead to inbreeding, resulting in a decrease in the effective population size.

Genomic Signatures of Selection in Yunnan Semi-Fine Wool Sheep

Genetic variation profoundly affects phenotypic variation. When either natural or artificial selection occurs, it will leave traces in the genome. Therefore, we adopted intrapopulation detection methods (Tajima’s D and iHS) to discover genome-wide footprints caused by natural and artificial selection in Yunnan semi-fine wool sheep. In this study, overlap of the selected genes simultaneously identified by the two approaches was observed in a few cases (108 genes, Figure 3E). The iHS approach shows higher efficiency in detecting partial sweeps; thus, it can detect an event in which a favourable allele increases rapidly from a low frequency but has not yet reached a fixed state (Pritchard et al., 2010). The Tajima’s D test is especially powerful for the detection of fixation signatures. According to the characteristics of the above two signature detection methods, the two approaches can be employed together to identify more traces of selection in the genome. Considering the short breeding history of Yunnan semi-fine wool sheep, positive selection may still be affecting some genes related to wool, meat quality and reproductive traits. Furthermore, some alleles for genes associated with domestication and adaptation may be fixed.

Domestication refers to the process by which wild animals come to be maintained under domestic conditions, in which they reproduce over generations and are used by humans. Domestic animals undergo fundamental changes in their phenotypes, morphology, and behaviour relative to those of wild animals (Naval-Sanchez et al., 2018). Some key traits were selected and fixed in most domestic sheep breeds in the early stages of domestication. The changes in these traits mainly manifest as reductions in brain volume and weight, which lead to dulled sensory organ function, more docile personalities and slower activities, such as vision, smell and motor abilities (Kruska, 1996). In this study, we identified some biological pathways associated with nervous system regulation and olfactory receptor activity, consistent with the results of previous studies (Axelsson et al., 2013; Li et al., 2020).

Due to the spread of domestic sheep with human migration activities worldwide, some traits associated with specific human needs have been fixed, leading to greater diversity of some phenotypes. In the Yunnan semi-fine wool sheep population investigated in this study, in addition to the findings related to domestication, we also found many selected regions and novel functional genes that may be responsible for traits such as reproduction and wool production. Some biological pathways associated with reproductive traits, including vascular smooth muscle contraction and folate biosynthesis, were identified, which may be related to the development of the endometrial vascular system and pregnancy (Cullinan-Bove et al., 1993).

In previous studies on a segregated flock based on QTL analysis and GWAS mapping, some mutations such as FecB (Mulsant et al., 2001), FecX (Galloway et al., 2000), and FecG (Nicol et al., 2009) that may affect ovulation in sheep were identified. Interestingly, we also found that the BMPR1B gene was located in a selected region in this study. The BMPR1B gene encodes a member of the bone morphogenetic protein (BMP) receptor family of transmembrane serine/threonine kinases and regulates animal cell growth and differentiation, embryonic development and reproductive performance (Wilson et al., 2001). Since the mutant FecB [nonsynonymous substitution (Q249R)] form of the sheep BMPR1B gene was proven to be the main gene regulating high-fertility performance in Booroola sheep, researchers have attempted to identify the FecB gene in some indigenous breeds such small-tail Han and Hu sheep (Chu et al., 2011). However, some researchers have reported no significant link between this gene and high fertility traits in breeds such as Tan sheep (Chong et al., 2019). In this study, the FecB mutation was not detected in Yunnan semi-fine wool sheep, suggesting that this BMPR1B gene mutation may not affect the reproductive performance of the breed. Of course, it is also possible that the examined sample size was too small for detection of this mutation in this study. Nevertheless, we detected 44 SNPs with significantly different frequencies between the single- and double-lambing groups. However, these SNPs were mainly located in intronic regions. Many mutations in intronic regions have previously been found to (Figure 5D) lead to functional changes through aberrant splicing (Busslinger et al., 1981; Vaz-Drago et al., 2017), so it is necessary to expand the studied population or conduct experimental verification.

Two other genes that participate in reproduction-related hormone secretion are FSHR and OXT. The OXT gene encodes a precursor protein that is processed to produce oxytocin and neurophysin I. This precursor seems to be activated as it is transported along the axon to the posterior pituitary. This hormone causes the contraction of smooth muscle during parturition and lactation and functions as a neurotransmitter in the central nervous system, playing a role in cognition, tolerance, adaptation, and complex sexual and maternal behaviours (Gimpl and Fahrenholz, 2001). In the present study, there was little evidence that this gene was associated with litter size, but we did find differences in SNPs in this gene region between the two litter size groups, suggesting that OXT may be related to litter size. Most relevant studies have shown that the FSHR gene is closely related to reproductive traits, such as ovarian function (Laan et al., 2012) and testicular development (Lend et al., 2010), in mammals. SNPs in the 10th exon of FSHR are significantly correlated with fertility traits of the Chinese native pig breed Xiaomeishan (Wu and Wang, 2012). Mutations in the 5′ flanking region of the ovine FSHR gene may also be associated with litter size (Pan et al., 2014). In this study, we identified specific missense mutations in the Yunnan semi-fine wool sheep population that differed from the reference genome. However, their functions still need to be further explored. On this basis, 23 SNPs with significantly different frequencies were also identified between the single- and double-lambing groups. All of these SNPs were located in intronic regions (Figure 5C), but their effect on lambing size needs to be further investigated.

Conclusion

In conclusion, using second-generation sequencing technology, the population structure of Yunnan semi-fine wool sheep was detected, and strategies for the conservation of this breed were proposed. In addition, some genes related to adaptation and reproductive traits were shown to have experienced strong selection. The results of this study increase our knowledge of the genetic basis of litter size in Yunnan semi-fine wool sheep, shed light on the changes in heritable phenotypes during the processes of adaptation and artificial breeding, and provide a theoretical basis for the breeding and conservation of Yunnan semi-fine wool sheep.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: NCBI, PRJNA783661.

Ethics Statement

The animal study was reviewed and approved by the ethical committee of the Yunnan Animal Science and Veterinary Institute (201909006).

Author Contributions

YG, GQ, and XD conceived the study. JL and CL collected the samples and recorded the phenotypes. YG and YW completed the DNA data analysis. JL, CL, and GW contributed to the visualization of the data. GQ and XD supervised the study and proposed revisions to the manuscript. YG, GQ, and XD wrote and revised the manuscript. All authors read and approved the manuscript.

Funding

This investigation was funded by the National Natural Science Foundation of China (Grant No. 31960663), the National Wool Caprine Industrial Technology System (Grant No. CARS-39), Yunnan Applied Basic Research Projects (Grant No. 202001AS070001), the Yunnan Young Academic Leaders Program (Grant No. 202005AC160004), and the Major Science and Technology Project of Yunnan Province (Grant No. 202102AE090039). This work was also supported by grants from the National Key Research and Development Project (2019YFE0106800).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

The authors gratefully acknowledge the reviewers for their constructive comments.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2022.812753/full#supplementary-material

References

Abied, A., Bagadi, A., Bordbar, F., Pu, Y., Augustino, S. M. A., and Xue, X. (2020). Genomic Diversity, Population Structure, and Signature of Selection in Five Chinese Native Sheep Breeds Adapted to Extreme Environments. Genes 11 (5), 494. doi:10.3390/genes11050494

PubMed Abstract | CrossRef Full Text | Google Scholar

Alexander, D. H., and Lange, K. (2011). Enhancements to the ADMIXTURE Algorithm for Individual Ancestry Estimation. BMC Bioinformatics 12, 246. doi:10.1186/1471-2105-12-246

PubMed Abstract | CrossRef Full Text | Google Scholar

Axelsson, E., Ratnakumar, A., Arendt, M.-L., Maqbool, K., Webster, M. T., Perloski, M., et al. (2013). The Genomic Signature of Dog Domestication Reveals Adaptation to a Starch-Rich Diet. Nature 495 (7441), 360–364. doi:10.1038/nature11837

PubMed Abstract | CrossRef Full Text | Google Scholar

Bolger, A. M., Lohse, M., and Usadel, B. (2014). Trimmomatic: a Flexible Trimmer for Illumina Sequence Data. Bioinformatics 30 (15), 2114–2120. doi:10.1093/bioinformatics/btu170

PubMed Abstract | CrossRef Full Text | Google Scholar

Busslinger, M., Moschonas, N., and Flavell, R. A. (1981). Beta + Thalassemia: Aberrant Splicing Results from a Single point Mutation in an Intron. Cell 27 (2 Pt 1), 289–298. doi:10.1016/0092-8674(81)90412-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Cheng, J., Zhao, H., Chen, N., Cao, X., Hanif, Q., Pi, L., et al. (2020). Population Structure, Genetic Diversity, and Selective Signature of Chaka Sheep Revealed by Whole Genome Sequencing. BMC Genomics 21 (1), 520. doi:10.1186/s12864-020-06925-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Chong, Y., Liu, G., and Jiang, X. (2019). Effect of BMPRIB Gene on Litter Size of Sheep in China: A Meta-Analysis. Anim. Reprod. Sci. 210, 106175. doi:10.1016/j.anireprosci.2019.106175

PubMed Abstract | CrossRef Full Text | Google Scholar

Chu, M., Jia, L., Zhang, Y., Jin, M., Chen, H., Fang, L., et al. (2011). Polymorphisms of Coding Region of BMPR-IB Gene and Their Relationship with Litter Size in Sheep. Mol. Biol. Rep. 38 (6), 4071–4076. doi:10.1007/s11033-010-0526-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Cingolani, P., Platts, A., Wang, L. L., Coon, M., Nguyen, T., and Wang, L. (2012). A Program for Annotating and Predicting the Effects of Single Nucleotide Polymorphisms, SnpEff. Fly 6 (2), 80–92. doi:10.4161/fly.19695

PubMed Abstract | CrossRef Full Text | Google Scholar

Cullinan-Bove, K., and Koos, R. D. (1993). Vascular Endothelial Growth Factor/vascular Permeability Factor Expression in the Rat Uterus: Rapid Stimulation by Estrogen Correlates with Estrogen-Induced Increases in Uterine Capillary Permeability and Growth. Endocrinology 133 (2), 829–837. doi:10.1210/endo.133.2.8344219

PubMed Abstract | CrossRef Full Text | Google Scholar

Danecek, P., Auton, A., Abecasis, G., Albers, C. A., Banks, E., and DePristo, M. A. (2011). The Variant Call Format and VCFtools. Bioinformatics 27 (15), 2156–2158. doi:10.1093/bioinformatics/btr330

PubMed Abstract | CrossRef Full Text | Google Scholar

Danecek, P., and McCarthy, S. A. (2017). BCFtools/Csq: Haplotype-Aware Variant Consequences. Bioinformatics 33 (13), 2037–2039. doi:10.1093/bioinformatics/btx100

PubMed Abstract | CrossRef Full Text | Google Scholar

Delaneau, O., Marchini, J., and Marchini, J. (2014). Integrating Sequence and Array Data to Create an Improved 1000 Genomes Project Haplotype Reference Panel. Nat. Commun. 5, 3934. doi:10.1038/ncomms4934

PubMed Abstract | CrossRef Full Text | Google Scholar

Do, C., Waples, R. S., Peel, D., Macbeth, G. M., Tillett, B. J., and Ovenden, J. R. (2014). NeEstimatorv2: Re-implementation of Software for the Estimation of Contemporary Effective Population Size (Ne) from Genetic Data. Mol. Ecol. Resour. 14 (1), 209–214. doi:10.1111/1755-0998.12157

PubMed Abstract | CrossRef Full Text | Google Scholar

E, G.-X., Duan, X.-H., Zhang, J.-H., Huang, Y.-F., Zhao, Y.-J., Na, R.-S., et al. (2019). Genome-wide Selection Signatures Analysis of Litter Size in Dazu Black Goats Using Single-Nucleotide Polymorphism. Biotech. 9 (9), 336. doi:10.1007/s13205-019-1869-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Fang, Y., Hao, X., Xu, Z., Sun, H., Zhao, Q., Cao, R., et al. (2021). Genome-Wide Detection of Runs of Homozygosity in Laiwu Pigs Revealed by Sequencing Data. Front. Genet. 12, 629966. doi:10.3389/fgene.2021.629966

PubMed Abstract | CrossRef Full Text | Google Scholar

Galloway, S. M., McNatty, K. P., Cambridge, L. M., Laitinen, M. P., Juengel, J. L., Jokiranta, T. S., et al. (2000). Mutations in an Oocyte-Derived Growth Factor Gene (BMP15) Cause Increased Ovulation Rate and Infertility in a Dosage-Sensitive Manner. Nat. Genet. 25 (3), 279–283. doi:10.1038/77033

PubMed Abstract | CrossRef Full Text | Google Scholar

Gimpl, G., and Fahrenholz, F. (2001). The Oxytocin Receptor System: Structure, Function, and Regulation. Physiol. Rev. 81 (2), 629–683. doi:10.1152/physrev.2001.81.2.629

PubMed Abstract | CrossRef Full Text | Google Scholar

Huang, D. W., Sherman, B. T., and Lempicki, R. A. (2009). Systematic and Integrative Analysis of Large Gene Lists Using DAVID Bioinformatics Resources. Nat. Protoc. 4 (1), 44–57. doi:10.1038/nprot.2008.211

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, E.-S., Sonstegard, T. S., Van Tassell, C. P., Wiggans, G., and Rothschild, M. F. (2015). The Relationship between Runs of Homozygosity and Inbreeding in Jersey Cattle under Selection. PLoS One 10 (7), e0129967. doi:10.1371/journal.pone.0129967

PubMed Abstract | CrossRef Full Text | Google Scholar

Kruska, D. (1996). The Effect of Domestication on Brain Size and Composition in the Mink ( Mustela Vison ). J. Zoolog. 239 (4), 645–661. doi:10.1111/j.1469-7998.1996.tb05468.x

CrossRef Full Text | Google Scholar

Laan, M., Grigorova, M., and Huhtaniemi, I. T. (2012). Pharmacogenetics of Follicle-Stimulating Hormone Action. Curr. Opin. Endocrinol. Diabetes Obes. 19 (3), 220–227. doi:10.1097/MED.0b013e3283534b11

PubMed Abstract | CrossRef Full Text | Google Scholar

Lend, A. K., Belousova, A., Haller-Kikkatalo, K., Punab, M., Poolamets, O., Peters, M., et al. (2010). Follicle-Stimulating Hormone Receptor Gene Haplotypes and Male Infertility in Estonian Population and Meta-Analysis. Syst. Biol. Reprod. Med. 56 (1), 84–90. doi:10.3109/19396360903456676

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, H., and Durbin, R. (2009). Fast and Accurate Short Read Alignment with Burrows-Wheeler Transform. Bioinformatics 25 (14), 1754–1760. doi:10.1093/bioinformatics/btp324

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, H., Handsaker, B., Wysoker, A., Fennell, T., Ruan, J., and Homer, N. (2009). The Sequence Alignment/Map Format and SAMtools. Bioinformatics 25 (16), 2078–2079. doi:10.1093/bioinformatics/btp352

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, X., Yang, J., Shen, M., Xie, X. L., Liu, G. J., Xu, Y. X., et al. (2020). Whole-genome Resequencing of Wild and Domestic Sheep Identifies Genes Associated with Morphological and Agronomic Traits. Nat. Commun. 11 (1), 2815. doi:10.1038/s41467-020-16485-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, Z., Sun, C., Qu, L., Wang, K., and Yang, N. (2016). Genome-Wide Detection of Selective Signatures in Chicken through High Density SNPs. PLoS One 11 (11), e166146. doi:10.1371/journal.pone.0166146

PubMed Abstract | CrossRef Full Text | Google Scholar

Mastrangelo, S., Ciani, E., Sardina, M. T., Sottile, G., Pilla, F., and Portolano, B. (2018). Runs of Homozygosity Reveal Genome-wide Autozygosity in Italian Sheep Breeds. Anim. Genet. 49 (1), 71–81. doi:10.1111/age.12634

PubMed Abstract | CrossRef Full Text | Google Scholar

McKenna, A., Hanna, M., Banks, E., Sivachenko, A., Cibulskis, K., Kernytsky, A., et al. (2010). The Genome Analysis Toolkit: A MapReduce Framework for Analyzing Next-Generation DNA Sequencing Data. Genome Res. 20 (9), 1297–1303. doi:10.1101/gr.107524.110

PubMed Abstract | CrossRef Full Text | Google Scholar

McQuillan, R., Leutenegger, A. L., Abdel-Rahman, R., Franklin, C. S., Pericic, M., Barac-Lauc, L., et al. (2008). Runs of Homozygosity in European Populations. Am. J. Hum. Genet. 83 (3), 359–372. doi:10.1016/j.ajhg.2008.08.007

CrossRef Full Text | Google Scholar

Meuwissen, T. H., and Woolliams, J. A. (1994). Effective Sizes of Livestock Populations to Prevent a Decline in Fitness. Theor. Appl. Genet. 89 (7-8), 1019–1026. doi:10.1007/BF00224533

PubMed Abstract | CrossRef Full Text | Google Scholar

Mulsant, P., Lecerf, F., Fabre, S., Schibler, L., Monget, P., Lanneluc, I., et al. (2001). Mutation in Bone Morphogenetic Protein Receptor-IB Is Associated with Increased Ovulation Rate in Booroola Merino Ewes. Proc. Natl. Acad. Sci. U S A. 98 (9), 5104–5109. doi:10.1073/pnas.091577598

PubMed Abstract | CrossRef Full Text | Google Scholar

National Commission On Livestock and Poultry Genetic Resources (2011). Animal Genetic Resources in China: Sheep and Goats. Beijing: China Agriculture Press.

Google Scholar

Naval-Sanchez, M., Nguyen, Q., McWilliam, S., Porto-Neto, L. R., Tellam, R., Vuocolo, T., et al. (2018). Sheep Genome Functional Annotation Reveals Proximal Regulatory Elements Contributed to the Evolution of Modern Breeds. Nat. Commun. 9 (1), 859. doi:10.1038/s41467-017-02809-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Nicol, L., Bishop, S. C., Pong-Wong, R., Bendixen, C., Holm, L. E., Rhind, S. M., et al. (2009). Homozygosity for a Single Base-Pair Mutation in the Oocyte-specific GDF9 Gene Results in Sterility in Thoka Sheep. Reproduction 138 (6), 921–933. doi:10.1530/REP-09-0193

PubMed Abstract | CrossRef Full Text | Google Scholar

Nielsen, R., Williamson, S., Kim, Y., Hubisz, M. J., Clark, A. G., and Bustamante, C. (2005). Genomic Scans for Selective Sweeps Using SNP Data. Genome Res. 15 (11), 1566–1575. doi:10.1101/gr.4252305

PubMed Abstract | CrossRef Full Text | Google Scholar

Pan, X., Liu, S., Li, F., Wang, W., Li, C., Ma, Y., et al. (2014). Molecular Characterization, Expression Profiles of the Ovine FSHR Gene and its Association with Litter Size. Mol. Biol. Rep. 41 (12), 7749–7754. doi:10.1007/s11033-014-3666-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Porto-Neto, L. R., Lee, S. H., Sonstegard, T. S., Van Tassell, C. P., Lee, H. K., Gibson, J. P., et al. (2014). Genome-wide Detection of Signatures of Selection in Korean Hanwoo Cattle. Anim. Genet. 45 (2), 180–190. doi:10.1111/age.12119

PubMed Abstract | CrossRef Full Text | Google Scholar

Pritchard, J. K., Pickrell, J. K., and Coop, G. (2010). The Genetics of Human Adaptation: Hard Sweeps, Soft Sweeps, and Polygenic Adaptation. Curr. Biol. 20 (4), R208–R215. doi:10.1016/j.cub.2009.11.055

PubMed Abstract | CrossRef Full Text | Google Scholar

Pudovkin, A. I., Zaykin, D. V., and Hedgecock, D. (1996). On the Potential for Estimating the Effective Number of Breeders from Heterozygote-Excess in Progeny. Genetics 144 (1), 383–387. doi:10.1093/genetics/144.1.383

PubMed Abstract | CrossRef Full Text | Google Scholar

Purcell, S., Neale, B., Todd-Brown, K., Thomas, L., Ferreira, M. A., Bender, D., et al. (2007). PLINK: A Tool Set for Whole-Genome Association and Population-Based Linkage Analyses. Am. J. Hum. Genet. 81 (3), 559–575. doi:10.1086/519795

CrossRef Full Text | Google Scholar

Sabeti, P. C., Reich, D. E., Higgins, J. M., Levine, H. Z., Richter, D. J., Schaffner, S. F., et al. (2002). Detecting Recent Positive Selection in the Human Genome from Haplotype Structure. Nature 419 (6909), 832–837. doi:10.1038/nature01140

PubMed Abstract | CrossRef Full Text | Google Scholar

Szpiech, Z. A., and Hernandez, R. D. (2014). Selscan: An Efficient Multithreaded Program to Perform EHH-Based Scans for Positive Selection. Mol. Biol. Evol. 31 (10), 2824–2827. doi:10.1093/molbev/msu211

PubMed Abstract | CrossRef Full Text | Google Scholar

Tajima, F. (1989). Statistical Method for Testing the Neutral Mutation Hypothesis by DNA Polymorphism. Genetics 123 (3), 585–595. doi:10.1093/genetics/123.3.585

PubMed Abstract | CrossRef Full Text | Google Scholar

Vaz-Drago, R., Custodio, N., and Carmo-Fonseca, M. (2017). Deep Intronic Mutations and Human Disease. Hum. Genet. 136 (9), 1093–1111. doi:10.1007/s00439-017-1809-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Voight, B. F., Kudaravalli, S., Wen, X., and Pritchard, J. K. (2006). A Map of Recent Positive Selection in the Human Genome. Plos Biol. 4 (3), e72. doi:10.1371/journal.pbio.0040072

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, W. (2009). Conservation and Utilization of Yunnan Semi-fine Wool Sheep. Yunnan J. Anim. Sci. Vet. Med. (S1), 42–43.

Google Scholar

Weir, B. S., and Cockerham, C. C. (1984). Estimating F-Statistics for the Analysis of Population Structure. Evolution 38 (6), 1358–1370. doi:10.1111/j.1558-5646.1984.tb05657.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Wilson, T., Wu, X. Y., Juengel, J. L., Ross, I. K., Lumsden, J. M., Lord, E. A., et al. (2001). Highly Prolific Booroola Sheep Have a Mutation in the Intracellular Kinase Domain of Bone Morphogenetic Protein IB Receptor (ALK-6) that Is Expressed in Both Oocytes and Granulosa Cells. Biol. Reprod. 64 (4), 1225–1235. doi:10.1095/biolreprod64.4.1225

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, J. S., and Wang, J. Y. (2012). Polymorphism of Exon10 of FSHR Gene and its Relationship with Litter Size in Xiaomeishan Pigs. Scientia Agricultura Sinica 45 (13), 2728–2736.

Google Scholar

Yang, J., Lee, S. H., Goddard, M. E., and Visscher, P. M. (2011). GCTA: A Tool for Genome-wide Complex Trait Analysis. Am. J. Hum. Genet. 88 (1), 76–82. doi:10.1016/j.ajhg.2010.11.011

CrossRef Full Text | Google Scholar

Yuan, Y. Y., and Sun, L. M. (2014). Animal Genetic Resources in Yunnan. Kunming: Yunnan Scientific Publishing Press.

Google Scholar

Zhang, C., Dong, S. S., Xu, J. Y., He, W. M., and Yang, T. L. (2019). PopLDdecay: A Fast and Effective Tool for Linkage Disequilibrium Decay Analysis Based on Variant Call Format Files. Bioinformatics 35 (10), 1786–1788. doi:10.1093/bioinformatics/bty875

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, J., Nie, C., Li, X., Ning, Z., Chen, Y., Jia, Y., et al. (2020). Genome-Wide Population Genetic Analysis of Commercial, Indigenous, Game, and Wild Chickens Using 600K SNP Microarray Data. Front. Genet. 11, 543294. doi:10.3389/fgene.2020.543294

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhao, F., Deng, T., Shi, L., Wang, W., Zhang, Q., Du, L., et al. (2020). Genomic Scan for Selection Signature Reveals Fat Deposition in Chinese Indigenous Sheep with Extreme Tail Types. Animals (Basel) 10 (5). doi:10.3390/ani10050773

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhong, X. U., Hao, S., Zhe, Z., Qing-bo, Z., Olasege, B. S., Qiu-meng, L. I., et al. (2020). Genome-wide Detection of Selective Signatures in a Jinhua Pig Population. J. Integr. Agr. 19 (5), 1314–1322. doi:10.1016/S2095-3119(19)62833-9

CrossRef Full Text | Google Scholar

Zhu, C., Li, M., Qin, S., Zhao, F., and Fang, S. (2020). Detection of Copy Number Variation and Selection Signatures on the X Chromosome in Chinese Indigenous Sheep with Different Types of Tail. Asian-australas J. Anim. Sci. 33 (9), 1378–1386. doi:10.5713/ajas.18.0661

CrossRef Full Text | Google Scholar