In current study, whole genome sequences of four Tibetan chicken populations reared in Tibetan plateau were analyzed by focusing on genetic diversity, run of homozygosity, genomic inbreeding and selection signature. A composite Tibetan local breed, Lhasa white was also included in the analysis to compare results among populations. Lhasa white is a synthetic breed generated by crossing male White Leghorn and female Tibetan chicken, which has more than sixty years of history in Tibetan plateau.
Observed (Ho) and expected (He) heterozygosity were effective indices to evaluate the genetic diversity within populations. We calculated Ho as mostly previous studies using SNPs with MAF ≥ 0.05, and found values close to 0.3 in all populations, which are similar to values estimated in modern chicken populations using sequence data [15] and are higher than values estimated in Italian local chickens [12]. Moreover, when we kept all SNPs in the analysis to avoid any bias [39], we found similarly lower Ho for all populations compared to Ho reported in Dutch local chickens evaluated by chicken 60K SNP arrays [10]. In our study, we observed slightly lower heterozygosity than expected heterozygosity in SH, NM, DZ and LW, suggesting subtle inbreeding in these populations. However, a little heterozygosity excess (Ho > He) existed in LZ population. This may indicate a recent bottleneck or an isolate-breaking effect [40] which may likely be due to the recent domestication and selection process. According to Wright’s interpretation of Fst [41], pair-wise Fst among SH, NM and DZ were less than 0.05, indicating a little genetic differentiation. LZ was moderately distant from SH, NM and DZ, suggesting that LZ population was bred with little or no admixture in a relatively isolated environment. These findings confirmed that two or more Tibetan chicken populations live in the plateau [2, 4].
In present study, number and the distribution of ROHs identified for Tibetan native chickens were comparable to that previously reported in broiler [14]. Most ROHs identified in our study belonged to the short class, indicating that little deleterious inbreeding happened in five chicken populations [42]. The relationship between total number of ROH and total length of the genome covered by ROH showed considerable variation among animals within and across populations. Similar distributions were also reported in other livestock species, such as indigenous sheep [17] and cattle [43]. Genomic data is the only reliable source to estimate inbreeding and relatedness of marginalized populations in the absence of other data sources, such as pedigree. It is commonly accepted that the proportion of the genome within ROH (FROH), genomic SNP-by-SNP inbreeding coefficient (FGRM), excess of homozygosity (FHOM) and correlation between uniting gametes (FUNI) were indicators for inbreeding assessment [44]. Herein, the ROH-based genomic inbreeding coefficients of Tibetan chicken were similar to estimates in other Chinese local chickens that were under conversation [11], while the values were much lower than those in modern chickens [15, 16] and other local chickens from Italy [12]. This suggests that Tibetan chickens maintained their natural diversity in the plateau. The correlation between FROH and FHOM, as well as FUNI were significantly high, similar to those estimated in cattle [45], pig [46], horse [47] and modern chicken [15], indicating that FROH can be used as an accurate estimate of the proportion of IBD (identity by descent). However, the correlation between FGRM and FROH and that between FGRM and FHOM close to 0 may probably be due to the number of SNPs used in this study, since sequences generated much more SNPs than SNP arrays, which significantly affected FROH estimates [42, 48]. Similar situation was reported in other domestic animals [49, 50]. Moreover, the FHOM measurements that employed homozygous sites in the genome were more sensitivity to allelic frequencies [42].
Whole genome sequencing allows ROH to be more reliably detected, and analyzing the effect of common ROHs may reveal the demographic history of animal populations and allows exploration and testing of new approaches to understand complex traits [5]. The ROH islands of Tibetan chickens and Lhasa white chickens harbored many QTLs and candidate genes controlling economically important traits, including conformation, production, egg and meat quality, digestion and absorption, reproduction and growth traits. Regarding common genes located on GGA5, leucine-rich repeat-containing G protein-coupled receptor 4 (LGR4), enriched in biological process of osteoblast in GO database, contributes to regulation of energy metabolism including food intake and energy expenditure [51]. Brain-derived neurotrophic factor (BDNF), enriched in positive regulation of synapse assembly is considered important for the temperature perception in chicken [52]. In rats, BDNF administration in the paraventricular nucleus reduced energy intake and decreased body weight [53]. STAT1 and STAT4 are members of Janus kinase (JAK)-signal transducer and activators of transcription (STAT) pathway that plays critical roles in facilitating various cellular reactions to diverse forms of cellular stress, including hypoxia/reperfusion, endotoxin, ultraviolet light, and hyperosmolarity [54]. Moreover, these genes were previously identified as ROH islands-associated genes in Italian [12] and Mexican native chickens [20], suggesting that ROH-related homeostasis and STATs pathway are important in highland chickens. Particularly, metal ion binding was enriched in 34 genes. Although the process of how metal ion binding affect animal’s physiology and production is rarely reported, some genes enriched in the term including VAV3, NOS2, COL3A1 and PRKD1 were putative candidate genes associated with highland adaptation [55], implying that metal ion binding may be associated with highland adaptation.
ROH islands might be indicative of genomic regions underwent natural and artificial selection [56]. The iHS approach appears to be the most powerful for detecting ongoing selection processes for which the target allele has a moderate to high frequency within a population. If the iHS method detects a genomic region, this region can contain several loci that may be undergoing selection within the breed. Therefore, the iHS method can detect breed-specific candidate genes under selection [57]. Our iHS analysis revealed that the common genomic region with different ROH islands on GGA8 were overlapped with putative selection signature in SH, NM, LZ and LW populations, indicating selective forces were undergoing in the regions. Commonly identified regions by both iHS and ROH analysis harbored AMY2A, NTNG1 and VAV3 genes. While AMY2A encodes a member of the alpha-amylase family of proteins, which is involved in carbohydrates and glycogen metabolism, affecting growth, carcass traits and feed intake efficiency in chicken [58]. Previous report shown that AMY2A was under selection for metabolism, energy availability and response to thermal stress in African chickens [59]. Similar to African village conditions, chicken feeding is mainly based on scavenging, household waste and some grain supplementation in the Tibetan plateau. Therefore, carbohydrate metabolism, energy generation and transport are important traits for feeding adaptation. NTNG1 is a responsible gene for axon and neurite growth [60] and has been linked to body weight gain of beef cattle in a meta-analysis [61], as well as birth weight of pigs in a GWAS research [62]. The gene was also differentially expressed in chicken hepatocellular cell line in response to stress [63]. VAV3 is a member of the VAV gene family that play vital roles as guanosine nucleotide exchange factors for Rho GTPases and signaling adaptors downstream of protein tyrosine kinases [64], and it regulates osteoclast function, bone mass, hypothyroidism and renal systems [65]. Specifically, VAV3 has been identified as candidate gene associated with highland adaptation in Ethiopians [66] and Ethiopian sheep [67]. This probably resulted from the role it plays in homeostasis of the cardiovascular [68]. We therefore suggest that VAV3 also functions putatively in the adaptation to high altitude of Tibetan chicken.