We performed single-trait GWAS, multi-trait GWAS, and LONG-GWAS for six body size traits on three growth stages in Simmental beef cattle. However, the three methods yielded different results with few shared loci. The reason for this discrepancy was likely due to the restricted dataset in single-trait GWAS and LONG-GWAS analyses. One universal phenomenon that cannot be ignored is that growth traits are controlled by multiple genes [29], and each method had its specific advantages in the identification of distinct loci. For example, single-trait GWAS is robust in detecting trait-specific QTLs and multi-trait GWAS is efficient for mapping pleiotropic QTLs [30], whereas LONG-GWAS can improve the detection power for time-dependent and consistent loci [31]. Thus, combining these three GWAS methods was expected to markedly improve the analysis the genetic mechanism of the body traits of beef cattle. In addition, since many complex traits have a similar architecture across diverse species [7], which prompted us to compare some of our significant genes with the previous reports that investigate the same genes and their involvements in growth. As a result, 21 suggestive genes were considered as candidate genes that were involved in the growth of cattle, swine, mice, and humans.
Candidate Genes Identified by Single‑trait GWAS
On BTA1, two SNPs were near SOX2 (SRY-Box 2), which encodes a transcription factor involved in the regulation of embryonic development [32, 33]. The paralog of this gene is SOX17, which positively affects the growth traits of cattle, and the conserved regions of this gene in human genome is closely related to body development [34]. On BTA2, two SNPs were near SNRPD1 (small nuclear ribonucleoprotein D1 polypeptide), which is a member of the ghrelin receptor family, and the encoded protein is involved in zinc-dependent signaling in epithelial tissue [35]. On BTA6, variations near RASGEF1B (RasGEF domain family member 1B) were associated with body height [36]. In addition, body height was positively correlated with calcium absorption, which is an important determinant of calcium balance [37]. On BTA7, a SNP near EFNA5 (ephrin A5) was associated with two traits (HS and AS) at the same stage. It was also identified as a candidate gene for growth traits in broiler chicken [38]. Another SNP near PTBP1 (polypyrimidine tract binding protein 1) was found to show genome-wide association with growth traits at 6 months by both single-trait and multi-trait GWAS. Its expression level determined the release of insulin, thereby affecting development [39]. On BTA9, a SNP (BovineHD0900027283) located in SNX9 (sorting nexin 9), as olfactory receptor, was associated with growth traits in Yorkshire pigs [40]. The SNP near SV2C (synaptic vesicle glycoprotein 2C) was associated with BH by both single-trait and multi-trait GWAS. This gene was reported to modulate dopamine release in neural and endocrine cells [41]. On BTA11, PKDCC (protein kinase domain containing, cytoplasmic) was associated with HS in both single-trait and multi-trait analyses. This gene was involved in the maintenance of bone density in humans [42]. On BTA13, SYNDIG1 (synapse differentiation inducing 1) has been reported as a factor influencing the final weight and backfat thickness of Landrace pigs [43], whereas the AKR1E2 (aldo-keto reductase family 1 member E2) variant was associated with body length and girth in cattle [44]. On BTA23, the PRIM2 (DNA primase subunit 2) was associated with body weight and trait changes in pigs [45, 46].
Single‑trait GWAS versus multi‑trait GWAS
Multiple-trait analysis of linkage experiments has been reported to significantly enhance the power to detect common SNPs across traits [47, 48]. Therefore, we used multi-trait analysis to complement single-trait GWAS rather than to replace it. In the single-trait GWAS, the minimum p values for the three stages were 3.90E-07, 9.92E-07 and 2.11E-08 respectively. These three values decreased to 8.23E-13, 1.73E-09 and 3.84E-10 in the multi-trait GWAS, respectively. We also identified several critical loci as follows. On BTA1, the SLC37A1 (solute carrier family 37, member A1) gene, which encodes a glucose-6-phosphate transporter that is involved in the homeostasis of blood glucose [49], was found to be the best candidate gene for modifying milk production traits [50]. On BTA6, LAP3 (leucine aminopeptidase 3) was reported to play critical roles in the regulation of hormone levels and protein maturation. Another study demonstrated that putative regulatory elements in the PCDH7 (protocadherin 7) gene may have roles in residual feed intake in Nelore cattle [51]. In addition, MANEA (mannosidase endo-alpha), which has roles in proteolysis, was associated with the birth weight of Canchim beef cattle [52]. On BTA11, a mutation in the LHCGR (luteinizing hormone/choriogonadotropin receptor) gene was as the cause of empty follicle syndrome [53].
Single‑trait GWAS versus LONG-GWAS
We used LONG-GWAS, which involved multiple phenotype measurements for each individual [24]. One disadvantage of this method was that incorporating all data may have overwhelmed significant signal, that is, if QTL effects varied during the different stages [54]. In this study, these time-specific expressed QTLs identified by the single-trait and multi-trait GWAS were not detected by LONG-GWAS. However LONG-GWAS also detected some significant functional loci as follows. On BTA1, P2RY1 (purinergic receptor P2Y1), a candidate gene that affects the serum Ca2+, encoded for a member of the family of G protein-coupled receptor family that works as receptor for extracellular ATP and ADP [55]. On BTA3, the MPZL1 (myelin protein zero like 1) gene could significantly enhance the migratory and metastatic potential of hepatocellular carcinoma cells by phosphorylating and activating the pro-metastatic protein [56]. Besides on BTA8, LINGO2 (leucine rich repeat and Ig domain containing 2), which is expressed in the central nervous system of mouse embryos, has been reported to associate with the body mass in a cohort of elderly Swedes [57]. On BTA18, CMIP (C-Maf inducing protein), a candidate gene for reading-related traits, was also associated with plasma lipoprotein levels [58]. Moreover, WSCD1 (WSC domain containing 1), which encodes a protein with sulfotransferase activity that participates in the metabolism of glucose, was a candidate gene for feed efficiency and feeding behaviors in the White Duroc × Erhualian F2 population [59].