Skip to main content

Breeding history and candidate genes responsible for black skin of Xichuan black-bone chicken

Abstract

Background

Domesticated chickens have a wide variety of phenotypes, in contrast with their wild progenitors. Unlike other chicken breeds, Xichuan black-bone chickens have blue-shelled eggs, and black meat, beaks, skin, bones, and legs. The breeding history and the economically important traits of this breed have not yet been explored at the genomic level. We therefore used whole genome resequencing to analyze the breeding history of the Xichuan black-bone chickens and to identify genes responsible for its unique phenotype.

Results

Principal component and population structure analysis showed that Xichuan black-bone chicken is in a distinct clade apart from eight other breeds. Linkage disequilibrium analysis showed that the selection intensity of Xichuan black-bone chickens is higher than for other chicken breeds. The estimated time of divergence between the Xichuan black-bone chickens and other breeds is 2.89 ka years ago. Fst analysis identified a selective sweep that contains genes related to melanogenesis. This region is probably associated with the black skin of the Xichuan black-bone chickens and may be the product of long-term artificial selection. A combined analysis of genomic and transcriptomic data suggests that the candidate gene related to the black-bone trait, EDN3, might interact with the upstream ncRNA LOC101747896 to generate black skin color during melanogenesis.

Conclusions

These findings help explain the unique genetic and phenotypic characteristics of Xichuan black-bone chickens, and provide basic research data for studying melanin deposition in animals.

Background

Domestication is a distinctive co-evolutionary, mutualistic relationship between humans and wild animals or plants that results in a range of genotypic and phenotypic impacts [1]. Animal domestication during the Neolithic period transformed the human lifestyle from hunting to farming, which enabled rapid changes in social organization and productivity [2]. Domesticated animals have spread to every region of the globe along with their human domesticators. The study of domestic animals contributes to our understanding of the evolution of animals under artificial selection, the influence of domestication or adaptive evolution on the animal genome, the genetic basis of evolution and phenotypic differentiation, and the optimization of animal breeding and diversity protection programs.

Chicken is one of the major sources of animal protein for humans. More recently, chickens have also become an important research model in fields such as physiology, disease, development, and aging, [3,4,5,6]. Chickens (Gallus Gallus domesticus) were the first domesticated bird species and were subjected for more than 8000 years to the combined effects of natural selection and human-driven artificial selection. Compared with their wild progenitors (red junglefowl, Gallus gallus), chickens present many characteristics associated with domestication that impact behavior, morphology, physiology, egg production, and skin color. A variety of studies have used whole-genome high-throughput DNA sequencing to reveal the genetic basis for traits acquired by natural and artificial selection in domesticated animals. For example, this approach has recently been applied to dogs [7], pigs [8, 9], chickens [10], sheep [11], rabbits [12], cattle [13, 14], and ducks [15, 16].

Studies suggest that northern China (alongside the Yellow River) is likely to have been one of several sites where chicken domestication occurred [17]. However, the history of breeding and the resulting phenotypic differentiation of indigenous chickens in this region have not been examined in sufficient detail. Both breeding history, and the genetic mechanisms underlying breed differentiation, have important theoretical implications for understanding the domestication, evolution, and phenotypic formation of chickens, and may also provide valuable insights for future breeding programs [18].

The Xichuan black-bone chicken (XBC), named for the Chinese prefecture of Xichuan, typically has five black parts (beak, skin, bones, legs, and meat) that distinguish it from other chicken breeds. Black-bone chickens are commonly believed to have medicinal properties and have been used as remedies to enhance the human immune system [19], prevent emaciation [20], treat diabetes [21], and cure conditions such as menstrual abnormalities and postpartum complications [22]. Xichuan black-bone chickens were primarily developed in a mountainous and inaccessible region in Xichuan County, China. Because transportation and trade were restricted, local resources were used for medical treatments, such as diet supplementation therapy. These incentives encouraged the selection of black skin in chickens. Xichuan black-bone chickens are highly prized, and became major income producers for farmers in the region. As demand increased, local farmers began to raise chickens at large scale and gradually developed this unique local chicken breed. However, the breeding history has not been examined in any detail, and the identity of the genetic factors responsible for the black body parts are unknown, particularly at the genomic level. In addition, the extent to which the Xichuan black-bone chickens has contributed germplasm to other breeds has not been studied.

In this study, we used whole genome resequencing to identify genomic features that illuminate the breeding history of the Xichuan black-bone chicken, and to correlate these features with the black characteristics of the breed. Genes with Xichuan black-bone chicken specific genomic variants were identified. Some genes in this class have undergone positive selection, and may be clues to the adaptive evolutionary history of the breed. By combining the genomic data with transcriptomic data, we were able to examine the genetic basis of adaptive evolution and breed differentiation more closely. Our results provide insight into the genetic factors underlying traits specific to Xichuan black-bone chickens. Furthermore, the data will provide a foundation for studying black coloration and modeling evolutionary selection mechanisms in this breed.

Results

Genetic variation in Xichuan black-bone chicken

To identify genetic variations, we resequenced the genomes of 5 Xichuan black-bone chickens. Using the Illumina sequencing platform, each animal yielded over 213 million clean reads, representing approximately 32 Gb per chicken (Supplementary Table S4). Q30 scores exceeded 94%, and the average GC content was 43.89%. The average sequencing depth was 28-fold per individual (a total of 160.90 Gb of high-quality paired-end sequence data). Around 98.45% of Xichuan black-bone chicken sequences were identical to those of red jungle fowl. The gene density per 200 kb and the number of SNPs, insertion/deletion polymorphisms (InDels), CNVs and SVs per 100 kb are shown in Fig. 1b. We identified 5,062,529 SNPs, including 247,054 homozygous SNPs, 830,606 InDels (372,903 insertions and 457,703 deletions), 1279 CNVs, and 11,433 SVs (Supplementary Fig. S1).

Fig. 1
figure1

Experimental design and variant statistics. a Geographical origins of the 9 chicken breeds used in this study. The map was created by the Adobe Illustrator (AI) 2019 software (http://adobe.e-bridge.com.cn/shop/index.html). b Summary of genomic resequencing data from 5 Xichuan black-bone chickens. The figure shows the distribution of SNPs, indels, and SVs. Chromosomes are represented in different colors in the outermost circle. The remaining circles, in order from the outside to the inside, are as follows: genes, SNPs, indels, and SVs. c Distribution of SNPs based on context. Different colors represent SNPs within various functional regions. One circle represents 1% of total SNPs. d Heat map depicting genetic relationships based on SNP data obtained for the 29 chickens that represented the 9 breeds. Colors represent pairwise genetic distances

To better understand the distribution of SNPs, we classified them according to their context into 13 categories (Fig. 1c). 39.23% of SNPs were found in intergenic regions. An even larger number were located in introns (47.23%). Smaller numbers of SNPs were within 5′ and 3′ UTRs (2.55% in aggregate), upstream or downstream of genes (3.14% in aggregate), splicing regions (0.01%), or were associated with ncRNA genes (6.88% in aggregate). Moreover, IS analysis indicated that chickens from the same habitat were more likely to have similar genetic distance and the clearest clusters (Fig. 1d).

Population structure and domestication

To examine genome-wide relationships and divergence between Xichuan black-bone chickens and other populations, we constructed a NJ tree with 1000 bootstrap replicates based on whole-genome polymorphic SNPs (Fig. 2a). The chickens clustered into three major branches that reflect geographic origins and breed. While YVC formed a distinct clade, the 5 RJFs fell into two clades (RJF-S1 and RJF-S2) that associated with other chickens (RJF-S1 and XGF, and RJF-S2 and TBC). This result aligns with their geographical distributions rather than their breed ascription, and is consistent with previous studies [23, 24], reinforcing the hypothesis that RJF was domesticated in several areas, as suggested by an analysis of mtDNA sequences [25, 26]. XBC clustered with WC, SK, DX, and LX, within the same branch. Results obtained using PCA were consistent with the NJ tree. The first three principal components accounted for 7.02% (PC1), 5.94% (PC2), and 5.27% (PC3) of total variability. The breeds formed distinct groupings except for RJF (Fig. 2b).

Fig. 2
figure2

Population genetics and LD decay. a NJ tree generated using polymorphisms detected in the 29 individual chickens. The scale bar represents the evolutionary distances measured by p-distance. Each of the 9 breeds has been assigned a distinct color. b Three-way PCA plots based on the 29 chickens. Symbol colors indicate breed (key on right). c Genetic structure of samples from 29 individuals for K groups using the ADMIXTURE program. K is the number of presumed ancestral groups, which was varied in the analysis from 2 to 10. The optimal K value was obtained with the least CV error value. d Decay of LD for XGF, YVC, RJF, TBC and XBC chickens measured by R2

Figure 2c shows a structure plot representing the 29 sampled chickens. At a low value of K (K = 2), XBC is clearly associated with a separate ancestor. When K was set to 4, 2 individuals from RJF, YVC and XGF (i.e., six individuals total) presented similar major components, while XBC and TBC appeared to cluster separately. XBC became almost distinguishable from all other breeds with increasing values of K. To estimate LD in the XGF, YVC, RJF, TBC, and XBC breeds, we calculated the squared correlations for two loci against the genome distance R2 between pairs of SNPs. Because the NJ tree analysis and PCA showed RJF-S2 to be relatively independent, we used RJF-S2 for the LD analysis. As shown in Fig. 2d, the most rapid attenuation was observed in RJF, followed by YVC, XGF, TBC, and XBC (Supplementary Table S5), indicating these breeds have relatively higher diversity and lower selection intensity. In contrast, the slowest attenuation was found in XBC. Thus, the XBC has experienced more domestication and selection intensity than have the other breeds.

Analysis of gene flow, time of divergence, and demographic history

We used TreeMix [27] to examine the topology of relationships and migration history among populations. We observed an early split between western (TBC), central (DX, LX, SK and XBC), and southern (XGF and YVC) populations (Fig. 3a & Fig. 3b). We detected a genetic contribution from RJF-S2 to XBC, and also found gene flow between XGF and LX, both of which are raised for cockfighting. Interestingly, a gene flow was observed from XBC to TBC (Fig. 3a & Supplementary Fig. S2).

Fig. 3
figure3

Gene flow analysis. a Maximum likelihood tree with 5 migration events. Migration events are shown as colored arrows, shaded according to their weight. Horizontal branch lengths are proportional to the amount of genetic drift that has occurred on each branch. The scale bar shows 10 times the average standard error of the entries in the sample covariance matrix. b Residual fit from the maximum likelihood tree in (a)

XBC differentiation time was analyzed by combining the available literature and fossil archeological data [17, 26, 28]. As shown in Fig. 4a (Supplementary Table S6), we found that XGF is much more closely related to YVC than to other chickens. The estimated time of divergence between XGF and YVC is 1.45 ka years ago, while the estimated divergence time of XBC is 2.89 ka years ago. Moreover, the results also demonstrated that the earliest differentiation occurred in RJF-S2 (5.78 ka years ago).

Fig. 4
figure4

Population genetics and demographic history. a Time of divergence between populations. The number at each node represents the time of divergence in thousands of years. b Demographic history of XBC and 4 other Chinese chicken breeds. Generation time (g) = 1 year and trans-version mutation rate (μ) = 1.91 × 10− 9 mutations per base pair per generation

As the unique genetic characteristics of XBC might be related to distinct divergence events, we conducted a PSMC for XBC and other Chinese populations to infer historical changes in effective population size (Ne). A tendency toward increased population size was detected in five of the populations 5 ka years ago (Fig. 4b), reaching a peak at the Last Glacial Maximum (20–26.5 ka), with a dramatic decrease following that peak [29]. The effective population size declined after the Last Glacial Maximum, with the increase in global temperature and development of human civilization [30].

Genome-wide selective sweep signals and functional analysis

In order to better detect genome-wide selection signals related to the unique black-skin trait, we divided the populations into black-skin and non-black-skin groups. Both the pi cut-off ratio (top 5%, pi ratio > 0.95 or < 0.05) and high Fst values (top 5%, Fst value> 0.17) were used as criteria for classifying selective sweeps. A total of 1469 candidate genes within these sweeps were associated with black-skin (Fig. 5a & Supplementary Table S7). Among them, one sweep located on chromosome 20 exhibited a high Fst value, indicating obvious genetic differentiation between black-skin and non-black-skin populations. In addition, a large difference in the pi value was also observed within this sweep between black and non-black groups. The pi value in the black-skin group was lower, suggesting that the sweep has been positively selected in the black-skin population (Fig. 5b). The candidate genes were then subjected to functional analysis. Top 30 of GO analysis revealed 10 GO terms enriched in biological processes, 8 terms in molecular functions, and 12 in cellular components (Supplementary Fig. S3). The KEGG pathway analysis results showed that candidate genes are mainly involved in the neuroactive ligand-receptor interaction, Jak-STAT signaling pathway, and so on (Fig. 6).

Fig. 5
figure5

Identification of genomic regions in Xichuan black-bone chickens with strong selective sweep signals. a Selective sweep signals are located to the left and right of the vertical dashed lines (representing Pi ratio values > 0.95 or < 0.05, respectively), and above the horizontal dashed line (representing an Fst value > 0.17). Regions selected for black skin are shown using blue points, while other skin colors are shown in green. The x-axis shows the pi ratio between black-skin and non-black-skin groups, and the y-axis shows the Fst values. b Genes with selective sweep signals in black and non-black skin. The shaded genomic regions contain selective signals for both skin types. pop1 and pop2 represent black skin non-black skin, respectively

Fig. 6
figure6

KEGG pathway enrichment analysis of candidate genes under selection in black skin chickens

Candidate genes for skin pigmentation

Skin color is an important domestication trait in chickens. We analyzed differences between black-bone chickens and other domesticated chickens to detect selection signatures associated with black skin. Genes associated with pigmentation that play important roles in the regulation of melanin deposition in mammals were identified in several genomic regions within selective sweeps. The candidate list includes solute carrier family 45 member 2 (SLC45A2), oxysterol binding protein like 2 (OSBPL2),

solute carrier family 24 member 2 (SLC24A2), PRELI domain containing 3B (SLMO2), ATP synthase, H+ transporting, mitochondrial F1 complex, epsilon subunit (ATP5e), cyclin dependent kinase inhibitor 2A (CDKN2A), GRAM domain containing 3 (GRAMD3), fibroblast growth factor 10 (FGF10) and Endothelin3 (EDN3).

To better understand the evolution of the selected genes, we hybridized Gushi chickens with Xichuan black-bone chickens to obtain F2 full-sibs with different skin colors. Yellow-skin and black-skin individuals with the same genetic background were then used for RNA-seq (Fig. 7). Four selected genes (SLC45A2, SLMO2, ATP5e, and EDN3) were differentially expressed between the black and yellow skin groups (Supplementary Table S8). A differentially expressed long noncoding RNA (TCONS-00054154) was also identified that is potentially related to black skin (Supplementary Table S9). TCONS-00054154 was used as a query using NCBI BLAST to identify related sequences, and was found to be identical to Gallus gallus uncharacterized LOC101747896, transcript variant X5, ncRNA. Both ncRNAs are located near EDN3, SLMO2 and ATP5e (Supplementary Table S10), suggesting that TCONS-00054154 might function in regulating these genes to produce different skin colors during melanogenesis. NJ analysis of EDN3 and TCONS-00054154 showed that they cluster to one branch in black-skin chickens (SK, DX, and XBC) (Fig. 8). This supports the hypothesis that these two genes might interact with the candidate genes and affect pigmentation in black skin.

Fig. 7
figure7

Histomorphological examination of Xichuan black-bone chickens with black (a) and yellow (b) skin. The picture was taken by Canon camera (5D Mark IV)

Fig. 8
figure8

Gene trees for LOC101747896 (a) and EDN3 (b), based on 29 chickens

Among the candidate genes, we found that fatty acid desaturase 6 (FADS6), Leptin receptor (LEPR), Lipoprotein lipase (LPL), melanin-concentrating hormone receptor 1 (MCHR1) and perilipin 2 (PLIN2) were associated with fatty acid-related pathways. Additionally, solute carrier organic anion transporter family member 1A2 (SLCO1A2), solute carrier organic anion transporter family member 1B1 (SLCO1B1), and solute carrier organic anion transporter family member 1C1 (SLCO1C1) belong to the organic ion transporter polypeptide (OATP) gene family, and may affect the transport of pigment in eggshells.

Identification of differentially expressed genes by qRT-PCR

As shown in Supplementary Fig. S4, gene expression levels were determined by RNA-seq sequencing and qRT-PCR. Our result suggested that the expression level for six genes obtained by qPCR were consistent with the high-throughput data, which showed the reliability of RNA-seq data.

Discussion

Xichuan black-bone chickens are rare in China and elsewhere, and have not been the focus for many studies. Due to the limited information available, we performed in-depth whole-genome sequencing of 5 black-bone chickens to explore the population structure, genetic diversity, and history of this breed. The data were analyzed along with sequences for other Chinese breeds that were downloaded from NCBI.

Sequence variations were distributed across different chromosome in numbers roughly proportional to chromosome size. Using the variations as the basis for our analysis, we found that Chinese chickens are divided into three large groups. Chickens from neighboring geographical locations and from similar altitudes are gathered together. Gene exchanges are likely to be correlated with geographical location, which is consistent with results from previous studies [14, 24, 31,32,33]. Population genomics analysis, including PCA, NJ tree, and structure, reveal that red jungle fowl can be separated into two branches, possibly because there have been different venues for domestication [17, 34, 35]. The LD decay analysis showed that XBC had the highest attenuation rate among the breeds tested, suggesting that XBC has been more domesticated and subjected to higher selection intensity than other chickens. This is in accordance with our expectations, since XBC has actually undergone a relative stronger artificial selection program than other breeds. Analysis of migration patterns among different varieties revealed that XBC evolved from the red jungle fowl, and another gene flow was observed between XBC and TBC. This is because breeds such as XBC were domesticated in the middle or lower reaches of the Yellow River, and then brought to many other regions by humans [36]. Some traits such as cockfighting obviously suffered positive selection, so it is reasonable that gene flow appeared between XGF and LX. Additionally, gene flow between chickens that were geographical neighbors has been proven. Obviously, the role of human activities cannot be ignored in the domestication history of native breeds. Human cultural exchange, trade, and migration have also promoted the domestication and spread of chickens.

Because the chicken is the most common domesticated animal in the world, their domestication has always been a topic of great interest in archeology and evolutionary biology. In 2014, Chinese and foreign scholars concluded that Northern China was a domestication center for chickens, dating back to the early Holocene, approximately 10,000 years ago [17]. However, this conclusion has been contested by other researchers [37, 38]. According to animal skeleton data and genomic information, the differentiation of chickens occurred 5.78 ka years ago. This finding is consistent with the most reliable archeological evidence [39]. Our results suggest that XBC diverged about 2.89 ka years ago, earlier than the appearance of cockfighting in China. This result also shows that human activities had a considerable influence on the domestication of chickens. Consistent with previous results, PSMC analysis revealed that the effective population of chickens fluctuated over time as the earth’s climate changed. The effective population peaked at the last glacial maximum (20–26 ka), and then decreased dramatically, as is the case for populations of other species [14, 31, 40].

Animals with distinct phenotypes in morphology, physiology, and behavior are excellent models for studying potential genetic changes and evolutionary mechanisms. Many traits in chickens, such as comb mass [41], aggressive behavior [42], reproductive performance [43], plateau adaptability [23], and vision [18] were influenced by selective sweeps. These studies do not include many other distinctive phenotypes found in various chicken breeds, such as black skin caused by pigmentation. Multiple genes affect skin and shank color. Many genes related to coat and skin color, such as SLC45A2, SLC24A2, SLMO2, ATP5e, CDKN2A, GRAMD3, FGF10, and EDN3, have been found in genomic regions under selective sweep [44,45,46,47,48,49].

Variations at the nucleotide level can affect gene transcription and the evolution of an organism. In order to better study the formation of melanin, we hybridized an exogenous chicken breed with XBC and collected F2 full-sibs with black and yellow dorsal skin. We used these individuals for high-throughput RNA-Seq and identified candidate genes for melanin deposition. We found that SLC45A2, SLMO2, ATP5e, and EDN3 were differentially expressed between the black and yellow dorsal skin types.

A large number of studies have established that the synthesis of melanin begins with the oxidation of tyrosine by tyrosinase, followed by a series of complex reactions. This process is regulated by a variety of signal molecules that affect the deposition and distribution of melanin in different animals. Membrane-associated transporter SLC45A2 is a member of the free vector family of SLC45A (solute carrier family 45) [50]. It is a transporter that regulates melanin synthesis and plays an important part in melanin production [51]. Mutations in SLC45A2 can affect the processing of tyrosinase, causing albinism in human eyes and skin [52]. SLC45A2 is differentially expressed in the skins of 1-month-old grey and white Lander geese [53]. In agreement with these studies, we found that SLC45A2 is differentially expressed in chicken skins of different color. The results suggest that the differential expression of SLC45A2 affects the deposition of melanin in black skin of the XBC, although the specific regulatory mechanism remained unknown.

Recent research has shown that the Fm phenotype (associated with the Fm region) is related to domestication and has been intensely artificially selected, with about 1.5-folder blacker level [48]. A previous study found that an inverted duplication affecting EDN3, a key gene in dermal melanin located on chromosome 20 in the Fm region, caused excessive accumulation of black pigment in chickens [54, 55]. We also detected an inversion at 10.97–11.52 Mb and a duplication at 10.74–12.55 Mb in EDN3, which overlaps the Fm region. EDN3 promotes the survival, proliferation, and differentiation of melanocytes [56, 57]. Differences in EDN3 copy number correspond with black and white pigmentation of chicken skin [58]. EDN3 is also the most likely candidate gene for coloration phenotypes in other domesticated animals, such as sheep [59], cats [60], and cattle [61]. The expression of lncRNA has been breed- and tissue- specific. In our study, no lnc-TMEM184C found by Hong et al. [62].

Green shanks are another typical characteristic of Xichuan black-bone chickens. GRAMD3 is highly expressed in human retinal pigment epithelial cells, and might be involved in the synthesis or transport of melanin. A mutation in the noncoding region of exon 4 of GRAMD3 may be linked to shank color in chickens [63]. Mutations in the flanking region of GRAMD3 cause abnormal expression of the gene in green shank [64]. Further research is needed to understand how GRAMD3 affects the shank color.

Our results suggest that FADS6, LEPR, LPL, MCHR1, PLIN2, and PYROXD1 might be associated with fatty acid-related pathways. Fatty acid desaturase 6 (FADS6) is related to fatty acid composition (FAC) in Hanwoo Cattle [65]. LPL produces fatty acids and monoglycerides for tissue use or storage. LEPR is an important mediator for leptin function. Mutations in the LEPR gene are associated with obesity in humans and fat deposition in animals such as cattle and pigs [66,67,68]. Most notably, MCHR1, when it is expressed in the hypothalamic ventromedial nucleus and paraventricular nucleus in the brain, is related to feeding behavior [69, 70]. When MCHR1 is distributed in the peripheral tissues of mice, it affects the synthesis of fat and is related to obesity [71, 72]. There are highly specific autoantibodies against MCHR1 in vitiligo patients that play an important role in the pathogenesis of vitiligo autoimmunity [73]. When MCHR1 is highly expressed in the dorsal skin of Misgurnus anguillicaudatus (pond loach), it helps the formation of dark spots. We suggest that MCHR1 might play an important role in the meat quality and pigmentation in XBC.

SLCO1A2, SLCO1B1 and SLCO1C1 function in pigment transportation in poultry eggshells. SLCO2B1 is a member of the OATP gene family. The protein encoded by SLCO2B1 might be related to biliverdin transport [74]. Brown shell traits are closely related to protoporphyrin transportation in the eggshell gland, a process that is regulated by OATP gene family members such as SLCO1A2 and SLCO1C1 [75].

Conclusions

In summary, the evolutionary history, molecular phylogeny, and selection evidence of several chicken breeds were explored in this study. The results provide new insight into Xichuan black-bone chickens and their unique qualities.

Methods

Sample collection and sequencing

Whole blood samples from five 30-week-old female Xichuan black-bone chickens were collected from Xichuan County (32°56′11″N, 111°24′46″E) in Henan Province, which come from different conservation populations and have no genetic relationship. The pentobarbital was used for euthanization by intraperitoneal injection at a dose of 40 mg/kg of body weight according our previous study [76]. We made all efforts to decrease animal suffering. Genomic DNA was isolated using a standard phenol/chloroform extraction method. Briefly, library construction involved the following steps: DNA shearing, purification, end-repairing, ligating adaptors, size selection, and DNA amplification. The libraries were sequenced using the Illumina HiSeq 4000 platform (PE300). Sequencing and base calling were performed using the manufacturer’s protocols to generate primary read data.

Sequence data for 24 individuals from other chicken breeds was obtained from published dataset (Supplementary Table S1) [23, 44]. The downloaded data were generated using Illumina HiSeq 2000 technology with 18X coverage per sample (Supplementary Table S2). Combining our data with the downloaded data, there were 29 chickens in total, including 5 XBC, 5 Tibetan chickens (TBC) originating from Tibet, 5 Xishuangbanna game fowl (XGF) originating from Dai Autonomous Prefecture of Xishuangbanna of Yunnan Province (raised for cockfighting), 5 Yunnan Village chickens (YVC, raised for food), 5 Red Junglefowl (RJF) originating from Yunnan Province, and one Dongxiang blue-eggshell chicken (DX) and one Silkie chicken (SK), both originating from Jiangxi Province, one Luxi fighting chicken (LX) originating from Shandong Province, and one Wenchang chicken (WC) originating from Hainan Province. The data were used to analyze the evolutionary history of Xichuan black-bone chickens and other populations (Fig. 1a). Breed origins were as described by the Animal Genetic Resources in China for Poultry.

Quality control processing and variant calling

To ensure high-quality data, the raw sequencing data for all breeds were filtered to remove adaptors, trim low-quality and unidentified nucleotides (Ns) at the 5′ and 3′ ends, trim nucleotides having an average quality score less than 30 within a window of 4 nucleotides, and removing reads containing more than 10% Ns. After these steps, reads with lengths less than 30 bp were removed. High-quality reads were aligned to the chicken reference genome (version: galGal5) using BWA-MEM (version: 0.7.12) [77]. Sequence variants were identified using a GATK Best Practices workflow to identify single nucleotide polymorphisms (SNPs) and insertion/deletion sequences (InDels) [78,79,80]. Detection of Structural variations (SVs) and copy number variants (CNVs) was performed using BreakDancer [81] and CNVnator [82]. Finally, gene-based annotation of the variants was accomplished using SnpEff [83], for which the corresponding chicken gene annotation file was downloaded from Ensembl. The distribution of InDels, SNPs, CNVs, and SVs was visualized using Circos [84] and arranged in chromosomal order. Whole-genome identity score (IS) were calculated with SNP frequencies [85]. To compare the Xichuan black-bone chicken with the other breeds, we adjusted the sequencing data for the Xichuan black-bone chicken so that the average sequencing depth was similar across breeds.

Population structure and linkage disequilibrium analysis

To reduce noise due to linkage disequilibrium (LD), SNPs with high pair-wise R2 values (R2 > 0.2) were pruned from the dataset using PLINK (v1.90, arguments: --indep-pairwise 50 5 0.2) [86]. The neighbor-joining tree (NJ) was constructed with Nei’s genetic distances using the phylogeny program MEGA v7.0 [87, 88]. Principal component analysis (PCA) was performed using Genome-wide Complex Trait Analysis (GCTA) [89] and visualized with R package (ggplot2). The population structure was investigated using the supervised ADMIXTURE program (Version 1.3) [90] using the maximum likelihood method. We predefined K, the number of genetic clusters (ranging from K = 2 to K = 10). LD for five breeds was calculated on the basis of the correlation coefficient R2 statistics of two loci using PopLDdecay [91].

Gene flow, estimation of divergence time, and demographic history

TreeMix [27] was used to construct a population phylogenetic tree for gene flow analysis of the 9 breeds. The divergence time between XBC and other groups was estimated using Beast2 and TreeAnnotator [92]. The alignment result was converted into Nexus file format and used as the input for Beast2. Parameters were set as follows: Partition, Unlink; Site Model, HYK; Gamma Category Count, 4; Clock model, Strict; Priors tree prior, Yule Model; MCMC, 10,000,000; Tracelog, Screenlog and treelog.t:tree, default. A tree file was one of the outputs generated by Beast2. The tree file was used as input to TreeAnnotator. We used a pairwise sequentially Markovian coalescent (PSMC) model [93] to detect changes in the effective ancestral population sizes of the XBC population at a high level of detail. We set g = 1 and a mutation rate of 1.91 × 10− 9 per generation to estimate the distribution over time.

Selective sweep analysis

To define candidate regions that had undergone directional selection in the XBC, we calculated the population differentiation index (Fst) between black-skin chickens (XBC, DX, SK) and non-black-skin chickens (TBC, XGF, RJF, YVC, LX, WC), as described by Weir & Cockerham [94]. We calculated the average Fst value in 10 kb sliding windows with a 5 kb sliding step between all black and all non-black populations using VCFtools (v0.1.13) [95]. The fixation index (Fst) measures the genetic differentiation between two populations. For each SNP, Fst is defined as:

\( {F}_{ST}=1-\frac{p_1\left(1-{p}_1\right)+{p}_2\left(1-{p}_2\right)}{\left({p}_1+{p}_2\right)\left[1-\left({p}_1+{p}_2\right)/2\right]} \), where p1 and p2 are the frequencies of an allele in all black and non-black populations, respectively.

Nucleotide diversity (pi) was also estimated for each sliding window across all sampled populations. The pi primarily reflects the nucleotide difference. The top 5% of Fst values were retained as candidate sweeps. The sweeps also exhibited large differences in pi values between black and non-black populations. Finally, we performed a functional enrichment analysis of the candidate genes, using Gene Ontology (GO) categories and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways. The analysis was conducted with the Database for kobas3.0 (http://kobas.cbi.pku.edu.cn/kobas3).

RNA-Seq and data processing

To determine whether novel allelic variants located in the top 5% of Fst regions of the genome affect gene expression, we utilized RNA-Seq to construct six cDNA libraries from black (BS) and yellow (YS) dorsal skin of a Xichuan black-bone chickens hybrid population. The hybrid F2 population was obtained by crossing the Gushi chickens with Xichuan black-bone chickens. All the chickens were raised in cages under the same environment with ad libitum conditions. Full-sibs with different skin colors of 1-day-old female chickens were collected for RNA-seq. RNA-seq libraries with insert sizes of approximately 300 bp were prepared using the Illumina standard RNA-seq library preparation pipeline and sequenced on the Illumina HiSeq 4000 platform.

Adaptor sequences and low-quality sequence reads were removed from the data sets. Clean reads were then mapped to the reference genome sequence (https://www.ncbi.nlm.nih.gov/genome/?term=chicken). Differentially expressed genes were detected with edgeR (v3.6) [96, 97] using a threshold of |log2 (fold change) | ≥1 and adjusted p < 0.05.

Quantitative RT-PCR and data analysis

To verify the accuracy of the data obtained from RNA-seq, the results were confirmed by qRT-PCR. Six gene-specific qRT-PCR primers were designed using the Primers-BLAST online program at NCBI and synthesized by Shenggong Biological Engineering (Wuhan, China) Limited by Share Ltd. Primer sequences are shown in Supplementary Table S3. All qPCR was conducted using a LightCycler® 96 Real-Time PCR system (Roche Applied Science). Each sample was assayed three times. Chicken β-actin gene was used as a reference gene for the qPCR assays. Relative expression was calculated using the 2-ΔΔCt method [98].

Availability of data and materials

Data for the 5 Xichuan black-bone chickens used in the whole-genome resequencing analysis, and for the 6 Xichuan black-bone chickens used in the RNA-Seq analysis are accessible at NCBI under BioProject accession numbers PRJNA642410 and PRJNA418694, respectively. Illumina paired-end sequences for the other 24 chickens used in this study were downloaded from NCBI with SRA accession numbers SRP034930 and SRP040477. Chicken reference genome sequence was downloaded from NCBI (https://www.ncbi.nlm.nih.gov/genome/?term=chicken). chicken gene annotation file was downloaded from Ensemble (ftp://ftp.ensembl.org/pub/release-100/gff3/gallus_gallus).

Abbreviations

XBC:

Xichuan black-bone chicken

TBC:

Tibetan chickens

WC:

Wenchang chicken

SK:

Silkie chicken

DX:

Dongxiang blue-eggshell chicken

LX:

Luxi fighting chicken

BS:

Black skin

YS:

Yellow dorsal skin

InDels:

Insertion/deletion polymorphisms

SNPs:

Single nucleotide polymorphisms

CNVs:

Copy numbervariations

SVs:

Structural variations

IS:

Identity score

LD:

Linkage disequilibriu

NJ:

Neighbor-joining

PCA:

Principal component analysis

GCTA:

Genome-wide Complex Trait Analysis

GO:

Gene Ontology

KEGG:

Kyoto Encyclopedia of Genes and Genomes

References

  1. 1.

    Zeder MA. Core questions in domestication research. Proc Natl Acad Sci U S A. 2015;112(11):3191–8.

    CAS  PubMed  PubMed Central  Google Scholar 

  2. 2.

    Jing L, Yaping Z. Advances in research of the origin and domestication of domestic animals. Biodivers Sci. 2009;17(4):319–29.

    Google Scholar 

  3. 3.

    Zhang M, Yan FB, Li F, Jiang KR, Li DH, Han RL, Li ZJ, Jiang RR, Liu XJ, Kang XT. Genome-wide DNA methylation profiles reveal novel candidate genes associated with meat quality at different age stages in hens. Sci Rep. 2017;7:45564.

    CAS  PubMed  PubMed Central  Google Scholar 

  4. 4.

    Zhang Y, Li D, Han R, Wang Y, Li G, Liu X, Tian Y, Kang X, Li Z, Yun Z. Transcriptome analysis of the pectoral muscles of local chickens and commercial broilers using Ribo-Zero ribonucleic acid sequencing. PLoS One. 2017;12(9):e0184115.

    PubMed  PubMed Central  Google Scholar 

  5. 5.

    Deist MS, Gallardo RA, Bunn DA, Kelly TR, Dekkers JCM, Zhou H, Lamont SJ. Novel analysis of the Harderian gland transcriptome response to Newcastle disease virus in two inbred chicken lines. Sci Rep. 2018;8(1):6558.

    PubMed  PubMed Central  Google Scholar 

  6. 6.

    Zhang J, Kaiser MG, Deist MS, Gallardo RA, Bunn DA, Kelly TR, Dekkers JCM, Zhou H, Lamont SJ. Transcriptome analysis in spleen reveals differential regulation of response to Newcastle disease virus in two chicken lines. Sci Rep. 2018;8(1):1278.

    PubMed  PubMed Central  Google Scholar 

  7. 7.

    Axelsson E, Ratnakumar A, Arendt ML, Maqbool K, Webster MT, Perloski M, Liberg O, Arnemo JM, Hedhammar A, Lindblad-Toh K. The genomic signature of dog domestication reveals adaptation to a starch-rich diet. Nature. 2013;495(7441):360–364,a363.

    CAS  PubMed  Google Scholar 

  8. 8.

    Rubin CJ, Megens HJ, Barrio AM, Maqbool K, Sayyab S, Schwochow D, Wang C, Carlborg O, Jern P, Jorgensen CB. Strong signatures of selection in the domestic pig genome. Proc Natl Acad Sci U S A. 2012;109(48):19529–36.

    CAS  PubMed  PubMed Central  Google Scholar 

  9. 9.

    Li M, Tian S, Jin L, Zhou G, Li Y, Zhang Y, Wang T, Yeung CKL, Chen L, Ma J. Genomic analyses identify distinct patterns of selection in domesticated pigs and Tibetan wild boars. Nat Genet. 2013;45(12):1431–8.

    CAS  PubMed  Google Scholar 

  10. 10.

    Keeling L, Andersson L, Schütz KE, Kerje S, Fredriksson R, Carlborg R, Cornwallis CK, Pizzari T, Jensen P. Chicken genomics: feather-pecking and victim pigmentation. Nature. 2004;431(7009):645–6.

    CAS  PubMed  Google Scholar 

  11. 11.

    Li X, Su R, Wan W, Zhang W, Jiang H, Qiao X, Fan Y, Zhang Y, Wang R, Liu Z. Identification of selection signals by large-scale whole-genome resequencing of cashmere goats. Sci Rep. 2017;7(1):15142.

    PubMed  PubMed Central  Google Scholar 

  12. 12.

    Carneiro M, Rubin C, Palma FD, Albert FW, Alfoldi J, Barrio AM, Pielberg G, Rafati N, Sayyab S, Turnermaier J. Rabbit genome analysis reveals a polygenic basis for phenotypic change during domestication. Science. 2014;345(6200):1074–9.

    CAS  PubMed  PubMed Central  Google Scholar 

  13. 13.

    Jiang J, Gao Y, Hou Y, Li W, Zhang S, Zhang Q, Sun D. Whole-genome resequencing of holstein bulls for indel discovery and identification of genes associated with milk composition traits in dairy cattle. PLoS One. 11(12):e0168946.

  14. 14.

    Lan D, Xiong X, Mipam TD, Fu C, Li Q, Ai Y, Hou D, Chai Z, Zhong J, Li J. Genetic diversity, molecular phylogeny and selection evidence of Jinchuan Yak revealed by whole-genome resequencing. G3. 2018;8(3):945–52.

    CAS  PubMed  Google Scholar 

  15. 15.

    Zhang Z, Jia Y, Pedro A, Mank JE, Marcel V, Wang Q, Jiang Z, Yu C, Kai Z, Hou S. Whole-genome resequencing reveals signatures of selection and timing of duck domestication. Gigascience. 2018;(4):4.

  16. 16.

    Zhou Z, Ming L, Hong C, Fan W, Yuan Z, Qiang G, Xu Y, Guo Z, Zhang Y, Jian H. An intercross population study reveals genes associated with body size and plumage color in ducks. Nat Commun. 2018;9(1):2648.

    PubMed  PubMed Central  Google Scholar 

  17. 17.

    Xiang H, Gao J, Yu B, Zhou H, Zhao X. Early Holocene chicken domestication in northern China. Proc Natl Acad Sci U S A. 2014;111(49):17564–9.

    CAS  PubMed  PubMed Central  Google Scholar 

  18. 18.

    Wang MS, Zhang RW, Su LY, Li Y, Peng MS, Liu HQ, Zeng L, Irwin DM, Du JL, Yao YG. Positive selection rather than relaxation of functional constraint drives the evolution of vision during chicken domestication. Cell Res. 2018;26(5):556.

    Google Scholar 

  19. 19.

    Bai F, Li YH, Huang FJ, Wu ZZ. Study on nutritional composition and functional properties of peptides from black-bone silky fowl (Gallus gallus domesticus Brisson). Food Res Dev (in China). 2019;40(002):17–21.

    Google Scholar 

  20. 20.

    Tu YG, Sun YZ, Tian YG, Xie MY, Chen J. Physicochemical characterisation and antioxidant activity of melanin from the muscles of Taihe Black-bone silky fowl (Gallus gallus domesticus Brisson). Food Chem. 2009;114(4):1345–50.

    CAS  Google Scholar 

  21. 21.

    Zhang FW. 100 cases of chronic prostatitis treated with Wuji Baifeng pill. J Anhui Univ Tradit Chin Med (in China). 2002;3(6):23.

    Google Scholar 

  22. 22.

    Tian Y, Xie M, Wang W, Wu H, Fu Z, Lin L. Determination of carnosine in black-bone silky fowl (Gallus gallus domesticus Brisson) and common chicken by HPLC. Eur Food Res Technol. 2007;226(1–2):311–4.

    CAS  Google Scholar 

  23. 23.

    Wang MS, Li Y, Peng MS, Zhong L, Wang ZJ, Li QY, Tu XL, Dong Y, Zhu CL, Wang L. Genomic analyses reveal potential independent adaptation to high altitude in Tibetan chickens. Mol Biol Evol. 2015;32(7):1880–9.

    CAS  PubMed  Google Scholar 

  24. 24.

    Li D, Che T, Chen B, Tian S, Zhou X, Zhang G, Li M, Gaur U, Li Y, Luo M. Genomic data for 78 chickens from 14 populations. GigaScience. 2017;6(6):1–5.

    PubMed  PubMed Central  Google Scholar 

  25. 25.

    Liu YP, Wu GS, Yao YG, Miao YW, Zhang YP. Multiple maternal origins of chickens: out of the Asian jungles. Mol Phylogenet Evol. 2006;38(1):12–9.

    CAS  PubMed  Google Scholar 

  26. 26.

    Miao YW, Peng MS, Wu GS, Ouyang YN, Yang ZY, Yu N, Liang JP, Pianchou G, Bejapereira A, Mitra B. Chicken domestication: an updated perspective based on mitochondrial genomes. Heredity. 2013;3(110):277–82.

    Google Scholar 

  27. 27.

    Pickrell JK, Pritchard JK, Tang H. Inference of population splits and mixtures from genome-wide allele frequency data. PLoS Genet. 2012;8(11):e1002967.

    CAS  PubMed  PubMed Central  Google Scholar 

  28. 28.

    Nam K, Mugal CF, Nabholz B, Schielzeth H, Wolf JBW, Backstrom N, Kunstner A, Balakrishnan CN, Heger A, Ponting CP. Molecular evolution of genes in avian genomes. Genome Biol. 2010;11(6):1–17.

    Google Scholar 

  29. 29.

    Clark PU, Dyke AS, Shakun JD, Carlson AE, Clark J, Wohlfarth B, Mitrovica JX, Hostetler SW, Mccabe AM. The last glacial maximum. Science. 2009;325(5941):710–4.

    CAS  Google Scholar 

  30. 30.

    Tallavaara M, Luoto M, Korhonen N, Järvinen H, Seppä H. Human population dynamics in Europe over the last glacial maximum. Proc Natl Acad Sci U S A. 2015;112(27):8232–7.

    CAS  PubMed  PubMed Central  Google Scholar 

  31. 31.

    Li L, Shi X, Zheng F, Wu D, Li A-A, Sun F-Y, Li C-C, Wu J-C, Li T. Transcriptome analysis of Dlm mutants reveals the potential formation mechanism of lesion mimic in wheat. Eur J Plant Pathol. 2016;146(4):987–97.

    CAS  Google Scholar 

  32. 32.

    Zhao P, Yu Y, Feng W, Du H, Yu J, Kang H, Zheng X, Wang Z, Liu GE, Ernst CW. Evidence of evolutionary history and selective sweeps in the genome of Meishan pig reveals its genetic and phenotypic characterization. Gigascience. 2018;5:5.

    Google Scholar 

  33. 33.

    Freedman AH, Lohmueller KE, Wayne RK. Evolutionary history, selective sweeps, and deleterious variation in the dog. Annu Rev Ecol Evol Syst. 2016;47(1):73–96.

    Google Scholar 

  34. 34.

    Tixier-Boichard M, Bed’Hom B, Rognon X. Chicken domestication: from archeology to genomics. C R Biol. 2011;334(3):0–204.

    Google Scholar 

  35. 35.

    Eda M, Lu P, Kikuchi H, Li Z, Li F, Yuan J. Reevaluation of early Holocene chicken domestication in northern China. J Archaeol Sci. 2016;67:25–31.

    Google Scholar 

  36. 36.

    Wang Z, Xiang H, Yuan J, Luo Y, Zhao X. Exploring the origin of domesticated pigs in the Yellow River area using information from ancient DNA. Sci Bull. 2012;57(12):1011.

    Google Scholar 

  37. 37.

    Peng MS, Shi NN, Yao YG, Zhang YP. Caveats about interpretation of ancient chicken mtDNAs from northern China. Proc Natl Acad Sci U S A. 2015;112(16):1970–1.

    Google Scholar 

  38. 38.

    Peters J, Lebrasseur O, Best J, Miller H, Fothergill T, Dobney K, Thomas RM, Maltby M, Sykes N, Hanotte O. Questioning new answers regarding Holocene chicken domestication in China. Proc Natl Acad Sci U S A. 2015;112(19):201503579.

    Google Scholar 

  39. 39.

    Yuan J. Zooarchaeological study on the domestic animals in ancient China. Quaternary Sci. 2010;30:298–306.

    Google Scholar 

  40. 40.

    Zhou X, Wang B, Pan Q, Zhang J, Kumar S, Sun X, Liu Z, Pan H, Lin Y, Liu G. Whole-genome sequencing of the snub-nosed monkey provides insights into folivory and evolutionary history. Nat Genet. 2014;46(12):1303–10.

    CAS  PubMed  Google Scholar 

  41. 41.

    Johnsson M, Gering E, Willis PM, Lopez S, Van Dorp L, Hellenthal G, Henriksen R, Friberg U, Wright D. Feralisation targets different genomic loci to domestication in the chicken. Nat Commun. 2016;7(1):12950.

    CAS  PubMed  PubMed Central  Google Scholar 

  42. 42.

    Guo X, Fang Q, Ma C, Zhou B, Wan Y, Jiang R. Whole-genome resequencing of Xishuangbanna fighting chicken to identify signatures of selection. Genet Sel Evol. 2016;48(1):62.

    PubMed  PubMed Central  Google Scholar 

  43. 43.

    Moges F, Mellesse A, Dessie T. Assessment of village chicken production system and evaluation of the productive and reproductive performance of local chicken ecotype in bure district, north West Ethiopia. Afr J Agric Res. 2010;5(13):1739–48.

    Google Scholar 

  44. 44.

    Yi G, Qu L, Liu J, Yan Y, Xu G, Yang N. Genome-wide patterns of copy number variation in the diversified chicken genomes using next-generation sequencing. BMC Genomics. 2014;15(1):962.

    PubMed  PubMed Central  Google Scholar 

  45. 45.

    Gao X, Qiao YE, Xing LI, Wang J, Guo Y, Jiguo XU, Nie Q. Identification of mutations and related frequency distributions involving appearance traits in Xuefeng black-boned chicken. China Poult. 2016;22(2):5–9.

    Google Scholar 

  46. 46.

    Kim JDF, Velie BD, Jeanette A, Rachel A, Hamilton NA, Leif A, Meadows JRS, Gabriella L. A potential regulatory region near the EDN3 gene may control both harness racing performance and coat color variation in horses. Phys Rep. 2018;6(10):e13700.

    Google Scholar 

  47. 47.

    Liu Y, Xue L, Chang L, Bai Y, Lixin LI, Jingwei LI, Xiuju YU, Xiaoyan HE, Wang H. Effects of EDN3 on the expression of melanin-related genes in mouse melanocytes. China Anim Husbandry Vet Med. 2018;45(8):581–9.

    Google Scholar 

  48. 48.

    Budhi DA, Yohei T, Sri S, Arifin ZMS, Toyoko A, Yoko S, Petr H. The origin and evolution of fibromelanosis in domesticated chickens: genomic comparison of Indonesian Cemani and Chinese Silkie breeds. PLoS One. 2017;12(4):e0173147.

    Google Scholar 

  49. 49.

    Darwish HYA, Zhang Y, Cui K, Yang Z, Han D, Dong X, Mao H, Deng W, Deng X. Molecular cloning and characterization of the endothelin 3 gene in black bone sheep. J Anim Sci Biotechnol. 2018;9(1):57.

    PubMed  PubMed Central  Google Scholar 

  50. 50.

    Vitavska O, Wieczorek H. The SLC45 gene family of putative sugar transporters. Mol Asp Med. 2013;34(2–3):655–60.

    CAS  Google Scholar 

  51. 51.

    Gunnarsson U, Hellstrom AR, Tixierboichard M, Minvielle F, Bedhom B, Ito S, Jensen P, Rattink A, Vereijken A, Andersson L. Mutations in SLC45A2 cause plumage color variation in chicken and Japanese quail. Genetics. 2007;175(2):867–77.

    CAS  PubMed  PubMed Central  Google Scholar 

  52. 52.

    Cullinane AR, Vilboux T, Obrien KJ, Curry JA, Maynard DM, Carlsondonohoe H, Ciccone C, Markello TC, Gunayaygun M, Huizing M. Homozygosity mapping and whole-exome sequencing to detect SLC45A2 and G6PC3 mutations in a single patient with oculocutaneous albinism and neutropenia. J Investig Dermatol. 2011;131(10):2017–25.

    CAS  PubMed  Google Scholar 

  53. 53.

    Wei A, Yang X, Lian S, Li W. Genetic analyses of Chinese patients with digenic oculocutaneous albinism. Chin Med J. 2013;126(2):226.

    CAS  PubMed  Google Scholar 

  54. 54.

    Dorshorst B, Molin A, Rubin C, Johansson AM, Stromstedt L, Pham M, Chen C, Hallbook F, Ashwell CM, Andersson L. A complex genomic rearrangement involving the Endothelin 3 locus causes dermal hyperpigmentation in the chicken. PLoS Genet. 2011;7(12):e1002412.

    CAS  PubMed  PubMed Central  Google Scholar 

  55. 55.

    Sohn J, Nam K, Hong H, Kim J, Lim D, Lee K, Do YJ, Cho CY, Kim N, Chai H. Whole genome and transcriptome maps of the entirely black native Korean chicken breed Yeonsan Ogye. GigaScience. 2018;7(7):1–14.

    CAS  Google Scholar 

  56. 56.

    Hutt FB. The genetics of the fowl. J Genet. 1930;22(1):109–27.

    Google Scholar 

  57. 57.

    Shinomiya A, Kayashima Y, Kinoshita K, Mizutani M, Namikawa T, Matsuda Y, Akiyama T. Gene duplication of endothelin 3 is closely correlated with the hyperpigmentation of the internal organs (Fibromelanosis) in silky chickens. Genetics. 2012;190(2):627–38.

    CAS  PubMed  PubMed Central  Google Scholar 

  58. 58.

    Han R, Yang P, Tian Y, Wang D, Zhang Z, Wang L, Li Z, Jiang R, Kang X. Identification and functional characterization of copy number variations in diverse chicken breeds. BMC Genomics. 2014;15(1):934.

    PubMed  PubMed Central  Google Scholar 

  59. 59.

    Yanan LI, Zhao BL, Wang HD, Chen TZ, Liu Y, Chang LC, Dong CS. Effect of EDN3 on of sheep skin melanocytes with different coat colors in vitro. Sci Agric Sin. 2017;50(6):1139–46.

    Google Scholar 

  60. 60.

    Kaelin CB, Xu X, Hong LZ, David VA, Mcgowan KA, Schmidtkuntzel A, Roelke ME, Pino J, Pontius J, Cooper GM. Specifying and sustaining pigmentation patterns in domestic and wild cats. Science. 2012;337(6101):1536–41.

    CAS  PubMed  PubMed Central  Google Scholar 

  61. 61.

    Qanbari S, Pausch H, Jansen S, Somel M, Strom TM, Fries R, Nielsen R, Simianer H. Classic selective sweeps revealed by massive sequencing in cattle. PLoS Genet. 2014;10(2):e1004148.

    PubMed  PubMed Central  Google Scholar 

  62. 62.

    Hong H, Chai HH, Nam K, Lim D, Nam JW. Non-coding transcriptome maps across twenty tissues of the Korean black chicken, Yeonsan Ogye. Int J Mol Sci. 2018;19(8):2359.

    PubMed Central  Google Scholar 

  63. 63.

    Clement P. Polymorphisms and associations of VCAN and GRAMD3 gene with black shank pigmentation in chicken. Master Degree Diss ( in China). 2015;9:24.

    Google Scholar 

  64. 64.

    Xu J, Lin S, Gao X, Nie Q, Zhang X. Mapping of Id locus for dermal shank melanin in a Chinese indigenous chicken breed. J Genet. 2017;96(6):977–83.

    CAS  PubMed  Google Scholar 

  65. 65.

    Srikanth K, Lee E, Kwan A, Lim Y, Chung H. Genetic variations in the bovine fatty acid desaturase 6 (FADS6) are associated with fatty acid composition in Hanwoo cattle. 2016;8(2):41–9.

  66. 66.

    Komisarek J. Impact of LEP and LEPR gene polymorphisms on functional traits in polish Holstein-Friesian cattle. Anim Sci Paper Rep. 2010;28(2):133–41.

    CAS  Google Scholar 

  67. 67.

    Mindekova S, Trakovicka A, Trandík J, Jr JB, Laszlo Z. Correlation of pig LEPR and H-FABP parental genotypes with fat content of meat in offsprings. Magyar Allatorvosok Lapja. 2010;132(1):14–21.

    CAS  Google Scholar 

  68. 68.

    Leońska-Duniec A, Jastrzbski Z, Jadewska A, Krzysztof F, Ciszczyk P. Leptin and leptin receptor genes are associated with obesity-related traits changes in response to aerobic training program. J Strength Cond Res. 2018;32(1):1036.

    PubMed  Google Scholar 

  69. 69.

    Lembo PMC, Grazzini E, Cao J, Hubatsch DA, Walker P. The receptor for the orexigenic peptide melanin-concentrating hormone is a G-protein-coupled receptor. Nat Cell Biol. 1999;1(5):267–71.

    CAS  PubMed  Google Scholar 

  70. 70.

    Hervieu GJ, Cluderay JE, Harrison D, Meakin J, Leslie R. The distribution of the mRNA and protein products of the melanin-concentrating hormone (MCH) receptor gene, slc-1, in the central nervous system of the rat. Eur J Neurosci. 2000;12(4):1194–216.

    CAS  PubMed  Google Scholar 

  71. 71.

    Bradley RL, Kokkotou EG, Maratos-Flier E, Cheatham B. Melanin-concentrating hormone regulates leptin synthesis and secretion in rat adipocytes. Diabetes. 2000;49(7):1073–7.

    CAS  PubMed  Google Scholar 

  72. 72.

    Bradley RL, Mansfield JPR, Eleftheria MF, Bentley C. Melanin-concentrating hormone activates signaling pathways in 3T3-L1 adipocytes. Am J Physiol Endocrinol Metab. 2002;283(3):584–92.

    Google Scholar 

  73. 73.

    Kemp EH, Waterman EA, Hawes BE, O'Neill K, Gottumukkala RV, Gawkrodger DJ, Weetman AP, Watson PF. The melanin-concentrating hormone receptor 1, a novel target of autoantibody responses in vitiligo. J Clin Investig. 2002;109(7):923–30.

    CAS  PubMed  Google Scholar 

  74. 74.

    Shitara Y, Takeuchi K, Nagamatsu Y, Wada S, Sugiyama Y, Horie T. Long-lasting inhibitory effects of cyclosporin a, but not tacrolimus, on OATP1B1- and OATP1B3-mediated uptake. Drug Metab Pharmacokinet. 2012;27(4):368–78.

    CAS  PubMed  Google Scholar 

  75. 75.

    Zheng C, Li Z, Yang N, Ning Z. Quantitative expression of candidate genes affecting eggshell color. Anim Sci J. 2014;85(5):506–10.

    CAS  PubMed  Google Scholar 

  76. 76.

    Zhang M, Li F, Ma X, Li W, Jiang R, Han R, Li G, Wang Y, Li Z, Tian Y. Identification of differentially expressed genes and pathways between intramuscular and abdominal fat-derived preadipocyte differentiation of chickens in vitro. BMC Genomics. 2019;20(1):743.

    PubMed  PubMed Central  Google Scholar 

  77. 77.

    Li H. Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM. arXiv: Genomics. 2013.

  78. 78.

    Mckenna A, Hanna M, Banks E, Sivachenko A, Cibulskis K, Kernytsky A, Garimella K, Altshuler D, Gabriel S, Daly M. The genome analysis toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Genome Res. 2010;20(9):1297–303.

    CAS  PubMed  PubMed Central  Google Scholar 

  79. 79.

    Der-Auwera GAV, Carneiro MO, Hartl C, Poplin R, Thibault J. From FastQ data to high-confidence variant calls: the genome analysis toolkit best practices pipeline. Curr Protoc Bioinformatics. 2013;11(1110):11–33.

    Google Scholar 

  80. 80.

    Depristo MA, Banks E, Poplin R, Garimella KV, Maguire J, Hartl C, Philippakis AA, del Angel G, Rivas MA, Hanna M, McKenna A, Fennel TJ, Kernytsky AM, Sivachenko AY, Cibulskis K, Gabriel SB, Altshuler D, Daly MJ. A framework for variation discovery and genotyping using next-generation DNA sequencing data. Nat Genet. 2011;43(5):491–8.

    CAS  PubMed  PubMed Central  Google Scholar 

  81. 81.

    Chen K, Wallis JW, Mclellan MD, Larson DE, Kalicki JM, Pohl CS, Mcgrath SD, Wendl MC, Zhang Q, Locke DP. BreakDancer: an algorithm for high-resolution mapping of genomic structural variation. Nat Methods. 2009;6(9):677–81.

    CAS  PubMed  PubMed Central  Google Scholar 

  82. 82.

    Abyzov A, Urban AE, Snyder M, Gerstein M. CNVnator: an approach to discover, genotype and characterize typical and atypical CNVs from family and population genome sequencing. Genome Res. 2011;21(6):974–84.

    CAS  PubMed  PubMed Central  Google Scholar 

  83. 83.

    Kai W, Li M, Hakon H. ANNOVAR: functional annotation of genetic variants from high-throughput sequencing data. Nucleic Acids Res. 2010;38(16):e164.

    Google Scholar 

  84. 84.

    Krzywinski M, Schein J, Birol I, Connors J, Gascoyne R, Horsman D, Jones SJ, Marra MA. Circos: an information aesthetic for comparative genomics. Genome Res. 2009;19:1639–45.

    CAS  PubMed  PubMed Central  Google Scholar 

  85. 85.

    Rubin C, Zody MC, Eriksson J, Meadows JRS, Sherwood E, Webster MT, Jiang L, Ingman M, Sharpe T, Ka S. Whole-genome resequencing reveals loci under selection during chicken domestication. Nature. 2010;464(7288):587–91.

    CAS  PubMed  Google Scholar 

  86. 86.

    Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira MAR, Bender D, Maller J, Sklar P, Bakker PIWD, Daly MJ. PLINK: a tool set for whole-genome association and population-based linkage analyses. Am J Hum Genet. 2007;81(3):559–75.

    CAS  PubMed  PubMed Central  Google Scholar 

  87. 87.

    Nei M. Genetic distance between populations. Am Nat. 1972;106(949):283–92.

    Google Scholar 

  88. 88.

    Kumar S, Stecher G, Tamura K. MEGA7: molecular evolutionary genetics analysis version 7.0 for bigger datasets. Mol Biol Evol. 2016;33(7):1870–4.

    CAS  Google Scholar 

  89. 89.

    Yang J, Lee SH, Goddard ME, Visscher PM. GCTA: a tool for genome-wide complex trait analysis. Am J Hum Genet. 2011;88(1):76–82.

    CAS  PubMed  PubMed Central  Google Scholar 

  90. 90.

    Alexander D, Novembre J, Lange K. Fast model-based estimation of ancestry in unrelated individuals. Genome Res. 2009;19(9):1655–64.

    CAS  PubMed  PubMed Central  Google Scholar 

  91. 91.

    Zhang C, Dong S, Xu J, He W, Yang T. PopLDdecay: a fast and effective tool for linkage disequilibrium decay analysis based on variant call format files. Bioinformatics. 2019;35(10):1786–8.

    CAS  PubMed  Google Scholar 

  92. 92.

    Bouckaert RR, Heled J, Kuhnert D, Vaughan TG, Wu C, Xie D, Suchard MA, Rambaut A, Drummond AJ. BEAST 2: a software platform for bayesian evolutionary analysis. PLoS Comput Biol. 2014;10(4):e1003537.

    PubMed  PubMed Central  Google Scholar 

  93. 93.

    Li H, Durbin R. Inference of human population history from individual whole-genome sequences. Nature. 2011;475(7357):493–6.

    CAS  PubMed  PubMed Central  Google Scholar 

  94. 94.

    Weir BS, Cockerham CC. Estimating F-statistics for the analysis of population structure. Evolution. 1984;38(6):1358–70.

    CAS  Google Scholar 

  95. 95.

    Danecek P, Auton A, Abecasis G, Albers CA, Banks E, Depristo MA, Handsaker RE, Lunter G, Marth GT, Sherry ST. The variant call format and VCFtools. Bioinformatics. 2011;27(15):2156–8.

    CAS  PubMed  PubMed Central  Google Scholar 

  96. 96.

    Mccarthy DJ, Chen Y, Smyth GK. Differential expression analysis of multifactor RNA-Seq experiments with respect to biological variation. Nucleic Acids Res. 2012;40(10):4288–97.

    CAS  PubMed  PubMed Central  Google Scholar 

  97. 97.

    Lun ATL, Chen Y, Smyth GK. It's DE-licious: a recipe for differential expression analyses of RNA-seq experiments using quasi-likelihood methods in edgeR. Methods Mol Biol. 2016;1418:391–416.

    PubMed  Google Scholar 

  98. 98.

    Schmittgen TD. Real-time quantitative PCR. Methods. 2001;25(4):383–5.

    CAS  PubMed  Google Scholar 

Download references

Acknowledgments

Not applicable.

Funding

This work was supported by the Earmarked Fund for China Agriculture Research System (No. CARS-40-K04), the Program for Innovation Research Team of the Ministry of Education (IRT16R23), the Key Science and Technology Research Project of Henan Province (112101110800), and the Young Backbone Teacher Training Program of Henan Higher Educational Institution (2017GGJS036). We are grateful to H.L. of Shanghai Majorbio Pharmaceutical Technology Co., Ltd. for their help with data analysis.

Author information

Affiliations

Authors

Contributions

XK and GS conceived and designed the experiments. YC, CZ, WL, FL, YF, YW and YT collected the samples and performed most of the experiments. GL, RH, RJ, ZL and XL analyzed the data and helped revise the manuscript. DL, MZ and WL wrote the paper. All authors reviewed and approved the final manuscript.

Corresponding authors

Correspondence to Wenting Li or Xiangtao Kang.

Ethics declarations

Ethics approval and consent to participate

All sample collection and treatment procedures were conducted in strict accordance with a protocol approved by the Institutional Animal Care and Use Committee (IACUC) of Henan Agricultural University, China (11–00929).

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary information

Additional file 1: Supplementary Figure S1.

Statistics of genetic variants in Xichuan black-bone chickens.

Additional file 2: Supplementary Figure S2.

Gene flow analysis of 29 chickens. Graphs were generated for one (A), two (B), three (C), and four (D) migration events. Migration events are shown by colored arrows. Migration weight is represented by color intensity.

Additional file 3: Supplementary Figure S3.

GO analysis of candidate genes under selection in Xichuan black-bone chickens.

Additional file 4: Supplementary Figure S4.

Verification of RNA-seq results by qRT-PCR.

Additional file 5: Supplementary Table S1.

Data sources for chickens used in this study.

Additional file 6: Supplementary Table S2.

Statistical summary for analysis of downloaded chicken sequence data (24 libraries).

Additional file 7: Supplementary Table S3.

Primers used for qRT-PCR in this study.

Additional file 8: Supplementary Table S4.

Statistical summary of sequencing data obtained for 5 XBC libraries by Illumina deep sequencing.

Additional file 9: Supplementary Table S5.

LD decay.

Additional file 10: Supplementary Table S6.

Time of divergence between populations.

Additional file 11: Supplementary Table S7.

Gene annotation of genomic regions with strong selective sweep signals in Xichuan black-bone chickens.

Additional file 12: Supplementary Table S8.

Difference in expression for dorsal skin between black skin (BS) and yellow skin (YS) in Xichuan black-bone chickens by RNA-Seq.

Additional file 13: Supplementary Table S9.

Coding potential of TCONS_00054154.

Additional file 14: Supplementary Table S10.

Bioinformatics analysis of TCONS_00054154

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Li, D., Sun, G., Zhang, M. et al. Breeding history and candidate genes responsible for black skin of Xichuan black-bone chicken. BMC Genomics 21, 511 (2020). https://doi.org/10.1186/s12864-020-06900-8

Download citation

Keywords

  • Xichuan black-bone chicken
  • Structural variants
  • Selective sweep
  • Black skin
  • Integration of whole genome and transcriptome