Genome-wide association study for carcass traits in Tsarskoye Selo chicken breed

Full Text

Introduction

Genotyping of farm animals is of great importance in modern selection programs, including poultry breeding. This is due to the fact that genotyping data are used to improve the accuracy of predicting the breeding value of an animal [1]. Furthermore, the presence of a range of genotypes allows researchers to perform a genome-wide association study (GWAS), which is designed to discover potential genomic regions that influence traits of interest [2]. GWASs reveal the relationship between genotype and phenotype by identifying single nucleotide polymorphisms (SNPs) with a high frequency of occurrence in individuals with similar trait values. Further study of the regions in which the identified SNPs are located makes it possible to identify candidate genes that influence the development of a particular trait. To date, a wide range of GWASs have been published in poultry, ranging from associations with various performance traits [3, 4], to plumage coloration [5], beak shape [6], and spur length [7].

 Carcass traits are an important group of chicken productive characteristics. They are mainly regulated by genetic factors such as genetic polymorphism and gene expression [8] and include a wide range of terms. Some researchers include in this group both quantitative and qualitative characteristics: weight values of muscle, fat and bone tissues, intramuscular fat content, the muscle proportion in the carcass, meat tenderness, amino acid composition, etc. [8]. However, an overwhelming number of studies refer to carcass traits as weight values of individual carcass parts [9-13]. This is probably due to the fact that weight values reflect the growth characteristics of an individual, i.e., they indicate its growth potential, which plays an important role in the economic efficiency of poultry meat production [8]. Although specialized broiler crosses are capable of providing the required amount of meat for the world’s population [14], there is a growing consumer demand for sustainable local products from alternative farming methods using local breeds [15, 16].

Tsarskoye Selo chicken breed has been developed by RRIFAGB researches on the basis of Center of Collective Use (CCU) «Genetic collection of rare and endangered chicken breeds» (RRIFAGB, Saint-Petersburg – Pushkin) [17]. As a dual-purpose breed, it has good meat characteristics, while maintaining optimal egg production. On average, the Tsarskoye Selo breed reaches a live weight of 3.5 kg for roosters and 2.4 kg for hens by 52 weeks of age [18]. In contrast to broilers that reach a live weight of 1.8 kg at 35 days of age [19], the Tsarskoye Selo breed weights 1.8 kg at 89 days of age [20]. Such a weight gain period has a favorable effect on the physiological homeostasis of the breed, because high productivity of broilers is often accompanied by an increased incidence of muscle myopathies [21]. The latter have a negative impact both on the commercial appearance of the carcass and on the technological properties of meat [22], but from an economic point of view, broiler poultry farming has better profitability. Nevertheless, the growing consumer interest in animal welfare focuses public attention on the ethical side of production [23]. Under the current circumstances, the need to optimize the production process is a matter of time, for which it is necessary to be prepared in advance.

In order to consider the Tsarskoye Selo breed as one of the parental forms for the production of specialized meat crosses, it is necessary to have a sufficient amount of genome data, as well as to identify genetic associations of genotype with phenotype. In the long term, these data may allow both effective selection in the Tsarskoye Selo breed and justify the feasibility of its use in local meat production. In this regard, the aim of the present study was to conduct a GWAS for carcass traits in Tsarskoye Selo chicken breed to establish the genetic determinants of meat productivity.

Materials and methods.

Healthy individuals of Tsarskoye Selo breed (n=96) from CCU «Genetic collection of rare and endangered chicken breeds» (RRIFAGB, Saint-Petersburg – Pushkin) were used as material for the study. At the age of 58 weeks all animals were slaughtered for postmortem collection of phenotypic data. During the slaughter, the blood sampling for DNA extraction was also performed. Postmortem collection of phenotypic data was done by weighting method and included the following parameters: live weight (LW), heart weight (H), liver weight (L), spleen weight, muscular and glandular stomaches weight, breast muscles weight (BM), thigh muscles weight (TM), shin muscles weight (SM), femur (FB) and sternum bone weight, and carcass without feathers weight (CWF).

DNA extraction from blood was performed according to the so-called “gold standard” - phenol-chloroform extraction method [24]. DNA concentration and purity of the samples were determined by spectrophotometry using the NanoDrop 2000c spectrophotometer (Thermofisher Scientific Inc., USA). The obtained samples were sent for whole-genome genotyping using Illumina Chicken 60K SNP iSelectBeadChip (Illumina Inc., USA) with a coverage density of 57,636 SNPs.  Genotyping quality control was performed using Genome Studio (Illumina Inc., USA). Samples with genotyping quality greater than 95% were used for further analysis. All samples (n=96) passed quality control. Editing of the obtained data to create adaptive extension files (.ped, .map, .fam, .bed, .bim) was performed using PLINK 1.9 software with minor allele frequency (MAF) > 0.05. This filter reduced the number of SNPs involved in the analysis to 43,123 SNPs.

The GWAS of genetic information with productive traits was performed using EMMAX statistical software. For this purpose, an «identity-by-state» kinship matrix was generated in EMMAX. The influence of SNPs on the trait was calculated according to the model:

Y=Xb+u+e,

where Y is the vector of phenotypes, b is the SNP effect, X is the estimated SNP genotype matrix, u is the vector of additive genetic effects, which is assumed to be normally distributed with mean equal to 0 and (co)variance σ2aG, where σ2a is the additive genetic variance and G is the genomic relationship matrix, e is the vector of random residual effects.

Bonferroni correction was applied to exclude false positives and to establish significance levels for the SNP effect. Significant and suggestive levels were set as 1.16E-06 (0.05/43.123) and 2.31E-05 (1.00/43.123), respectively. Genome-wide significance was assessed using the simple method in R, and the calculation of the effective number of independent tests was calculated using the Meff program.

Manhattan and quantile-quantile (Q-Q) plots were constructed based on GWAS results using the qqman and ggplot2 packages in the R programming environment; a correlation matrix for the studied traits was constructed using the corrplot package. Genes matching or close to the genomic region of the candidate SNP were annotated in the ENSEMBL genomic browser based on the Chicken (Red Jungle Fowl) GRCg6a genome assembly. Information on functional characterization of candidate genes was retrieved from the NCBI database.

Results.

As a result, the GWAS for carcass traits in Tsarskoye Selo breed identified associations for 8 traits out of 12. A total of 11 suggestive SNPs (2,31E-05) for 8 studied traits were identified on chromosomes (GGA) 1,3,11,12,15,22,23 and 27 (Table 1). The largest number of SNPs was identified for the TM trait - 3 SNPs on chromosomes GGA3, 12 and 27, and for BM - 2 SNPs on GGA1. One SNP each was detected for the remaining traits (Table 1). No genomic associations were found for such traits as spleen weight, sternum bone weight, and glandular and muscular stomach weight.  Additionally, correlations between traits with identified associations were analyzed (Fig. 1).

Figure 1. Correlation matrix for carcass traits of Tsarskoye Selo chicken breed: LW – live weight, CWF – carcass without feathers, H – heart, L – liver, TM – thigh muscles, BM – breast muscles, SM – shin muscles, FB – femur bone.

Table 1. Suggestive SNPs associated with carcass traits of Tsarskoye Selo chicken breed

Trait	SNP	Chromosome: Position	P-value	Allele	Localization	Candidate gene
Live weight (LW)	rs15204278	23: 4373046	6,20E-06	A/G	Intron variant	CLSPN
					Nearby:	AGO3
Carcass without feathers (CWF)	rs15204278	23: 4373046	8,23E-06	A/G	Intron variant	CLSPN
Heart weight (H)	rs14690448	22: 2002179	2,33E-05	T/C	Intron variant	UNC5D
Liver weight (L)	rs14023503	11: 9628291	7,53E-06	A/G	Intergenic variant	-
					Nearby:	GPATCH1
Thigh muscles (TM)	rs14406990	3: 103521235	1,10E-05	A/C	Intergenic variant	-
	rs313220886	27: 6441863	1,53E-05	A/G	Intron variant	SKAP1
	rs317772775	12:2661479	2,22E-05	T/C	Intron variant	DCAF1
Shin muscles (SM)	rs15773720	15: 6536930	1,91E-05	A/C	intron variant	ISCU
					Nearby:	TRAFD1
						RPL6
						CMKLR1
Breast muscles (BM)	rs15236831	1: 36634355	8,58E-06	T/C	Intron variant	TPH2
						TBC1D15
					Nearby:	RAB21
	rs316883857	1: 36606104	1,49E-05	A/G	intron variant	TBC1D15
					Nearby:	RAB21
Femur Bone (FB)	rs14882816	1: 124033301	2,71E-06	A/C	Intergenic variant	-
					Nearby:	GPM6B
						RAB9A
						TRAPPC2

Discussion.

The SNP rs15204278 (6,20E-06) on GGA23, localized in the intron of the CLSPN gene, whose product plays an important role in cell cycle control in response to replicative stress or DNA damage, was associated with LW and CWF traits [25] (Figure 2A, B). Correlation analysis also revealed a strong positive relationship between these two traits (Fig. 1). Previously, CLSPN was proposed as a candidate gene affecting egg production in a study on chicken population obtained by reciprocal crossing of white Leghorn and Dongxiang blue-shelled chicken [26]. Within 0.5 Mb of the target SNP, the AGO3 gene is localized, which plays an important role in RNA-mediated gene silencing. In a recent study, overexpression of the circular RNAs originating from exons of AGO3, circAGO3s, in chicken skeletal muscle markedly impairs myogenesis by upregulating genes associated with muscle atrophy [27], suggesting a potential effect of AGO3 on myogenesis and meat productivity. Accordingly, the identification of SNPs in the CLSPN gene associated with LW and CWF is most likely due to the proximity of the latter to the AGO3 gene.

 

Figure 2. – Manhattan plots for LW (A), CWF (B), H (C), L (D) traits in Tsarskoye Selo chicken breed

An association with SNP rs14690448 (2.33E-05) on GGA22, localized in an intron of the UNC5D gene encoding the netrin receptor, was identified for the H trait (Figure 2C). Netrins are proteins associated with extracellular laminin that play a role in cell migration processes, as well as angiogenesis and morphogenesis in animals and humans [28, 29]. Another association was found for SNP rs14023503 (7.53E-06) on GGA11 with the L trait (Figure 2D). This SNP was localized in the intergenic space, and within 0.5 Mb of it was located the GPATCH1 gene, which presumably provides RNA binding activity. Notably, in a 2021 study on the population of native African chickens from Rwanda, the GPATCH1 gene was found in one of the promising genomic regions controlling body weight [30].

Three associations were identified for the TM trait, namely rs14406990 (1.10E-05), rs313220886 (1.53E-05) and rs317772775 (2.22E-05) on GGA3, 27, and 12, respectively (Figure 3A). SNP rs14406990 was localized in the intergenic space, but no genes associated with exterior or growth parameters were found near it, whereas SNP rs313220886 was located in the intron of the SKAP1 gene, and rs317772775 - in the intron of DCAF1.

 The SKAP1 gene plays multiple roles in the regulation of integrin activation, “stop-signaling” and cell cycle optimization of proliferating T cells through interaction with various mediators [31], whereas DCAF1 is involved in the regulation of cellular immune responses [32].

Figure 3. – Manhattan plots for TM (A), SM (B), BM (C), FB (D) traits in Tsarskoye Selo chicken breed

In a recent GWAS for bone quality in Nonghua ducks, a genomic region containing the SKAP1 gene was associated with bone parameters such as metatarsus, tibia, and femur weight, as well as tibia and femur length [33]. The existing data are consistent with our study, where an association with TM was found for the SNP in SKAP1 gene. The explanation for this could be that the increase of femur area caused by bone growth in length and width creates a larger area for muscle attachment, thereby determining the potential for growth [34]. This assumption is additionally confirmed by the results of correlation analysis, where a strong positive correlation was found for TM and FB traits (Fig. 1)  The SNP rs1577373720 (1.91E-05) on GGA15, localized in the intron of ISCU gene, was identified for the SM trait (Figure 3B). ISCU encodes a component of the iron-sulfur cluster (Fe-S) framework, which play a role in the function of a diverse set of enzymes, including those that regulate metabolism, iron homeostasis, and response to oxidative stress. Human studies have identified a number of mutations in ISCU leading to myopathies accompanied by episodes of muscle weakness, rhabdomyolysis and lactoacidosis [35, 36]. This is explained by the fact that mutations in ISCU gene lead to an increase in the level of iron and reactive oxygen species within mitochondria, which cause oxidative stress and thus impair mitochondrial function [37]. It’s also worth mentioning that ISCU was identified as one of the 111 differentially expressed genes in a study of molecular mechanisms of growth depression in broilers in response to immunologic stress [38]. Taken together, these data suggest that ISCU is capable of mediating effects on meat performance in chicken, particularly on SM trait. The TRAFD1, RPL6, and CMKLR1 genes are located within 0.3 Mb of the target SNP rs1577373720. TRAFD1 is a negative feedback regulator controlling excessive innate immune responses. Thus, in monocytes and macrophages, TRAFD1 suppresses inflammatory responses of innate immunity by inhibiting NF-κB and MAPK activation [39]. In recent studies, TRAFD1 is annotated as a candidate gene associated with intramuscular fat content in Beijing black pigs [40]. RPL6 encodes a component of the large ribosomal subunit, the main function of which is protein synthesis [41]. RPL6 was proposed as a reference gene for RT-qPCR in adipose tissue and muscle tissue of the longissimus dorsi muscle in buffalo [42], and was also associated with brain deformities in Danio Rerio during embryonic development [43]. CMKLR1 encodes a protein that enables adipokinetic hormones binding activity, as well as the activity of their receptors. In particular, CMKLR1 is the main receptor of chemerin, a hormone produced mainly by adipose tissue and liver [44]. In existing studies on chickens, a potential role of CMKLR1 in stimulating the growth of chicken embryos was noted [45], which made the chemerin system a potential candidate to explain the differences in growth and metabolism of laying hens and broilers [46]. The importance of chemerin and its receptors for embryogenesis was also observed in Beijing ducks [47], which may indicate that this system is highly conserved. The assumption about the influence of the system of chemerin and its receptors on the growth and metabolism of animals is also supported by the study on mice. Thus, mice with CMKLR1 knockout under both low- and high-fat diets were characterized by reduced food intake, total body weight, and fat percentage, compared to wild-type control individuals [48].

 

Two SNPs on GGA1, rs15236831 (8.58E-06) and rs316883857 (1.49E-05), were identified for the BM trait (Figure 3C). SNP rs15236831 is localized in the intron part of TPH2 and TBC1D15 genes, whereas rs316883857 is localized only in the intron of TBC1D15 gene.  TPH2 encodes the tryptophan hydroxylase and is involved in serotonin biosynthesis. Verterbates exhibit two serotonin systems that are regulated by different enzymes (TPH1 and TPH2) and, accordingly, have different functions. For example, TPH2 is involved in central nervous system effects such as aggression, anxiety, food intake and sleep [49]. This is supported by the results of a study on the effect of gut microbiota on thermogenesis in broilers, where an increase in body weight and THP2 gene expression in the hypothalamus, as well as improved feed conversion were observed in the group reared at high temperatures using antibiotics [50]. Similar data were obtained in a study of the hypothalamic transcriptome profile in broilers, where TPH2 was noted as one of the differentially expressed genes associated with appetite [51]. The existing data are consistent with our study and suggest that the association between BM and the TPH2 gene is due to the latter's influence on eating behavior. The TBC1D15 gene encodes a protein acting as a GTPase activator. A recent study noted that overexpression of TBC1D15 disrupts lysosomal morphology and blocks growth factor withdrawal-induced cell death [52]. Another study provided convincing evidence that TBC1D15 serves as the main regulator of GLUT4 (Glucose transporter type 4) translocation and further affects GLUT4-mediated glucose uptake [53]. The latter indicates that TBC1D15, affecting glucose uptake, influences metabolic processes in the organism. Within 0.1 Mb of both SNPs is located the RAB21 gene, which encodes a GTPase involved in the control of cell membrane trafficking. A number of candidate genes associated with body mass index (BMI) have been obtained in a whole-genome study in humans, among which RAB21 has been noted [54]. Thus, the genes located near both SNPs form a haplogroup influencing body metabolic processes, which explains their association with the BM trait.

In addition to the above, SNP rs14882816 (2.71E-06) associated with the FB trait was obtained on GGA1 (Figure 3D). The GPM6B, RAB9A, and TRAPPC2 genes were localized within 0.4 Mb of the target SNP located in the intergenic space. GPM6B encodes a neuronal membrane glycoprotein putatively involved in protein transport, neurogenesis, and osteogenesis, particularly in osteoblast differentiation. Osteoblasts are cells of mesenchymal origin, responsible for bone formation and mineralization. Mineralization is ensured due to the release of special organelles - matrix vesicles - from the membrane of osteoblasts [55]. Drabek et al. established that GPM6B expression is upregulated in primary human mesenchymal stem cells during osteogenic differentiation and in human fetal preosteoblast cell line, while GPM6B silencing leads to inhibition of mineralization and alteration of the actin cytoskeleton [56].  This is not only consistent with the results of our study, but also allows to identify GPM6B as a novel regulator of osteoblast function and bone formation. RAB9A encodes a gene that shows increased expression during osteoblast differentiation [57] and also induces osteoclast formation and bone resorption [58]. Finally, the last gene near the target SNP, TRAPPC2, encodes a protein that is part of a large multi-subunit complex involved in vesicular transport between the endoplasmic reticulum (ER) and the Golgi. Mutations in this gene are known to cause spondyloepiphyseal dysplasia tarda (SEDT) [59, 60]. This disease is a genetically heterogeneous hereditary osteochondropathy characterized by malformations of vertebrae, epiphyses of tubular bones, and joints. Thus, all genes near SNP rs14882816 associated with FB trait are involved in the processes of osteogenesis and bone homeostasis.

Conclusion.

The present GWAS analysis of Tsarskoye Selo chicken breed identified associations for 8 out of 12 carcass traits. A total of 11 suggestive SNPs (2,31E-05) were identified on GGA1,3,11,12,15,22,23 and 27. The highest number of SNPs was found for thigh muscles weight (TM) with 3 SNPs and for breast muscles weight (BM) with 2 SNPs. One SNP each was detected for the remaining traits. A total of 16 genes associated with immunity (SKAP1, DCAF1, ISCU, TRAFD1), metabolism (GPATCH1, CMKLR1, TBC1D15, RAB21), osteogenesis (GPM6B, RAB9A, TRAPPC2), protein synthesis (RPL6), serotonin biosynthesis and eating behavior (TPH2), myogenesis (AGO3), morphogenesis (UNC5D), and DNA damage response (CLSPN) were identified. The results obtained can be used in the selection of Tsarskoye Selo chicken breed. Recommendation of these candidate genes in selection programs for other chicken breeds requires preliminary approbation of the results obtained on other meat and dual-purpose chicken breeds.

About the authors

Anastasiia I. Azovtseva

Russian Research Institute of Farm Animal Genetics and Breeding — Branch of the L.K. Ernst Federal Research Center for Animal Husbandry (RRIFAGB)

Author for correspondence.
Email: ase4ica15@mail.ru
ORCID iD: 0000-0002-2963-378X
SPIN-code: 5784-2786

Junior Researcher, Laboratory of Molecular Genetics

Russian Federation, Saint-Petersburg, Tyarlevo settlement, Moskovskoe shosse, 55a

Anna E. Ryabova

Russian Research Institute of Farm Animal Genetics and Breeding — Branch of the L.K. Ernst Federal Research Center for Animal Husbandry (RRIFAGB)

Email: aniuta.riabova2016@yandex.ru
ORCID iD: 0000-0003-2362-2892
SPIN-code: 4336-0310

Junior Researcher, Laboratory of Molecular Genetics

Russian Federation, Saint-Petersburg, Tyarlevo settlement, Moskovskoe shosse, 55a

Natalia V. Dementieva

Russian Research Institute of Farm Animal Genetics and Breeding — Branch of the L.K. Ernst Federal Research Center for Animal Husbandry (RRIFAGB)

Email: dementevan@mail.ru
ORCID iD: 0000-0003-0210-9344
SPIN-code: 8768-8906

Head of the Laboratory of Molecular Genetics, Leading Researcher, Cand. Biol. Sci.

Russian Federation, Saint-Petersburg, Tyarlevo settlement, Moskovskoe shosse, 55a

References

Supplementary files