A significant quantitative trait locus on chromosome Z and its impact on egg production traits in seven maternal lines of meat-type chicken

Jiqiang Ding, Fan Ying, Qinghe Li, Gaomeng Zhang, Jin Zhang, Ranran Liu, Maiqing Zheng, Jie Wen and Guiping Zhao*

Abstract Background: Egg production is economically important in the meat-type chicken industry. To better understand the molecular genetic mechanism of egg production in meat-type chicken, genetic parameter estimation, genome-wide association analyses combined with meta-analyses, Bayesian analyses, and selective sweep analyses were performed to screen single nucleotide polymorphisms (SNPs) and other genetic loci that were significantly associated with egg number traits in 11,279 chickens from seven material lines.Results: Yellow-feathered meat-type chickens laid 115 eggs at 43 weeks of age and white-feathered chickens laid 143 eggs at 60 weeks of age, with heritability ranging from 0.034–0.258. Based on meta-analyses and selective sweep analyses, one region (10.81–13.05 Mb) on chromosome Z was associated with egg number in all lines. Further analyses using the W2 line was also associated with the same region, and 29 SNPs were identified that significantly affected estimation of breeding value of egg numbers. The 29 SNPs were identified as having a significant effect on the egg number EBV in 3194 birds in line W2. There are 36 genes in the region, with glial cell derived neurotrophic factor, DAB adaptor protein 2, protein kinase AMP-activated catalytic subunit alpha 1, NAD kinase 2, mitochondrial, WD repeat domain 70, leukemia inhibitory factor receptor alpha, complement C6, and complement C7 identified as being potentially affecting to egg number. In addition, three SNPs (rs318154184, rs13769886, and rs313325646) associated with egg number were located on or near the prolactin receptor gene.Conclusion: Our study used genomic information from different chicken lines and populations to identify a genomic region (spanning 2.24 Mb) associated with egg number. Nine genes and 29 SNPs were identified as the most likely candidate genes and variations for egg production. These results contribute to the identification of candidate genes and variants for egg traits in poultry.

Keywords: Candidate gene, Egg production, Genome-wide association study, Meat-type chicken, Selective sweep analysis

Introduction

Throughout the twentieth century, specialized commercial chicken populations were established for egg laying and meat production [1]. Currently, layers display high levels of egg production and have been developed by intensive artificial selection in the egg-laying industry to yield more than 320 eggs during a 52-week laying period [2, 3]. Meat-type chickens have received substantial attention for their meat production, and to avoid selection for both growth traits and egg production traits in the same bird [1]. For meat-type chickens, the characteristics of body weight, growth rate, and feed conversion rate have been widely studied [4].

Although significant genetic changes may have occurred during selective breeding practices, the relationships between genotypes and phenotypes are not well characterized [2]. Many studies have shown that the heritability of egg number production. For example, such heritability was 0.13–0.45 and 0.10–0.40 for egg number production in two commercial lines of White Leghorn hens [5]. In the past three decades, many studies have explored egg production at the genetic level. Recent studies of reproductive hormones have identified 31 candidate genes significantly associated with egg production, and there are 64 new candidate genes and 108 SNPs associated with oviposition performance based on genome sequencing analyses [6]. Over 440 quantitative trait loci (QTL) on theGallus galluschromosome (GGA) have been reported to be associated with egg number in chickens [7]. With the deeply sequencing and the mining of increasing SNPs, genome-wide association studies (GWAS) and selective sweep analyses have become powerful tools for detecting QTLs in many breeds of livestock and poultry [8].

The growth and feed conversion traits of meat-type chickens have been widely reported, however, egg-laying traits have been less studied. For layers, multiple SNPs and QTLs have been found [6], but due to the genetic differences between layers and meat-type chickens, the application of those findings in meat-type chickens has been limited. Mapping of QTLs and SNPs with potential impact on egg production in meat-type chickens to provide more reference information for subsequent genome selection. Our experimental design involved GWAS and selective sweep analyses to detect genetic variations associated with egg number based on a 55 K SNP Chicken chip in a population of 11,279 meat-type chickens.

Materials and methods

Animals and phenotypes

A total of 11,279 meat-type chickens were obtained from seven lines of two breeding companies in China. The W1, W2, W3, and W4 were four commercial white-feathered pure lines in Yunnan Province, China. Body weight exceeded 2.4 kg at 40 days of age in lines W1, W2, and W4. Body weight exceeded 1.6 kg at 40 days of age in line W3. The feed conversion ratio was under 1.69 between 28 and 40 days of age, and the age at first egg was 175–180 d.

The Y1, Y2, and Y3 were three commercial yellow-feathered pure lines in Jiangsu Province, China. The feed conversion ratio was ~ 3.13, which was considered medium growth rate meat-type chickens.

Line W2 was selected for the egg production more than eight generations. We used the individuals of three generations (3159 chickens: 2532 females; 627 males), including the 6th, 7th and 8th generations. All other lines were used for individual data of one generation. All chickens had pedigree, genomic information and were housed in identical individual cages. There were 11,279 individuals in all lines, of which 9279 were recorded for egg numbers. Blood samples were collected from the brachial veins of each chicken using the standard procedure of the breeding program, which was approved by the Animal Welfare and Ethics Committee of the Institute of Animal Sciences, Chinese Academy of Agricultural Sciences (IAS-CAAS, Beijing, China).

Staged egg numbers statistics

We made a detailed division of the laying period using 627 hens from the 7th generation of line W2. We recorded the age at first egg (AFE), and egg production for each hen from first egg to 354 days of age. We divided the egg-laying period into four stages according to the laying curve in Additional file 1: Fig. S1 and counted the egg numbers (EN) at each stage. Four stages, including the stage of rapidly increasing egg production (from the onset of laying eggs to 195 days of age; EN1); the stage of peak egg production (from 196 to 227 days of age; EN2); the stage at which the laying rate was between 60% and 75% (from 228 to 307 days of age; EN3); and the stage of decreasing egg production (from 308 to 354 days of age; EN4). Moreover, we counted the total egg number from first egg to 354 d (TEN).

Genotyping and quality control

Genomic DNA was extracted from blood samples using the phenol-chloroform method. All the birds were genotyped using a customized 55 K SNP chicken array obtained from Beijing Compass Biotechnology Co., Ltd. (Beijing, China) [9]. We updated the SNP positions using the newest release from the University of California-Santa Cruz (GRCg6a/galGal6 genome version). The following quality control criteria were applied to the target panel: individual call rate ≥ 90%, SNP call rate ≥ 90%, and minor allele frequency ≥ 0.05. Ultimately, 37,847–45,346 SNPs and 11,279 birds remained for further analyses.

Population structure testing

Twenty individuals randomly selected from each line for structure using ADMIXTURE software (version1.3.0) [10]. A principal component analysis (PCA) plot was constructed using the “scatterplot3d” package in R (version 4.0.2).

Genome?wide association study

GWAS was performed using 9778 individuals with egg number phenotypes from the all lines. The GWAS for egg number traits was performed using the univariate linear mixed model implemented in GEMMA software, version 0.98.1 [11]. The statistical model was as follows:

whereyis the vector of phenotypic values;Wis the vector of covariates, including a column of 1 s;αis the vector of the corresponding coefficients, including the intercept;xis the vector of marker genotypes;βis the effect size of the marker;uis the vector of random polygenic effects;εis the vector of errors;τ?1is the variance of the residual errors;λis the ratio between the two variance components;Kis the centered relatedness matrix estimated from SNPs, and Ιnis the identity matrix.MVNnis then-dimensional multivariate normal distribution. The Wald test was used to select SNPs associated with egg production traits.

We used the parameters of plink software --indeppairwise 25 5 0.2 to infer effective independent tests. Considering the over-conservation of the 5% Bonferroni correction method, we adjusted the threshold of genome-wide significantP-values based on the number of independent SNPs [12]. Manhattan and Q-Q plots were constructed for each trait using the “qqman” package in R (version 4.0.2). Boxplots were produced by the “ggplot2” package in R (version 4.0.2).

Meta?analyses

Meta-analyses of the seven lines were performed using standard method by Metal software (version 2011-03-25). Meta-analyses can combine either (a) test statistics and standard errors or (b)P-values across multiple GWAS for a single trait (taking sample size and direction of effect into account) [13]. The detailed calculation formula is introduced in reference [13]. Meta-analysis results in little or no loss of efficiency compared to analysis of a combined dataset including data from all individual studies.

Association study on GGA Z

Due to the particularity of GGA Z, we used “-method threshold” function of SNPTEST v2.5.6 software [14], a software for sex chromosome association study, for GGA Z association study. For chickens, female genotypes are encoded as 0/1 in SNPTAST. In addition to association test statistics, SNPTEST can output expected genotype and allele counts for diploid and haploid samples. Computation of allele frequencies and info statistics also take into account ploidy.

Selective sweep analysis

Selective sweep analysis was performed on all hens from the line W1 and 782 hens from the 8th generation of line W2, which had genomic information obtained by 55 K SNP chips. Line W2 was selected for the egg number more than eight generations, Line W1 was not selected for egg number. We used the “--weir-fst-pop” function of the Vcftools v0.1.14 software to calculate the population differentiation index (Fst), which were a window length of 40 kb and a step size of 2 kb [15].

Bayesian analysis

Two Bayesian approaches, Bayes B and Bayes R, were used to obtain the estimated marker effect explained by each SNP. Bayes B and Bayes R analyses were performed using the BGLR and hibayes package in R, respectively. The number of iterations after the burn-in phase and that of the burn-in period were 20,000 and 10,000, respectively. The purpose of the Bayesian analysis was to estimate whether the QTL we mapped had significantly larger effects. The Bayesian analysis’ result was the effect size of the SNPs, while the result of GWAS was thePvalue and effect size of SNPs associated with phenotype. Their results could be verified by each other.

Genome prediction models

The genome prediction models, genomic best line prediction (GBLUP) and single-step GBLUP (ssGBLUP), were utilized in this study and were as follows:

Whereyis the vector of phenotypic value,XandZare the association matrices of fixed and additive genetic effects,bis the vector of the fixed effect,αis the vector of the random additive genetic effect, andeis the vector of the random residual error. For GBLUP, when anαobeyed the following normal distributionα~N(0,Gσ2u), G was the consanguinity matrix (G matrix). In ssGBLUP, the H matrix replaces the G matrix.

H matrix construction:

WhereHis an inter-individual relational matrix using genome-wide markers and pedigrees,Ais thepedigree-based relatedness matrix,Gis the genomebased relationship matrix, andA22is the individual genealogical relation matrix of sequencing individuals. The relative weights ofGandA22in the H matrix were set asGω= 0.95 ×G+ 0.05 ×A22. The level of theGmatrix was corrected to theA22(matrix)

G?represents the adjustedGmatrix, andaandbare as follows:

The formula of theHmatrix after merging is:

The animal models of the two traits are as follows:

Wherey,X,b,Z,α, andeare the same formula (2).

The genetic parameters of line W2 for egg production traits were calculated using the ssGBLUP model, and the heritability of the seven lines was estimated by the GBLUP model. GBLUP and ssGBLUP were calculated using ASReml v4.1 software [16].

Results

Basic statistics and genetic parameters

A total of 9778 hens from seven lines were subjected to genome-wide association analyses. The egg number of each line was evaluated at a different time. The heritability of egg number varied between the seven lines. The results showed that the egg number traits exhibited low heritability (0.034–0.258) (Table 1).

Table 1 Basic statistics for egg production traits in the seven lines

The basic statistics for the distribution of egg production traits from the 627 hens used in the GWAS and genetic parameters are shown in Table 2. EN1 and EN4 had the largest phenotypic variation, as some chickens were late to laying and laid eggs at inconsistent intervals. Additionally, some birds died or were too old to lay eggs in EN4 because of their low rate of follicle development[17]. A total of 7771 independent tests of all chromosomal SNPs were performed. Therefore, the thresholdP-values for genome-wide suggestive and genome-wide significance were calculated as 1.29E-04 (1.00/7771) and 6.43E-06 (0.05/7771), respectively.

Table 3 Estimates of genetic parameters for egg production traits in 627 hens of line W2

The heritability and correlation of egg production traits was estimated using a relationship matrix composed of pedigrees and genotypes. AFE and EN1 had medium heritability (0.24 and 0.27), while the heritability of ENs in other egg-laying stages was low, ranging from 0.09 to 0.17. AFE had a high negative correlation with EN1, a moderately negative correlation with TEN, and a low correlation with ENs in other egg-laying stages. It was shown that AFE affected the laying performance by affecting the EN1. EN1 had a moderately positive correlation with TEN, and low positive correlations with ENs in other egg-laying stages. EN2, EN3, EN4, and TEN exhibited high correlations with each other (> 0.8) (Table 3).

Population structure testing

PCA using the first three principal components showed obvious population stratification between the seven lines and that four lines (W1, W2, W3, and W4) of whitefeathered meat-type chickens were mixed with each other (Fig. 1).

Figure S2 displays a bar plot based of the cross-validation error rate. When K = 2, the fast-growing whitefeathered and the yellow-feathered meat-type chickens appeared as two differentiated clusters. When K = 3, the four white-feathered meat-type chicken lines were divided into three groups. When K = 4, the three yellowfeathered populations were separated (Additional file 3: Fig. S3).

Fig. 1 PCA plots for the seven lines. The first three PCs were used to calculate the genetic relationships between the seven lines. Y1, Y2, and Y3 represent three lines of yellow-feathered meat-type chickens, W1, W2, W3 and W4 represent four lines of white-feathered meat-type chickens. The four white-feathered meat-type chicken lines were closely related, and three yellow-feathered meat-type chicken lines 101, 401 and 501 were closely related

GWAS for single?line egg number traits

GWAS was performed on each line, and 35 SNPs were associated with egg number. Using the Ensembl database to annotate SNPs, we found 31 genes around significant SNPs (Additional file 6: Table S1). For the GWAS of the seven lines, the genomic inflation factors (λ, 0.993–1.069) were calculated to evaluate the accuracy and reliability of the GWAS results (Additional file 4: Fig. S4).

Meta?analysis based on GWAS results of single?line egg number traits

Through the above analysis, the SNP effect of egg numbers in a single line could be obtained, but it was unknown whether the same region exists for multiple lines with large genetic background differences. Therefore, using the above GWAS results, a meta-analysis was performed on all individuals of the seven lines, and an obvious bulge region was found on the Z chromosome (Fig. 2A).

Significant region was differentiated between lines with different breeding directions

The Fst and major allele frequency analyses were performed on the line W1 and W2 chickens (Fig. 2B and C). It was a significant signal of selection that a QTL region overlapped with the result of meta-analysis located between 10.81 Mb and 13.05 Mb on the GGA Z. The region contains 36 genes. Eight genes were identified as candidate genes after functional annotation. The genes included glial cell derived neurotrophic factor (GDNF), DAB adaptor protein 2 (DAB2), protein kinase AMPactivated catalytic subunit alpha 1 (PRKAA1), NAD kinase 2, mitochondrial (NADK2), WD repeat domain 70 (WDR70), leukemia inhibitory factor receptor alpha (LIFR), complement C6 (C6), and complement C7 (C7).

Fig. 2 Selective sweep analyses and meta-analyses were located in a region of the GGA Z. A Manhattan plots of the meta-analyses for egg numbers of the seven lines; each dot represents one SNP. B The visualization of Fst in lines W1 and W2; each dot represents a calculation step size of 2 kb. C The major allele frequency analyses of SNPs located in the candidate region for lines W1 and W2

Phased egg number traits to verify significant QTL region

(1) GWAS for staged egg number traits

Combined with the results of GWAS and selective sweep analysis, we further analyzed the region between 10.81 Mb and 13.05 Mb on the GGA Z in line W2. For the GWAS of 627 hens in line W2, Manhattan and quantile–quantile (Q-Q) plots for egg production traits were generated, as shown in Fig. 3 and Additional file 5: Fig. S5. Detailed information on the SNPs significantly associated with EN2, EN3, EN4, and TEN is shown in Table 4. For the EN2, 26 genome-wide significant SNPs were identified on GGA Z, and 38 genome-wide suggestive SNPs were detected on GGA Z, GGA1, GGA3, GGA4, GGA10, and GGA13. The genomic inflation factor (λ) was 1.00 for EN2, which suggested that the population structure was well controlled.

For the TEN, one significant SNP and 38 suggestive SNPs were detected on GGA Z, GGA1, GGA3, and GGA23. The SNPs on the GGA Z were also associated with EN2. Three SNPs (rs318154184, rs13769886, rs313325646) associated with EN2, EN3 and TEN, were located on or near the prolactin receptor (PRLR) gene on the GGA Z. There were 33 coincidently associated SNPs for EN2 and TEN between 10.81 Mb and 13.05 Mb on GGA Z, among which 29 SNPs were associated with significantly different egg number EBVs of the 3159 individuals in line W2. The 29 SNPs were located in four LD blocks (Fig. 4), and the visualization for analysis of variance of one SNP selected each Block.

(2) Association tests on GGA Z

Considering the particularity of GGA Z, we used SNPTEST software to analyze GGA Z of EN2 and TEN traits of 627 hens in line W2 specially. The results (Fig. 5) were consistent with those in (1), which could indicate the reliability of the region between 10.81 Mb and 13.05 Mb on GGA Z.

(3) The effect of each SNP

We used Bayes B and Bayes R to estimate the effect value of each SNP using the EN2 and TEN phenotypes of 627 hens from line W2. There were obvious bulges in the significant QTL region, which was similar to GEMMA in that the calculated effect size also had obvious bulges in the region (Additional file 6: Fig. S6). In the effect estimated for this region, GEMMA, Bayes B and Bayes R estimated EN2 and TEN accounted for 0.89%, 1.01%, 0.47%, 0.37%, 0.52%, and 0.52% of total effects, respectively. This result indicated that there were indeed strong effects in the region that could affect phenotypic variation.

Significance analysis between phenotypes and SNPs of the QTL region

To further verify the correlation between egg numbers and 29 SNPs, we analyzed the SNPs’ change of dominant genotypes with breeding. We found that the frequency of dominant genotypes gradually increased with the selection of W2 for egg production, which indicated that the SNPs were also selected with the selection of phenotype, which could well explain these SNPs and egg production is significantly correlated (Fig. 6). Therefore, the 29 SNPs were identified as important candidate SNPs (Table 5).

Fig. 3 Manhattan and quantile–quantile plots of the GWAS for EN2 and TEN traits. Each dot represents a SNP in the dataset. The horizontal red and blue lines indicate the genome-wide significance (P-value = 6.43e-6) and suggestive thresholds (P-value = 1.29e?4), respectively. A, C Manhattan and quantile–quantile (Q-Q) plots of EN2 (egg number from 196 to 227 d); λ was the genomic inflation factor. B, D Manhattan and quantile–quantile (Q-Q) plots of TEN (egg number from onset of laying eggs to 354 d)

Table 4 Annotation of SNPs with significance associated with egg number traits in 627 of line W2 hens

Fig. 4 LD plots of the 29 significant SNPs from 3159 chickens in line W2. Box plots were structured using four SNPs selected from four blocks. The SNPs effect was tested using Kruskal-Wallis univariate ANOVA analyses for egg number EBV

Fig. 5 Manhattan plots of GGA Z association study for EN2 and TEN traits. Each dot represents a SNP in the dataset. The horizontal red and blue lines indicate the genome-wide significance (P-value = 2.14e-5) and suggestive thresholds (P-value = 4.26e? 4), respectively. A Manhattan plot of EN2 (egg number from 196 to 227 d); B Manhattan plot of TEN (egg number from onset of laying eggs to 354 d)

Fig. 6 Frequency variations of the 29 SNPs associated with higher egg number EBVs. G6, G7, and G8 represent the three generations of the line W2 chickens, and each line represents the differences in the frequency change curve of a SNP’ genotype with higher egg number EBV

Table 5 Annotation of the 29 SNPs that were significantly different in EBV of 3159 individuals

Discussion

In the twentieth century, layers and meat-type chickens were developed in different directions to avoid conflicts between egg-producing and meat-producing traits [18]. Such activities were extremely successful in improving productivity; laying hens produced more than 320 eggs during their 52-week laying period, while broilers gained 50–60 times their weight from hatching to sale [19]. Those activities also led to substantial genomic differences between layers and meat-type chickens. To explore the genetic mechanisms of laying traits in meat-type chickens, 11,279 chickens were used in this experiment.

The genetic contributions of egg number have low or medium levels of heritability, typically ranging from 0.16 to 0.64 [6]. However, among the seven lines evaluated in this study, the heritability of egg number was different (0.034–0.258), but it was not the color of their feathers that caused the difference. The range of heritability observed in this study was consistent with those reported in White Leghorn hens (0.05–0.27) [20]. The main reasons for other differences in heritability included environmental effects, population size, animal breeds, differences in assessment methods, and the reliability of the individual information records.

Further analysis the EN of partial periods and the total period in Line W2, Liu et al. predicted heritability range between 0.05–0.24 for White Leghorn hens [12]. Our results showed that AFE and EN1 had a high level (0.24–0.27) of heritability, with the heritability of the other stages being 0.093–0.17, which was slightly lower than those previously reported (0.17–0.20) [21]. AFE had a strong negative genetic correlation with EN1 (? 0.96), and a moderate negative correlation with TEN (? 0.37), which was consistent with the findings of Yuan et al. [21] and Liu et al. [12].

Meta-analyses have been widely used in multi-population GWAS studies. Falker-Gieske et al. used 1306 White Leghorn chickens and a meta-analysis of GWASs to identify the variations and genes related to feather pecking, and found 15 significant SNPs that may be related to feather pecking [22].

We used GWAS, meta-analyses, selective sweep analyses, and Bayesian analyses to scan the SNPs and QTLs associated with laying traits in meat-type chickens. A 10.81–13.05 Mb region on GGA Z was analyzed using these methods and 36 genes were annotated in this region. However, no studies have detailed the genes related to egg number traits in chickens. Notably, several regions associated with egg number traits have been identified on the GGA Z (https://www. anima lgeno me.org/cgi-bin/QTLdb/GG/index). The regions identified in our study were 1.9 Mb apart from those reported by Zhao et al., which were carried out in more than a dozen breeds, including White Leghorn and Rhode Island Red chickens [23]. Although the results were different between the present study and that of Zhao et al., they suggest that GGA Z is an important genetic region involved in egg production.

Referring to prior studies, we found that the regions between 10.81 Mb and 13.05 Mb on the GGA Z identified in this study were inconsistent with those reported in layers, which may be due to the differences in the genomes between layers and meat-type chickens. In the regions, genes that may be involved in reproductive traits were identified, including theGDNFgene,which has been reported to affect the in vitro maturation of human oocytes, increasing the prevalence of MII stage oocytes from 34% to 49% [24].GDNFhas also been reported to affect the development of oocytes in a variety of animals, such as humans, pigs and mice [25–27].

DAB2 is an endocytic adaptor protein in several NPXY motif-containing cell surface receptors, including the lipoprotein receptor and integrin. Budna et al. also suggested that theDAB2gene may be involve in the regulation of MII phase oocyte formation and other processes critical to porcine fertility [28]. ThePRKAA1gene encodes the catalytic α-subunit of 5′ AMP-activated protein kinase has an opposite effect--maintenance of the meiotic block in porcine and bovine oocytes [29, 30] AMPK can control steroidogenesis in ovarian cells (granulosa, theca cells and corpus luteum cells) and germ cell maturation [31]. TheNADK2gene in this region, the NAD kinase (NADK) is the sole NADP biosynthetic enzyme [32]. NADPH, which is the reduced form of NADP, regenerates cellular oxidative defense systems via glutathione reductase to counteract oxidative damage [33, 34]. Therefore, NADK2 may affect egg production by affecting reactive oxygen species levels [34]. TheWDR70gene was identified as a candidate gene for milk production and reproductive traits in Chinese and Northern European Holstein dairy cows [35]. The gene is homologous in chickens and cattle and may have the same function [36].

TheLIFRgene can form polymeric complexes with other receptors, such as glycoprotein 130 (GP130), to stimulate the JAK/STAT, MAPK and PI3K signaling pathways [37]. JAK/STAT and PI3K signaling pathways influence the hypothalamic-pituitary-ovarian (HPO) axis to regulate the breeding cycle of laying hens. Two genes,C6andC7,were identified in this region. The complement system plays an important role in the process of inflammation and infection, and also provides a link between innate and adaptive immunities. Immune capacity canaffect the reproductive performance, health, and welfare of hens [38].

In addition, three SNPs (rs318154184, rs13769886, rs313325646) associated with EN2, EN3, and TEN were located on or near the prolactin receptor (PRLR) gene on GGA Z. ThePRLRgene is associated with laying traits and consists of 15 exons and 14 introns in chicken [39, 40]. Liu et al. hypothesized thatPRLRmight be the primary gene responsible for sexual maturation and is considered a candidate genetic marker for reproductive traits [41]. Studies showed that the prolactin gene is expressed in the hypothalamus, pituitary, oviduct, and ovary of chickens. Some studies also showed that polymorphisms occur at different locations of thePRLRgene, including SNP in exons 2, 5 [40], and 10 [39], which were mainly associated with egg number traits.

Conclusion

In conclusion, our results indicated that the heritability of egg numbers was different in different chicken populations, and egg number was a low-heritability trait. Based on GWAS analyses of 627 individuals in one generation of line W2, meta-analyses of seven lines, and selective sweep analyses of lines W1 and W2, we found an interesting region on GGA Z. Where multiple SNPs and genes may be significantly related to egg production. Our findings contributed to a better understanding of the genetic basis of egg production.

Abbreviations

GGA:Gallus galluschromosome; GWAS: Genome-wide association study; LD: Linkage disequilibrium; QTL: Quantitative trait locus; Animal QTLdb: Animal QTL Database; SNP: Single nucleotide polymorphism; AFE: The age at first egg; EN: Egg numbers; EN1: The egg number from the onset of laying eggs to 195 days of age; EN2: The egg number from 196 to 227 days of age; EN3: The egg number from 228 to 307 days of age; EN4: The egg number from 308 to 354 days of age; TEN: The total egg number from first egg to 354 d; PCA: Principal component analysis; Q-Q: Quantile–quantile; Fst: The population differentiation index; GBLUP: Genomic best line prediction; ssGBLUP: Single-step GBLUP; GDNF: Glial cell derived neurotrophic factor; DAB2: DAB adaptor protein 2;

PRKAA1: Protein kinase AMP-activated catalytic subunit alpha 1; NADK2: NAD kinase 2, mitochondrial; WDR70: WD repeat domain 70; LIFR: Leukemia inhibitory factor receptor alpha; C6: Complement C6; C7: Complement C7; PRLR: Prolactin receptor; AMPK: Adenosine 5′-monophosphate (AMP)-activated protein kinase.

Supplementary Information

The online version contains supplementary material available at https://doi.org/10. 1186/s40104-022-00744-w.

Additional file 1: Fig. S1.Egg production curve of the complete laying period.

Additional file 2: Fig. S2.Line chart of cross validation error.

Additional file 3: Fig. S3.Admixture plot.

Additional file 4: Fig. S4.Manhattan and quantile–quantile (Q-Q) plots of the GWAS for egg number of 7 lines.

Additional file 5: Fig. S5.Manhattan and quantile–quantile (Q-Q) plots of the GWAS for AFE, EN1, EN3, and EN4 traits in line W2 chickens.

Additional file 6: Fig. S6.Manhattan plots of the estimated squaredmarker effect from GWAS, Bayes B and Bayes R.

Additional file 7: Table S1.Annotation of significant SNPs associated with the seven chicken lines.

Additional file 8.The phenotypic data used in this study.

Authors’ contributions

JQD contributed to data analysis and manuscript preparation. GMZ and YF contributed to data collection and manuscript revision. QHL, JZ, JW, RRL and MQZ contributed to the study design and revised the manuscript. DZ raised and sampled chickens, and collected the data. GPZ contributed to the design of the study, interpretation of the data, and manuscript preparation. All authors read and approved the final manuscript.

Funding

This work was supported by the Guangdong Provincial Key Research and Development Program (2020B02022002), the Shandong Agricultural Seed Improvement Project (2020LZGC013) and the Technological Innovation Project of Chinese Academy of Agricultural Sciences (CAAS-ZDRW202005).

Availability of data and materials

The whole genome sequencing data reported in this paper were deposited in the Figshare data center under accession number https://doi.org/10. 6084/m9. figsh are. 19178 906 that can be publicly accessed at https://figshare. com/.

Declarations

Ethics approval and consent to participate

Ethical approval on animal survival was provided by the Animal Welfare and Ethics Committee of the Institute of Animal Sciences (IAS), Chinese Academy of Agricultural Sciences (CAAS, Beijing, China), with the following reference number: IAS2019–44.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Author details

1Institute of Animal Sciences, Chinese Academy of Agricultural Sciences, Beijing, China.2MiLe Xinguang Agricultural and Animal Industrials Corporation, MiLe, China.

Received: 20 February 2022 Accepted: 9 June 2022