Genome-wide association study reveals novel candidate genes for litter size in Markhoz goats

Introduction: The Markhoz goat is the only breed that can produce high-quality fiber called mohair in Iran; however, the size of its population has faced a dramatic decline during the last decades, mainly due to the reluctance of farmers to rear Markhoz goats caused by a reduction in goat production income. Litter size at birth (LSB) and weaning (LSW) are two economically important reproductive traits for local goat breeders and have the potential of increasing the population growth rate. The present study was aimed to identify possible genomic regions that are associated with LSB and LSW in Markhoz goats using a genome-wide association study (GWAS).

Methods: To this end, 136 Markhoz goats with record(s) of kidding were selected for GWAS using the Illumina Caprine 50K bead chip. The individual breeding values (BV) of available LSB and LSW records estimated under an animal mixed model were used as the dependent variable in the GWAS, thereby incorporating repeated categorical variables of litter size.

Results: Four SNPs on chromosomes 2, 20 and 21 were identified to be significantly associated (FDR p<0.05) with LSB after multiple testing correction under a Bayesian-information and Linkage-disequilibrium Iteratively Nested Keyway (BLINK) model. Least-square analysis was performed to investigate the effects of detected genotypes on LSB. Ultimately, the GWAS results introduced six candidate genes, including GABRA5, AKAP13, SV2B, PPP1R1C, SSFA2 and TRNAS-GCU in a 100 kb adjacent region of the identified SNPs. Previous studies proposed functional roles of GABRA5 and AKAP13 genes in reproductive processes; however, the role of other candidate genes in reproduction is not clear.

Conclusion: These findings warrant further investigation for use in marker-assisted selection programs in Markhoz goats.

Introduction

Throughout history, goats have been a primary production species for mankind due to their ability to withstand variable and harsh environmental conditions and their desirable production of meat, milk, fiber, and skin. Improving economic traits in goats, such as prolificacy and viability, can be very profitable for rural people in countries with low-quality grazing lands, where goat farming is one of the primary sources of income.

Markhoz goats are the only mohair producing breed in Iran with different coat colors such as white, black, and various shades of brown (1). The population size of the Markhoz goat underwent a considerable decrease during the last two decades (2), such that only about 1,000 head remain in their main native habitat in western and northwestern Iran. The main reasons for the decreased population of the indigenous Markhoz goat include changes in the management system for how animals are reared and the generally low income gained by goat farming in local regions. Hence, increasing the number of kids born per kidding and subsequently increasing the total income through the selling of kids, fiber, and meat may encourage ranchers toward goat production.

Litter size at birth (LSB) and litter size at weaning (LSW) are two reproductive traits that are known to be controlled by several underlying genes in goats (3). Specifically, genetic variants with significant effects on LSB have been identified in genes including, GDF9, BMP15, GnRH1, KISS1, KITLG, NGF, POU1F1, PRLR and promoter of miR-9 gene in various breeds of goats (4–11). However, there may be other genes that affect LSB and LSW that have remained unknown. Nowadays, with the development and availability of SNP genotyping technologies, conducting genome-wide association studies (GWASs) and detecting candidate genes and genetic variants that may have a significant association with economic traits, have become much easier and faster. In this regard, many GWASs have been conducted using the Illumina Caprine 50K beadchip [Illumina Inc., San Diego, CA (12)] for various economically important traits in different goat breeds including coat color and mohair traits (13), body morphological traits (14), conformation and milk yield (15), and resistance to nematodes (16). Similarly, multiple researchers have conducted GWAS on liter size at birth in sheep (17–19). However, there is only one GWAS for the number of kids alive per kidding in goat (20) and no GWAS for the number of kids alive till weaning in sheep or goat.

Detection of significantly associated SNPs with LSB and LSW in goat would lead to enhanced efficacy of animal selection in breeding strategies by reducing the cost and time required to raise and phenotypically characterize animals as they mature. As complex traits, LSB and LSW are known to be controlled by many genes with variants having a small effect. Hence, the possibility of incorporating genetic variants into selection strategies will likely accelerate the rate of improvement of such reproductive traits as compared to traditional phenotypic selection. The purpose of this study was to conduct GWAS to detect possible genomic regions and variants associated with LSB and LSW in Markhoz goats, with the potential of applying results in genomic selection.

Materials and methods

Animals and phenotypes

All female goats (n = 184) existing at the Markhoz Goat Performance Testing Station in Sanandaj, Kurdistan, Iran, were selected for inclusion in the study. The herd is reared under a semi-intensive management system, in which animals graze on natural pastures from spring to early autumn and fed a diet consisting of alfalfa and wheat straw for the rest of the year. At the age of 16–18 months, does are mated for the first time. The kidding season starts in late winter and ends in early spring. Litter size at birth (LSB) and litter size at weaning (LSW) were the two reproductive traits evaluated. LSB described the number of live kids born to the doe. LSW described the number of live kids at weaning for each doe, typically evaluating kids at 22–27 weeks old. LSB and LSW were categorical variables potentially repeated per doe as they aged and had subsequent litters.

Statistical analyses

Prediction of breeding values

Predicted breeding values (PBV) for LSB and LSW were generated for use in the GWAS to capture the repeated categorical values of LSB and LSW. A total of 3,410 litter size records for the Markhoz goats were collected from 1994 to 2019 for use in predicting breeding values. Accuracy of EBVs is based on the amount of performance information available on the animal and its close relatives. Selection using EBVs is more accurate, especially for low heritable traits like litter size, which allows for faster genetic gain compared to mass selection using phenotypes. Breeding values have the advantage that they are free of systematic environmental effects on measured phenotypes, as these effects are considered in the statistical model used for the estimation of EBVs. Additionally, they reflect the genetic makeup more accurately because they do not solely rely on own records but include information from all measured relatives. The pedigree file included 5,396 animals with 1,533 dams and 252 sires. The number of founders, individuals with progeny, and individuals without progeny were 343, 1,785 and 3,611, respectively. Breeding values for each individual was estimated, applying a repeatability threshold animal model using ASReml 2.0 (21) as follows:

$\begin{array}{l}y=Xb+Za+Wpe+e\end{array}$

where y is a vector of phenotypic value for LSB/LSW, b, a, pe and e are vectors of fixed effects including year of kidding (2010–2019), age of dam (2–9 years) and parity (1–7), random animal effects, random permanent environmental effects and random residual effects, respectively. X, Z, and W are design matrices that relate records to fixed, animal and permanent environmental effects, respectively.

Genotyping and quality control

All animal procedures were approved by the Cornell University Institutional Animal Care and Use Committee prior to sampling (protocol #2014-0121). Vacutainer tubes containing K₂EDTA as an anticoagulant were used to collect whole blood (5 ml) samples from the jugular vein of goats. Samples were immediately stored at −20°C until DNA extraction. A standard Phenol-Chloroform DNA extraction method was utilized for extracting genomic DNA (22). The Illumina Caprine 50K beadchip (Illumina, Inc., San Diego, CA, United States), including 53,353 SNPs, was used for genotyping samples (VHL Genetics, Wageningen, Netherlands). Golden Helix SVS v8.3.4 (Golden Helix, Bozeman, MT, United States) software was used for quality control process as follows: (1) 624 SNPs were removed for a call rate less than 0.9; (2) 2,540 SNPs with a minor allele frequency less than 0.03 were excluded; (3) 810 SNPs were not assigned to a genomic location; thus they were removed; (4) five samples were removed for a genotyping call rate less than 0.9; and (5) three samples with an estimated identity-by-state (IBS) score greater than 0.9 were removed to eliminate the possible effects of substantial relatedness between individuals on the overall results. Furthermore, 39 samples were excluded because they had no history of kidding in the data set. After the quality control process, 49,764 SNPs and 136 animals remained for GWAS.

Genome-wide association studies

The GAPIT v3.0 R package was used for investigating the association between genomic regions and PBVs for LSB and LSW as phenotypes. Several models including General Linear Model (GLM), Mixed Linear Model (MLM), Multi Locus Mixed Model (MLMM), Compressed Mixed Linear Model (CMLM), Enriched Compressed Mixed Linear Model (ECMLM), Factored spectrally transformed Linear Mixed Model (Fast-LMM), Settlement of MLM Under Progressively Exclusive Relationship (SUPER), Fixed and random model Circulating Probability Unification (FarmCPU), Efficient Mixed-Model Association (EMMA), Efficient Mixed Model Association eXpedited (EMMAX) and BLINK were performed to find the best model for fitting PBV LSB and LSW data (Figure 2). The BLINK (23) and FarmCPU (24) models showed less deviation from expectation in the Q-Q plots than other models. Despite FarmCPU identifying more significantly associated SNPs in the model comparison, BLINK was selected for the final GWAS because it controlled both false positives and false negatives effectively, showing a sharp upward deviated tail and a straight line close to the 1:1 line. In the BLINK method, markers in LD (r² > 0.7) with the most significant marker are excluded from the analysis and the maximum likelihood of a random effect model is approximated by using the Bayesian Information Content (BIC) of a fixed-effects model. By applying that, the most significant markers will be selected among all markers that remained after LD exclusion and then used as cofactors in the model to test all markers across the genome. Because of the abovementioned reasons, the BLINK method had improved statistical power and better control on false positives in comparison with kinship-based methods (23). In the present study, the first 10 PCs explained about 23% of stratification, of which 12% was explained by the first three PCs. For considering population structure and avoiding biases due to population stratification in the present GWAS, only the first three PCs were included as covariates in the model, because no difference observed in the results when 4–10 top PCs (explained about 11% of stratification) were included in the analyses.

Significance of marker association was determined using a false discovery rate (FDR) adjusted p-value of less than 0.05. The ggplot R package was used to generate Q–Q and Manhattan plots (25). The genomic heritability was estimated using LDAK software v5.1 (26).

Least-square analysis and correlations between studied traits

Least-square analyses were conducted to investigate the negative/positive effects of detected genotypes on studied traits using PROC MIXED of SAS v8.2 software (27). Phenotypic and genotypic correlations between LSB and LSW were estimated using CORR procedure of SAS v8.2 and GCTA software tool (28), respectively.

Identifying candidate genes

A region with a distance of 100 kbp up-stream and 100 kbp down-stream of significant SNPs was explored to detect the nearest gene using Capra hircus ARS1 assembly (29) in the NCBI database. The LD decay pattern of the Markhoz goat population was estimated and the explored distance was selected based on the intersection point between the LD line and the r² threshold determined the LD decay value (r² < 0.1). Finally, gene databases, such as Genecards, National Center for Biotechnology Information (NCBI), Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) were scrutinized to find the functions and pathway of identified genes.

Results

Descriptive statistics and predicted breeding values

The descriptive statistics and predicted breeding values for LSB and LSW traits are presented in Table 1. The studied population had an average LSB of 1.16, while only 0.99 of them remained alive at the age of weaning (LSW). The mean PBV for LSB and LSW was 0.0031 and 0.0073, respectively. Furthermore, the distribution of studied traits is depicted in Figure 1.

TABLE 1

Table 1. Descriptive statistics and estimated breeding values for litter size at birth and litter size at weaning in the Markhoz goat.

FIGURE 1

Figure 1. The distribution of litter size at birth (A); litter size at weaning (B); PBVs for litter size at birth (C); and PBVs for litter size at weaning (D).

Genome-wide association studies

The Manhattan and Q-Q plots for PBV of LSB and LSW are depicted in Figure 2. Furthermore, the Q-Q plot for litter size at birth fitting other models is indicated in Figure 3. Four significantly associated SNPs on chromosome 2 (rs268267345, unadjusted p-value = 1.35e−07), chromosome 20 (rs268258357, unadjusted p-value = 1.33 e−08) and chromosome 21 (rs268288690, unadjusted p-value = 5.41e−13; rs268256209, unadjusted p-value = 2.56e−07) were identified for PBV LSB (Figure 2A). The FDR adjusted p-values for the identified SNPs were 0.002, 0.0003, 2.69e−08 and 0.003, respectively. The Q–Q plot for LSB showed a very slight genomic inflation factor with λ_GC = 1.01. In contrast, no significantly associated SNPs were found for PBV LSW after multiple testing correction under the various models tested for this trait. The detailed information, including location, alleles, and p-values, for significantly associated SNPs on the PBV LSB trait is provided in Table 2. Considering that the distance between the significantly associated SNPs from their neighboring SNPs was farther than 500 kb and the BLINK model uses a minimum distance of 300 kb for exclusion of markers, linkage disequilibrium analysis was not performed.

FIGURE 2

Figure 2. The Manhattan and quantile–quantile (Q–Q) plots of -log 10 (p-value) for litter size at birth (LSB) (A); and litter size at weaning (LSW) (B) in the Markhoz goats using the BLINK GWAS model. The green horizontal line indicates the FDR adjusted p-value threshold of 0.05.

FIGURE 3

Figure 3. The Q–Q plot for litter size at birth (LSB) fitting Bayesian-information and Linkage-disequilibrium iteratively nested keyway (BLINK) (A); general linear model (GLM) (B); mixed linear model (MLM) and multi locus mixed model (MLMM) (C); compressed mixed linear model (CMLM) (D); enriched compressed mixed linear model (ECMLM) (E), factored spectrally transformed linear mixed model (Fast-LMM) (F); settlement of MLM under progressively exclusive relationship (SUPER) (G); fixed and random model circulating probability unification (FarmCPU) (H); and efficient mixed-model association (EMMA) and efficient mixed model association eXpedited (EMMAX) (I) GWAS models.

TABLE 2

Table 2. Genome-wide association study identifies four SNPs significantly associated with the predicted breeding value of litter size at birth (LSB) in Markhoz goats.

Least-square analysis of identified genotypes for LSB and phenotypic and genotypic correlations

The results of the least-square analyses for LSB are provided in Table 3. Results show that markers rs268267345, rs268258357, and rs268288690 in the genome of the Markhoz goats significantly resulted in increased litter size in goats having one or two mutated alleles. In contrast, rs268256209 SNP has a significant negative effect on LSB in goats that carry two G alleles.

TABLE 3

Table 3. Genotypic frequency and least-square mean ± standard error of litter size at birth for identified SNPs in the Markhoz goats.

The phenotypic correlation between LSB and LSW was 0.697 (P < 0.01). Furthermore, analysis of the genotypic correlation between these traits showed a strong genetic relationship of 0.725 (P < 0.01).

Estimated genomic heritability

The estimated heritability using recorded data and pedigree for LSB and LSW of the entire population was 0.018 and 0.32e−6, respectively. However, the genomic heritability for studied traits was 0.011 and 0.21e−7, respectively. The detailed information for variance components of LSB and LSW traits is presented in Table 4.

TABLE 4

Table 4. Estimation of variance components and genetic parameters for litter size at birth (LSB) and litter size at weaning (LSW) in the Markhoz goats.

Detected candidate genes

Six candidate genes were found within the 100 kb windows up- and down-stream of the identified SNPs located on the genome of Capra hircus in the NCBI database (ARS1 assembly, accession number: GCF_001704415.1). The detected genes harboring SNPs, their distance from SNPs and the roles of each gene are presented in Table 5.

TABLE 5

Table 5. Genes within 100 kb distance from identified significant SNPs and their description.

Discussion

In Iran, the Markhoz goat is the only breed that can produce mohair in black, gray, white and varying shades of brown color. The fiber obtained from the Markhoz goat is both culturally and socially important for Kurdish people costumes (30). Thus, the decreasing population size of this valuable breed could be mitigated by using breeding programs incorporating genetic markers affecting reproductive traits such as litter size at birth (LSB). This is the first study on litter size traits in the Iranian goats at a genome-wide scale which can be employed as a practical tool to identify novel genetic markers that may influence litter size in the Markhoz and other goat breeds.

For both traits, the first three principal components (PCs) were included in the model to correct for population structure. A total of four SNPs, found on chromosomes 2, 20, and 21, were significantly associated with the PBV of LSB (Figure 2A). The Q–Q plots showed that observed versus expected data was well aligned, indicating minimal population stratification affecting the model (λ_GC = 1.01). In contrast, the results of the GWAS failed to reveal significant association with the PBV of LSW (Figure 2B).

The most significantly associated variant (rs268288690) with the PBV of LSB is located within the GABRA5 gene. GABRA5 encodes a 462 amino acid protein of the GABA receptors family known as GABA-A. It has been suggested that GABA acts at GABA-A receptors in the central and peripheral nervous systems as the major inhibitory neurotransmitter (31). Watanabe et al. (32) investigated the role of GABA in the regulation of GnRH neurons. They reported that in the course of fetal development, GABA is involved in the regulation of GnRH neuron migration from the olfactory placodes into the forebrain. They also stated that negative and positive feedback of estradiol are mediated by GABA, and there is a significant correlation between these feedbacks and frequency of GABA transmission to GnRH neurons. In a recent study, Di Giorgio et al. (33) demonstrated that GABA-A and GABABRs interact with kisspeptin (a protein encoded by the Kiss1 gene) in the regulation of reproductive processes. For instance, GABA increases Kiss1 expression by affecting GABA-A receptors in early embryo development. In addition, at the time of ovulation in adults, a main double excitatory input to (GnRH) neurons is provided by the AVPV/PeN neuron population, leading to the expression of GABA and kisspeptin.

Another significantly associated SNP identified in the present study (rs268256209) is located 99,644 bp upstream of the AKAP13 gene. AKAP13 is a member of the AKAP family which is a structurally diverse protein and is involved in the binding process to the regulatory subunit of protein kinase A (PKA) and confining the holoenzyme to discrete locations within the cell. Luconi et al. (34) reported that AKAP proteins are expressed in both female and male reproductive systems, especially during gametogenesis. It has been suggested that AKAP–PKA interactions control the maturation of oocytes (35).

Based on the functional role of PPP1R1C, SV2B, TRNAS-GCU and SSFA2 genes, there is no evidence to suspect a causative association with the LSB trait in our study. Whereas, both GABRA5 and AKAP13, are the more likely genes potentially influencing LSB based on their influence of reproductive processes.

The results of least-square analyses showed that rs268288690 SNP leads to a significant increase in litter size in Markhoz goats (p < 0.01) so that goats having two copies of the mutated allele (GG genotype) had more kids within the litter than those carrying one or no copy of the mutated allele (AG and AA genotypes). The GG genotype of rs268258357 SNP indicated the highest litter size (1.822 ± 0.161) among all identified genotypes. Similarly, rs268267345 SNP also had a positive impact on litter size. These findings revealed that only rs268256209 SNP negatively affected litter size, while the other three SNPs identified from GWAS positively influenced the number of kids in Markhoz goats.

Estimated heritability for LSB and LSW was 0.018 and 0.32e−6, respectively. These values are lower than the estimated genomic heritability of 0.011 for LSB and 0.21e−7 for LSW. One of the main reasons for the observed differences could be different sample sizes used for estimating heritability via BLUP and LDAK models. In BLUP, we used all available data (3,410) and pedigree (5,396) records, while only 136 genotyped individuals were used in the LDAK method to predict genomic heritability. Besides, much of the heritability of traits may not be accounted for by rare, low-frequency genetic variants, known as the missing heritability problem (36).

There are many genes that have been identified as associated with litter size in goat and sheep. Among them, GDF9, BMP15, BMPR1B, and IGF1 genes have been widely discussed in literature acknowledging their effects on litter size. However, to the best of our knowledge, the associations between candidate genes detected in the present study and litter size have not been reported previously. Therefore, according to their vital functions in reproductive processes, the GABRA5 and AKAP13 genes could be important novel candidate genes for litter size in small ruminants. However, our study had some limitations, including small sample size and relatively low heritability of studied traits. Thus, more genotyped animals are required to validate the impact of these potential candidate genes on litter size in goats.

Additionally, due to the high genomic correlation between LSB and LSW traits, we expected to find some common regions for the two traits, but the GWAS for the PBV for LSW failed to detect any significantly associated genomic regions. One possible reason could be that LSW is affected by environmental conditions or the pattern of effects may be altered by environment and genetic interactions. Furthermore, LSW may be influenced by rare causal variants that are not included in the current goat medium-density SNP-chip and not captured by linkage disequilibrium. Additionally, the complexity of gene interactions on this trait, such as epistasis, was not considered in this study and may play a more prominent role in regulation. The next possible reason could be the low frequency of the causal variants, which require a larger sample size to capture effect.

Despite the significant markers found for LSB, the present study has some limitations regarding LSW trait. It should be noted that LSB and LSW had extremely low heritabilities, suggesting a low possibility to achieve rapid genetic progress through phenotypic selection for LSB. Furthermore, LSW is generally connected to the mothering ability of the doe and environmental factors such as farm management. In addition, the sample size used in this study was limited, due to the low population size of Markhoz goats. Thus, caution must be taken when interpreting the results of the present GWAS, especially for LSW trait.

Conclusion

To conclude, we found plausible candidate genes based on SNPs associated with the EBV of LSB in the Markhoz goats using the 50K Caprine SNP-chip for the first time. The significant SNPs and genes identified in the present study can be beneficial for future molecular-based breeding for increased litter size at birth in goats; however, due to the low sample size used in this study, the results should be interpreted with caution. It is noteworthy that a breeding program focused on the major variations for LSB would not necessarily increase the number of surviving progenies due to the extremely low heritability of LSW. There may be no net impact on LSW from the slight increase in litter size caused by the substantial variations due to reduced viability.

Data availability statement

The SNP genotype data and EBVs for LSB and LSW are available in the Zenodo repository (https://zenodo.org/record/5824843).

Ethics statement

The animal study was reviewed and approved by Cornell University Institutional Animal Care and Use Committee (protocol #2014-0121). Written informed consent was obtained from the owners for the participation of their animals in this study.

Author contributions

PM performed data curation, formal analysis, visualization, methodology, and writing—original draft. AR and JR managed the project and contributed to writing— review and editing manuscript. AN-G conducted data curation, interpreted the results and contributed to review and editing manuscript. MR was responsible for data visualization. HH contributed to review and editing manuscript and provided SNP data and financial support. All authors contributed to the article and approved the submitted version.

Funding

Funding for the Markhoz goats' genotyping and part of the analysis has been supported by the laboratory of HH at Cornell University.

Acknowledgments

We would like to thank the Markhoz Goat Performance Testing Station Management for their cooperation by providing data used in this study.

Conflict of interest

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

Publisher's note

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

References

1. Mahmoudi P, Rashidi A, Razmkabir M. Inbreeding effects on some reproductive traits in Markhoz goats. Anim Prod Sci. (2018) 58:2178–83. doi: 10.1071/AN17043

CrossRef Full Text | Google Scholar

2. Bahmani HR, Tahmoorespur M, Aslaminejad AS, Vatankhah M, Rashidi A. Simulating past dynamics and assessing current status of Markhoz Goat population on its habitat. Iran J Appl Anim Sci. (2015) 5:347–54.

Google Scholar

3. Ahlawat S, Sharma R, Maitra A, Tantia MS. Current status of molecular genetics research of goat fecundity. Small Rumin Res. (2015) 125:34–42. doi: 10.1016/j.smallrumres.2015.01.027

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Feng T, Geng CX, Lang XZ, Chu MX, Cao GL, Di R, et al. Polymorphisms of caprine GDF9 gene and their association with litter size in Jining Grey goats. Mol Biol Rep. (2011) 38:5189–97. doi: 10.1007/s11033-010-0669-y

PubMed Abstract | CrossRef Full Text | Google Scholar

5. An XP, Hou JX, Zhao HB, Li G, Bai L, Peng JY, et al. Polymorphism identification in goat GNRH1 and GDF9 genes and their association analysis with litter size. Anim Genet. (2013) 44:234–8. doi: 10.1111/j.1365-2052.2012.02394.x

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Maitra A, Sharma R, Ahlawat S, Tantia MS, Roy M, Prakash V, et al. Association analysis of polymorphisms in caprine KiSS1 gene with reproductive traits. Anim Reprod Sci. (2014) 151:71–7. doi: 10.1016/j.anireprosci.2014.09.013

PubMed Abstract | CrossRef Full Text | Google Scholar

7. An, X. P., Hou, J. X., Gao, T. Y., Lei, Y. N., Song, Y. X., Wang, J. G., et al. Cao, and B.Y. (2015). Association analysis between variants in KITLG gene and litter size in goats. Gene. 558, 126–130. doi: 10.1016/j.gene.2014.12.058

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Hou JX, An XP, Han P, Peng JY, Cao BY. Two missense mutations in exon 9 of caprine PRLR gene were associated with litter size. Anim Genet. (2015) 46:87–90. doi: 10.1111/age.12223

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Ghoreishi H, Fathi-Yosefabad S, Shayegh J, Barzegari A. Identification of mutations in BMP15 and GDF9 genes associated with prolificacy of Markhoz goats. Arch Anim Breed. (2019) 62:565–70. doi: 10.5194/aab-62-565-2019

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Mahmoudi P, Rashidi A, Rostamzadeh J, Razmkabir M. Association between c.1189G>A single nucleotide polymorphism of GDF9 gene and litter size in goats: a meta-analysis. Anim Reprod Sci. (2019) 209, 106140. doi: 10.1016/j.anireprosci.2019.106140

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Mahmoudi P, Rashidi A, Rostamzadeh J, Razmkabir M. A novel variant in the promoter region of miR-9 gene strongly affects litter size in Markhoz goats. Theriogenology. (2020) 158:50–7. doi: 10.1016/j.theriogenology.2020.09.008

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Tosser-Klopp G, Bardou P, Bouchez O, Cabau C, Crooijmans R, Dong Y, et al. Design and characterization of a 52K SNP chip for goats. PLoS ONE. (2014) 9:e86227. doi: 10.1371/journal.pone.0086227

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Nazari-Ghadikolaei A, Mehrabani-Yeganeh H, Miarei-Aashtiani R, Staiger EA, Rashidi A, Huson HJ, et al. Genome-wide association studies identify candidate genes for coat color and mohair traits in the Iranian Markhoz goat. Front Genet. (2018) 9:105. doi: 10.3389/fgene.2018.00105

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Rahmatalla SA, Arends D, Reissmann M, Wimmers K, Reyer H, Brockmann GA, et al. Genome-wide association study of body morphological traits in Sudanese goats. Anim Genet. (2018) 49:478–82. doi: 10.1111/age.12686

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Mucha S, Mrode R, Coffey M, Kizilaslan M, Desire S, Conington J, et al. Genome-wide association study of conformation and milk yield in mixed-breed dairy goats. J Dairy Sci. (2018) 101:2213–25. doi: 10.3168/jds.2017-12919

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Silva FF, Bambou JC, Oliveira JA, Barbier C, Fleury J, Machado T, et al. Genome wide association study reveals new candidate genes for resistance to nematodes in Creole goat. Small Rumin Res. (2018) 166:109–14. doi: 10.1016/j.smallrumres.2018.06.004

CrossRef Full Text | Google Scholar

17. Gholizadeh M, Rahimi-Mianji G, Nejati-Javaremi A, Koning DJD, Jonas E. Genomewide association study to detect QTL for twinning rate in Baluchi sheep. J Genet. (2014) 93:489–93. doi: 10.1007/s12041-014-0372-1

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Xu S, Gao L, Xie X, Ren Y, Shen Z, Wang F, et al. Genome-wide association analyses highlight the potential for different genetic mechanisms for litter size among sheep breeds. Front Genet. (2018) 9:118. doi: 10.3389/fgene.2018.00118

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Hernández-Montiel W, Martínez-Núñez MA, Ramón-Ugalde JP, Román-Ponce SI, Calderón-Chagoya R, Zamora-Bustillos R, et al. Genome-wide association study reveals candidate genes for litter size traits in Pelibuey sheep. Animals. (2020) 10:434. doi: 10.3390/ani10030434

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Islam R, Liu X, Gebreselassie G, Abied A, Ma Q, Ma Y, et al. Genome-wide association analysis reveals the genetic locus for high reproduction trait in Chinese Arbas Cashmere goat. Genes Genom. (2020) 42:893–9. doi: 10.1007/s13258-020-00937-5

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Gilmour AR, Gogel BJ, Cullis BR, Thompson R. ASReml User Guide Release 2, 0. Hemel Hempstead: VSN International Ltd (2006).

Google Scholar

22. Sambrook J, Russell DW. Purification of nucleic acids by extraction with phenol:chloroform. CSH Protoc. (2006). doi: 10.1101/pdb.prot4455

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Huang M, Liu X, Zhou Y, Summers RM, Zhang Z. BLINK: a package for the next level of genome-wide association studies with both individuals and markers in the millions. Gigasci Giy. (2019) 154. doi: 10.1093/gigascience/giy154

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Liu X, Huang M, Fan B, Buckler ES, Zhang Z. Iterative usage of fixed and random effect models for powerful and efficient genome-wide association studies. PLoS Genet. (2016). 12:e1005767. doi: 10.1371/journal.pgen.1005767

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Wickham H. ggplot2: Elegant Graphics for Data Analysis. New York, NY: Springer-Verlag. (2016). doi: 10.1007/978-3-319-24277-4_9

CrossRef Full Text | Google Scholar

26. Speed D, Balding DJ. SumHer better estimates the SNP heritability of complex traits from summary statistics. Nat Genet. (2019) 51:277–84. doi: 10.1038/s41588-018-0279-5

PubMed Abstract | CrossRef Full Text | Google Scholar

27. SAS Institute. Users Guide, Version 8, 2. Statistics. Cary, NC: SAS Institute (2001).

Google Scholar

28. Yang J, Lee SH, Goddard ME, Visscher PM. GCTA: a tool for genome-wide complex trait analysis. Am J Hum Genet. (2011) 88:76–82. doi: 10.1016/j.ajhg.2010.11.011

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Bickhart D, Rosen B, Koren S. Single-molecule sequencing and chromatin conformation capture enable de novo reference assembly of the domestic goat genome. Nat Genet. (2017) 49:643–50. doi: 10.1038/ng.3802

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Rashidi A, Bishop SC, Matika O. Genetic parameter estimates for pre-weaning performance and reproduction traits in Markhoz goats. Small Rumin Res. (2011) 100:100–6. doi: 10.1016/j.smallrumres.2011.05.013

CrossRef Full Text | Google Scholar

31. Lee C, de Silva AJ. Interaction of neuromuscular blocking effects of neomycin and polymyxin B. Anesthesiology. (1979) 50:218–20. doi: 10.1097/00000542-197903000-00010

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Watanabe M, Fukuda A, Nabekura J. The role of GABA in the regulation of GnRH neurons. Front Neurosci. (2014) 8:387. doi: 10.3389/fnins.2014.00387

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Giorgio DI, Bizzozzero-Hiriart NP, Libertun MC, Lux-Lantos V. Unraveling the connection between GABA and kisspeptin in the control of reproduction. Reproduction. (2019) 157, 1741–7899. doi: 10.1530/REP-18-0527

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Luconi M, Cantini G, Baldi E, Forti G. Role of a-kinase anchoring proteins (AKAPs) in reproduction. Front Biol. (2011) 16:1315–30. doi: 10.2741/3791

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Newhall KJ, Criniti AR, Cheah CS, Smith KC, Kafer KE, Burkart AD, et al. Dynamic anchoring of PKA is essential during oocyte maturation. Curr Biol. (2006) 16:321–7. doi: 10.1016/j.cub.2005.12.031

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Young AI. Solving the missing heritability problem. PLoS Genet. (2019) 15:e1008222. doi: 10.1371/journal.pgen.1008222

PubMed Abstract | CrossRef Full Text | Google Scholar