Genetic Foundations of Nellore Traits: A Gene Prioritization and Functional Analyses of Genome-Wide Association Study Results

1. Introduction

Zebu cattle, also known as Bos taurus indicus, are prominently recognized by their pronounced hump, dewlap, heat tolerance, and tick resistance [1]. Over 8 000 years ago indicine cattle were known to have been domesticated in the Indus Valley [2]. Historically, Bos taurus indicus have been mainly selected for their ability to grow under harsh tropical environments and their capacity to maximize low-quality pasture; however, this resulted in a later age at puberty for Nellore cattle [3,4]. Nellore cattle are relevant in the global beef market and different breeding programs have been put in place to focus on genomic selection for the Nellore breed, to accelerate the rate of genetic progress, reduce generational interval, and increase the accuracy of selection [5,6,7].

Genome-wide association studies (GWAS) have been used to identify genomic regions associated with important economic traits in Nellore cattle. The GWAS uses the population linkage disequilibrium levels between markers and quantitative trait loci to identify genomic regions associated with a trait, thereby providing insights into their genetic architecture [8,9,10]. Examples of GWAS performed in Nellore cattle are: Irano [11], who uncovered chromosomal regions associated with indicator traits of sexual precocity [11]; Silva [12], who identified regions that have a significant impact on stayability [12]; Carvalho [13], who identified genomic regions influencing growth traits [13]; and Reis [14], who identified the possible genomic regions and candidate genes associated with carcass traits [14].

Although GWAS are a powerful tool for identifying genes and genetic variants associated with complex traits [10], the large number of GWAS has resulted in an extensive list of positional candidate genes for several important traits, with only a few validated through expression studies. This has created a gap in the comprehension of the genetic mechanism regulating these traits and made it difficult to determine which candidate genes are truly associated with the traits of interest. Oliveira [15] aimed to bridge this gap by fine-mapping genomic regions to pinpoint functional candidate genes specifically associated with reproductive traits in Nellore cattle [15]. By using these fine mapping techniques, the authors narrowed down the highly causal mutations and their associations with reproduction traits such as oocyte maturation and embryo development. Gene prioritization is a process through which the large list of candidate genes can be narrowed down and ranked to a fewer list of candidate genes, which could be potential causal variants for traits of interest [16]. Prioritizing candidate genes helps to optimize the fine-mapping process, as it ranks the long list of candidate genes before using fine-mapping. This is time efficient and potentializes the validation process [17].

In the absence of appropriate gene prioritization analysis, a long list of candidate genes can hamper the identification of functional candidate genes [17]. Therefore, gene prioritization and functional analyses of GWAS results are paramount to pinpoint potential functional candidate genes associated with traits currently evaluated in Nellore cattle, which can then be further validated through gene expression studies and/or fine mapping. Recent studies have performed gene prioritization and functional analyses of GWAS results in beef and dairy cattle [18,19], but to the best of our knowledge, only two studies has specifically focused on the Nellore breed [20,21]. Nonetheless, Silva [20,21] have only investigated candidate genes associated with testicular hypoplasia and feet and leg malformations, as no other traits were included in their analysis. Moreover, no previous studies have identified the shared genetic background for a group of traits (such as growth, reproduction, and carcass traits), which may be paramount in identifying the key candidate genes impacting several traits together [20,21]. Consequently, the objectives of this study were to: 1) understand the main genetic factors controlling the main trait groups in Nellore cattle (i.e., reproduction, growth, carcass, and morphological traits); and 2) pinpoint the functional candidate genes associated with these groups of traits within the Nellore cattle breed.

2. Materials and Methods

Approval from the Welfare and Animal Use Committee was not required for this study as the data was obtained from the literature.

2.1. Data Gathering and Editing

This systematic review was guided by the standards of the Preferred Reporting Items for Systematic Review and Meta-Analysis (PRISMA) Statement [22,23]. A search for genomic regions associated with the Nellore breed was performed using a combination of keywords, as detailed in the Supplementary Table S1. The keywords were created using the following criteria: terms related to economically important traits (e.g., reproduction, morphological, carcass, growth, scrotal circumference, age at puberty); type of analysis (e.g., GWAS, candidate genes, genome-wide association, genome-wide association studies); and specific cattle breed (e.g., Nellore cattle). The Keyword combination was used to search the Animal QTL database [24], National Center for Biotechnology Information [25], Web of Science (www.webofscience.com/), and Elicit (Elicit; Elicit: The AI Research Assistant; www.elicit.com; accessed November and December 2023). Covidence (www.covidence.org/) was used to extract and screen the papers. Additionally, the GALLO package [26] available in the R software [27] was also used to search for papers by focusing on the gene annotation database recorded on the Nellore breed.

A total of 561 scientific papers were first identified, genomic regions and single nucleotide polymorphism (SNP) markers reported as significant in association with target traits were recorded, and a quality control procedure was performed. The first control step consisted of removing duplicate papers (N=368) and conference proceeding abstracts that were not full-text publications (N=3); subsequently, papers with undefined traits and/or papers that do not report the genomic positions were also removed (N=14). Finally, papers that do not use purebred Nellore and/or perform genome-wide association analysis were also removed (N=62). A scheme of the quality control steps used in this study is shown in Figure 1. In the end, 1 830 genomic regions, reported in 52 scientific papers were available for the analysis.

When studies did not provide specific genomic coordinates for SNP markers and candidate genes, we retrieved them using the QTL ID, SNP ID, Ensembl gene ID, or gene name information provided in the literature, which were converted using the QTL database (QTLdb) [24], SNPchimp [28], and/or Ensembl genome browser [29]. The physical position of all genome coordinates was updated to the latest reference genome ARS-UCD1.2 [30] using Lift Genome Annotations [31], and the center of the updated genomic coordinates was used for gene annotation. Positional candidate genes were annotated within a 1MB interval (with 500Kb upstream and 500Kb downstream) using the GALLO package [26]. The final complete dataset used in this study is available in the Supplementary Table S2.

2.2. Trait Groups

The recorded traits were identified and categorized into four groups, each comprising a subgroup or indicator trait to further understand the genetic background of the traits. The groups are defined as:

Reproduction traits: Scrotal circumference, age at puberty, early pregnancy, early puberty, precocity, age at first calving, heifer rebreeding, number of calving, testicular hypoplasia, gestation length, preweaning calf mortality rate, antral follicle count, stayability, anti mullerian hormone, gestation length, calving interval, calving ease, post-natal mortality.

∙

Subgroup: Sexual precocity

Included traits: Scrotal circumference, early pregnancy, early puberty, age at first calving, antral follicle count, calving interval, calving ease.

2.

Carcass traits: Backfat thickness, rib eye area, rump fat thickness, hot carcass weight, subcutaneous fat thickness, intramuscular fat content, longissimus muscle area, marbling, tenderness, shear force tenderness, meat color.

∙

Subgroup: Meat quality

Included traits: Marbling, tenderness, shear force tenderness, meat color.

3.

Growth-related traits: Residual body weight gain, average daily gain, accumulated productivity, birth weight, weaning weight, yearling weight, adult cow weight, weight gain from birth to weaning, weight gain from weaning to yearling, yearling height, residual feed intake, dry matter intake, feed efficiency, feed conversion ratio.

∙

Subgroup: Efficiency

Included traits: Residual feed intake, dry matter intake, feed efficiency, feed conversion ratio.

4.

Morphological traits: Body conformation, muscularity, precocity, feet and leg deformation, feet and leg conformation.

∙

Subgroup: Visual scores

∙

Included traits: Conformation, precocity, muscling.

A final summary of the number of studies, genomic regions, and positional candidate genes retrieved per trait group is showed in Table 1.

2.3. Overrepresentation and Prioritization Analyses

Positional candidate genes annotated for each trait group were used in an overrepresentation analysis (ORA). The ORA was performed using the Medical Subject Headings (MeSH) vocabulary to retrieve annotations that appear more frequently in a selected gene group relative to their occurrence by chance [32]. The statistical significance (p-value) of ORA was assessed using the hypergeometric test [33]:

$$\mathrm{p}=\sum _{\mathrm{i}=\mathrm{k}}^{\mathrm{m}\mathrm{i}\mathrm{n}(\mathrm{M},\mathrm{n})}{\displaystyle \frac{\left(\begin{array}{c}\mathrm{M}\\ \mathrm{i}\end{array}\right)\left(\begin{array}{c}\mathrm{N}-\mathrm{M}\\ \mathrm{n}-\mathrm{i}\end{array}\right)}{\left(\begin{array}{c}\mathrm{N}\\ \mathrm{n}\end{array}\right)}}$$

(1)

where N is the total number of genes that were analyzed for each trait (reference genes), M is the total number of selected genes, n is the total number of genes in the MeSH term under study, k is the number of selected genes that belong to the MeSH term under study, and $\left(\begin{array}{c}\mathrm{y}\\ \mathrm{x}\end{array}\right)={\displaystyle \frac{\mathrm{y}!}{\mathrm{x}!}}\left(\mathrm{y}-\mathrm{x}\right)!$ is the binomial coefficient where y is the total number of genes analyzed and x is the number of genes selected from the set of y genes [34]. To increase the power of test, all traits within the same trait group were analyzed together in the ORA. The ORA was performed using the MESHR package [35,36]. Significant genes (i.e., adjusted p-value < 0.05) from the ORA were then considered “test genes” for the prioritization analysis.

The “guilt by association” prioritization analysis was performed using the software GUILDify and ToppGene [18,37,38,39]. The “guilt by association” principle suggests that genes with similar biological functions tend to be associated, thereby allowing the statistical inference of a gene’s function based on the association and prior knowledge of other genes [18,37,38]. First, a “trained list” of the top 200 ranked genes with a guild score > 0.6, which is a likelihood score to assess the relevance of associated preselected relevant trait keywords, was obtained from GUILDify v2.0 Web Server [19,40]. The selected keywords are as follows:

Reproduction and sexual precocity traits: “Fertility”, “implantation”, “preimplantation”, “endometrium”, “embryonic development”, “primordial follicles”, “uterus”, “luteal”, “gestation”, “embryo”, “ovulation”, “estrogen”, “estradiol”, “endocrine hormone”, “gamete, mammary glands”, “lactation”, “pregnancy”, “oocyte”, “fetus”, “zygote”, “ovary”, “amniotic”, “cervical mucosa follicle”, “Gonad”, “antral follicles”, “pre-eclampsia”, “placenta”, “testes”, “sperm”, “scrotal circumference”, “scrotal, testicular”, “testis”, “semen”, “spermatozoa”, “spermatogenesis”, “testicular hypoplasia”, “hypogonadism”.

Carcass and meat quality traits: “Longissimus muscle,” “Backfat thickness,” “Myosin,” “Actin,” “Fibroblast,” “Tissue deposition,” “Fat deposition,” “Muscle, Subcutaneous fat,” “Biceps femoris,” “Carcass weight,” “Skeletal muscle”.

Morphological and visual score traits: “body development,” “muscle mass,” “skeletal muscle,” “subcutaneous fat,” “conformational structure,” “adipose tissue,” “body depth,” “feet and leg,” “foot angle,” “tendon,” “joints,” “mobility aplomb,” “osteogenic differentiation,” “articular cartilage,” “fat deposition”.

Growth-related and Efficiency traits: “body weight”, “birth weight”, “adult weight”, “weight gain”, “body size”, “average daily gain”, “body weight gain”, “daily gain”, “stature”, “growth”, “adipose tissue”, “dry matter intake”, “metabolism”, “maternal behavior”, “maternal ability”, “feed conversion”, “feed efficiency”, “dry matter”, “milk composition”, “muscle development”.

The GUILDify software uses a Biologic Interaction and Network Analysis (BIANA) knowledge database to query gene products associated with keywords and rank the genes in this analysis [40]. This software uses a selected list of genes and species-specific (Homo sapiens) protein interaction network and applies graph theory algorithms to prioritize genes [19,40]. A network-based prioritization algorithm “Netscore” with three repetitions and two iterations was used to prioritize gene associations and phenotype protein interaction [38,39].

The ToppGene software was used to perform an annotation-based prioritization analysis through a fuzzy-based multivariate approach [37]. The functional information shared between the “trained” gene list and the “test” gene list was used to perform the multivariate analysis [37]. This functional information was retrieved from the following sources: Gene Ontology (GO) terms for molecular function (MF), biological process (BP), and cellular component (CC); as well as human and mouse phenotypes; metabolic pathways; PubMed publications; and diseases. Using a statistical meta-analysis, the p-values obtained in a random sampling of 5,000 genes from the whole genome for each annotation information were combined in an overall p-value. Gene ranking was performed by applying the PPIN-based candidate gene prioritization and the K-Step Markov method [37,40]. A total of 200 training genes were used for all trait groups and subgroups, and the number of test genes was 379 for reproduction, 277 for sexual precocity, 2 016 for carcass, 835 for meat quality, 187 for growth traits, 224 for efficiency, 91 for morphological, and 30 for visual score traits. A false discovery rate (FDR) of 5% multiple testing correction (p-value ≤ 10^-3) was used, and the significant prioritized genes were used in further analysis [18,40].

2.4. Functional Analysis

Enrichment analysis of the significant prioritized genes was performed using the Clusterprofiler to investigate their roles in the BP, MF, and CC [41]. The Bioconductor annotation package “org.Hs.eg.db” was imported for the genome-wide annotation using mapping Entrez gene identifiers for humans, as the prioritization analysis also used a Homo sapiens interaction network [41]. Significant prioritized genes have the same functional profile as the genes on the “trained” list [37].

2.5. Venn Diagrams and Gene Network Integration

The official gene names of the prioritized genes were used to integrate the available information and create the Venn diagrams. The Venn diagram tool available at www.bioinformatics.psb.ugent.be/webtools/Venn/ was used to calculate intersections between the gene lists thereby depicting what genes overlap or are unique to each trait group. To investigate possible interactions among the prioritized genes in each categorized trait, the GeneMANIA software (implemented as a plug-in on the cystoscape platform) was used [42]. This network analysis searches for related genes on publicly available biological data sets and classifies the link in the network based on their relationships, such as co-expression, physical interaction, genetic interaction, shared protein domains, co-localization, and pathway. These analyses give a better understanding of the genetic architecture of the complex polygenic traits, in order to show functions associated with genes in the network and their false discovery rate and coverage [43].

3. Results and Discussion

After quality control, a total of 52 articles were used in the analyses, as summarized in Figure 1. Of these, 23 articles focused on reproduction traits, 14 articles on carcass traits, 12 on growth traits, and five articles on morphological traits.

3.1. Gene annotation and Prioritization Analysis

The 52 articles selected contained a total of 1 761 unique genomic regions, being 605 genomic regions associated with reproduction, 387 with sexual precocity, 607 with carcass traits, 294 with meat quality, 360 with growth, 216 with efficiency, 189 with morphological traits, and 132 with visual scores. These genomic regions were annotated to identify candidate genes resulting in 8 569 positional candidate genes for reproduction traits, 5 412 for sexual precocity, 11 195 for carcass traits, 6 646 for meat quality, 5 239 for growth traits, 3 319 for efficiency traits, 3 483 for morphological traits, and 2 114 for visual score traits.

The candidate genes annotated for each trait group and subgroup were used in the overrepresentation analysis, where the numbers of significant genes identified in each category were: 379 for reproduction traits, 277 for sexual precocity, 2 016 for carcass, 835 for meat quality, 187 for growth, 224 for efficiency, 91 for morphological traits, and 30 for visual scores. These overrepresented significant genes were considered “test genes” for the prioritization analysis. For the prioritization analysis, a list of 200 ranked genes (Supplementary Table S3), derived from keywords related to the traits of interest, was used. This annotation-based prioritization analysis resulted in six prioritized genes for reproduction traits, two for sexual precocity, 10 for carcass, two for meat quality, two for growth, two for efficiency, and no prioritized genes for morphological and visual score traits. The list of prioritized genes for each trait group is included in the Supplementary Table S4.

The variation in the number of candidate genes prioritized across trait groups may be attributed to the differential focus of research efforts on specific traits in Nellore cattle, in addition to the different genetic background of the traits. For instance, out of 52 total articles used in this study, 23 focused on reproduction traits, 14 on carcass traits, 12 on growth traits, and only five on morphological traits. In this context, the relatively few studies on growth and morphological traits likely contributed to the identification of only two prioritized genes for growth and none for morphological traits. Consequently, there is a pressing need for more GWAS that focus on understanding the genomic regions associated with growth and morphological traits in Nellore cattle.

3.2. Enrichment and Gene Integration Analysis of the Functional Candidate Genes Identified within Groups

3.2.1. Reproduction and Sexual Precocity Traits

The six candidate genes prioritized for reproductive traits were GFRA4, RFWD3, SERTAD2, KIZ, REM2, and ANKRD34B. Among these, the ANKRD34B gene is common to both reproductive traits and sexual precocity, while the ACSM4 gene is uniquely prioritized for the sexual precocity trait. Gene functions of the six prioritized genes for reproductive traits and the two for sexual precocity traits were identified using GeneMANIA. The results indicated significant interactions among the genes associated with both traits: 80.30% and 77.64% physical interactions, suggesting that the genes are linked through protein-protein interactions; 9.48% and 5.37% predicted interactions, indicating that the genes likely participate in the same pathway reactions; 4.28% and 2.87% genetic interactions, showing that the genes are functionally associated, with the effect of one gene being modified by another; 3.72% and 8.01% co-expression, demonstrating similar patterns across gene expression studies; 0.74% and 3.63% co-localization, indicating that the genes are expressed in the same tissues; and 0.37% and 0.60% shared protein domains, showing common protein domain data among the genes.

Related genes, such as CEP72, TUBGCP2, and ESPNL (Figure 2a), and ACSM2A and SLC27A1 (Figure 2b), are associated with the prioritized genes through various biological functions. These genes are involved in the organization of the mitotic cell cycle, binding of actin filaments, and expression of reproductive traits such as the bipolar spindle, endometrium, and ovary [44,45,46]. Additionally, they encode mitochondrial acyl-coenzyme A and other fatty acids needed in the pathway [47,48].

The functional enrichment analysis of the reproductive traits highlighted several key genes. The KIZ gene plays a crucial role in spindle organization in Nellore cattle, essential for forming a robust mitotic centrosome architecture capable of withstanding the forces exerted on the centrosome during spindle formation [49]. The SERTAD2 gene is involved in the regulation of cell growth, which is important in the reproduction process. For instance, Zhang et al. [50] observed that patients with lower levels of lnc-SERTAD2-3, a member of the SERTAD family, had an unfavorable prognosis for osteosarcoma [50]. In contrast, patients with higher levels of lnc-SERTAD2-3 experienced inhibited cell growth and spreading, confirming the role of the SERTAD2 gene in cell growth regulation. Additionally, Darwish [51] found that inhibiting SERTAD3, another member of the SERTAD gene family, resulted in a two to four-fold reduction in cell growth rate [51].

The RFWD3 gene is associated with DNA recombination, DNA replication, DNA damage checkpoint signaling, and mitotic DNA damage regulation of the cell cycle phase. Fu et al. [52] investigated the role of the RFWD3–Mdm2 ubiquitin ligase complex in regulating p53 protein stability following DNA damage [52]. They identified that RFWD3 forms a complex with p53 to synergistically facilitate its ubiquitination, thereby playing a crucial role in cellular protection by regulating the response to DNA damage. The REM2 gene is associated with guanyl nucleotide (GTP) binding, a protein crucial for the development and maturation of reproductive cells. Casalotti et al. [53] investigated the expression of GTP binding proteins during the development of rat testes and discovered that the cellular expression of G proteins is temporally linked to testicular development in rats [53].

The ACSM4 gene plays a significant role in the mitochondrial matrix, which is crucial for cellular energy metabolism [54]. The mitochondrial matrix houses the citric acid cycle and enzymes such as C-acyltransferase and CoA-ligase, which are important for improving steroidogenesis, a process essential for the secretion of steroid hormones such as testosterone [55,56].

3.2.2. Carcass and Meat Quality Traits

The ten candidate genes prioritized for carcass traits were RFWD3, TMEM120A, MIEF2, FOXRED2, DUSP29, CARHSP1, OBI1, JOSD1, NOP58, and LOXL1-AS1. Among these, CARHSP1 and LOXL1-AS1 were also prioritized for meat quality traits. However, the genes DUSP29 and LOXL1-AS1 were not recognized in the GeneMANIA network. The LOXL1-AS1 gene is a non-protein coding antisense RNA gene, and DUSP29 was likely not found in GeneMANIA due to the lack of existing data sources in the interaction database.

The eight recognized prioritized genes for carcass traits (i.e., RFWD3, TMEM120A, MIEF2, FOXRED2, OBI1, JOSD1, NOP58) showed 71.22% co-expression, indicating similar expression levels across conditions, and 28.78% co-localization, indicating that these genes are expressed in the same tissues associated with carcass traits. The one prioritized gene for meat quality traits, CARHSP1, showed 77.64% physical interaction, suggesting that the genes are linked through protein-protein interactions; 8.01% co-expression, suggesting similar expression levels across conditions; 5.37% predicted interaction and 1.88% pathway participation, suggesting that the genes likely participate in the same pathway reactions; 3.63% co-localization, indicating expression in the same tissues; 2.87% genetic interaction, showing functional association with the effect of one gene being modified by another; and 0.60% shared protein domain, indicating common protein domain data among the genes.

Related genes such as ALKBH7, CDAN1, DYRK2, HSPA1A, and SNAP47 are involved in stabilizing existing proteins, cellular response to DNA damage, and biological functions related to carcass or meat quality, such as the nuclear envelope, which surrounds the nucleus, and intracellular membrane fusion [57,58,59,60,61].

The OBI1 gene is involved in the regulation of DNA replication and protein autoubiquitination and/or deubiquitination [62]. Nassar et al. [63] studied the mechanism behind OBI1 function using Xenopus laevis eggs and confirmed that the OBI1 gene encodes an E3 ubiquitin ligase, a protein essential for DNA replication in the genome. This gene plays a critical role in regulating DNA replication and maintaining cellular homeostasis through protein autoubiquitination or deubiquitination processes, which are important for ensuring good meat quality [63].

The TMEM120A gene is associated with sensory perception, which is important in the detection of mechanical external stimuli and fat cell differentiation [64]. Batrakou et al. [65] investigated the roles of two nuclear envelope transmembrane proteins, TMEM120A and TMEM120B, in adipocyte differentiation and metabolism. They confirmed that the nuclear envelope transmembrane protein TMEM120A is preferentially expressed in fat and is required for adipocyte differentiation and metabolism [65].

The DUSP29 gene is involved in protein dephosphorylation, which activates or inactivates enzymes and is also involved in muscle cell differentiation, an important characteristic of meat quality [66]. Cooper et al. [67] determined how DUSP29 is transcriptionally regulated in skeletal muscle and confirmed the importance of DUSP29 in muscle cell differentiation, providing insight into molecular and cellular mechanisms for skeletal muscle atrophy [67].

The JOSD1 gene is involved in protein modification and deubiquitination, reversing the ubiquitination process and preventing protein degradation, thereby improving meat quality [68]. Seki et al. [69] provided additional insights into the properties of the Josephin domain (JOSD), and confirmed that JOSD1 is activated by ubiquitination and helps regulate membrane dynamics [69].

The two candidate genes prioritized for meat quality traits, which were also prioritized for carcass traits, are LOXL1-AS1 and CARHSP1. LOXL1-AS1 is a non-protein coding antisense RNA [70], while CARHSP1 is involved in RNA binding to the 3’-UTR ends of histones [71], a process recently linked to adipose tissue accumulation in pigs [72]. Kociucka et al. [72] investigated the association between histone modifications and metabolism in pig adipose tissues, helping to understand adipose tissue development, which is essential for improving meat quality in pigs [72].

3.2.3. Growth-Related and Efficiency Traits

The candidate genes prioritized for growth-related and efficiency traits were the same, indicating that the genes identified are key genes regulating the expression of both growth-related and efficiency traits. The two prioritized genes for these trait groups are ANKRD34B and JOSD1, and their interaction with other genes is showed in Figure 4. A total of 77.64% physical interaction, 8.01% co-expression, 5.37% predicted interaction, 2.87% genetic interactions, 1.88% pathway, and 0.60% shared protein domains were observed. Related genes, such as TIMM8A, EIF3F, PJA1, and TENT4A are involved in insertion of hydrophobic membrane in the mitochondrial membrane and associated to enabling protein deubiquitination and DNA repair [73,74].

In addition to its role in reproduction processes [75], such as germ and pole plasm, ANKRD34B is also associated with ribonucleoprotein granules—a complex of RNA and RNA-binding proteins that play essential roles in regulating cell growth and development [76]. To address limitations related to ribonucleoproteins, Tartaglia [76] explored the crucial role of protein-RNA interactions in cellular regulation and the gaps in understanding how ribonucleoprotein complexes assemble and function. Their study found that these ribonucleoproteins are essential for various cellular processes, including growth [76]. Regarding the JOSD1 gene, this gene is involved in protein deubiquitination, which reverses the ubiquitination process and prevents protein degradation thereby improving growth [68]. This gene also plays a significant role in cellular processes, as it participates in the localization of plasma membranes, affecting endocytosis [69].

3.2.4. Morphological and Visual Score Traits

No candidate genes were prioritized for morphological and visual score traits; therefore, no gene was recognized by the network integration algorithm. The limited number of studies on morphological and visual score traits likely contributed to identifying no prioritized genes for both trait groups, as only five and three scientific publications were found for morphological and visual score traits, respectively.

3.3. Genes Shared among Groups

3.3.1. Overrepresented Genes

After assessing the overrepresentation of positional candidate genes, the intersection between gene lists was evaluated (Figure 5). Of the genes analyzed, 137 out of 379 were unique to reproduction traits, 348 out of 2 016 were unique to carcass traits, 79 out of 187 were unique to growth traits, and 34 out of 91 were unique to morphological traits.

In terms of shared genes, one gene was common among all four trait groups. Eight genes were common to reproduction, carcass, and growth traits, while two genes were common to reproduction, growth, and morphological traits. Only one gene was shared among carcass, growth, and morphological traits. Additionally, 30 genes were shared between reproduction and carcass traits, 14 genes between reproduction and growth traits, and seven genes between reproduction and morphological traits. The carcass and growth traits shared 15 genes, the carcass and morphological traits shared four genes, and the growth and morphological traits shared two genes. The complete list of overrepresented genes identified per trait group and subgroup is included in the Supplementary Table S5.

3.3.2. Prioritized Genes

Four genes (i.e., GFRA4, KIZ, REM2, and SERTAD2) were found to be unique to the reproduction trait, while eight genes (i.e., OBI1, NOP58, DUSP29, TMEM120A, LOXL1-AS1, MIEF2, CARHSP1, and FOXRED2) were unique to the carcass trait. Notably, no gene was found to be unique to the growth traits (Figure 6). The reproduction and carcass traits shared the gene RFWD3, the reproduction and growth traits shared the gene ANKRD34B, and the carcass and growth traits shared the gene JOSD1.

3.4. Genes Shared between Groups and Subgroups

The intersection of genes between trait groups and their respective subgroup are showed in Figure 7. Six prioritized genes—GFRA4, RFWD3, SERTAD2, KIZ, REM2, and ANKRD34B—are unique to reproduction traits. Among these, ANKRD34B is also common to both reproduction and sexual precocity traits, with ACSM4 being uniquely prioritized for sexual precocity.

The intersection between carcass and meat quality traits revealed two commonly prioritized genes: LOXL1-AS1 and CARHSP1. Eight prioritized genes—RFWD3, DUSP29, TMEM120A, CARHSP1, OBI1, LOXL1-AS1, FOXRED2, and NOP58—are unique to carcass traits, with no genes unique to meat quality traits. In the intersection between growth and efficiency traits, two prioritized genes—ANKRD34B and JOSD1—are shared, with no genes unique to either trait. No prioritized genes were identified for morphological and visual score traits.

The Venn diagram was used to calculate intersections between gene lists and illustrate the overlap of candidate genes among reproduction, carcass, and growth traits in Nellore cattle. By examining these intersections, we identify shared genes that may be pleiotropic, meaning that one gene affects the expression of multiple phenotypes [77].

The RFWD3 gene, common to both reproduction and carcass traits, plays a role in DNA repair, DNA replication, DNA recombination, and ubiquitination [78]. Efficient DNA repair mechanisms are crucial for maintaining genetic integrity during gametogenesis and embryonic development—a key aspect of reproduction [79]. Inano et al. [80] studied RFWD3's role as an E3 ubiquitin ligase in the DNA damage response. Their research revealed that RFWD3 facilitates homologous recombination and proper DNA repair by ubiquitinating and removing RPA and RAD51 proteins, preventing interference with the repair process [80]. The gene's involvement in DNA repair, replication, and recombination contributes to meat quality and yield, highlighting its dual role in both reproduction and carcass traits [78,79,80].

The ANKRD34B gene is shared between reproduction and growth traits. Known for its association with ribonucleoprotein granules—a complex of RNA and RNA-binding proteins—ANKRD34B plays essential roles in regulating cellular processes such as cell regulation, growth, and development [75]. This may influence gametogenesis or early embryo development [76]. Further research is needed to understand how ANKRD34B influences the reproduction process. The JOSD1 gene is common to both carcass and growth traits. This gene, involved in protein modification and deubiquitination, helps prevent protein degradation, thereby improving growth-related and meat quality traits. JOSD1's role in maintaining cellular protein balance through deubiquitination is crucial for muscle quality and growth regulation, impacting muscle composition and carcass quality [69].

4. Limitations of this Study and Future Directions

A primary limitation of this study was the absence of a cattle-specific database in the GUILDify software, which required the use of a human (Homo sapiens) protein interaction network for gene prioritization. Although cattle share approximately 80% of human orthologs [59,60], using a human database instead of a cattle-specific one may have affected the results. To address this limitation, future efforts should focus on developing and integrating a comprehensive cattle-specific protein interaction database into GUILDify. Additionally, enhancing cross-species mapping algorithms and improving cattle genome annotation could ensure that the genetic and protein interaction data used in studies are more accurate and relevant to cattle, thereby mitigating the impact of using human databases. Furthermore, GUILDify employs multiple algorithms (e.g., NetScore, NetZcore, NetShort) to prioritize genes based on protein-protein interactions [38]. However, this study utilized the default NetScore algorithm with a repetition of three and an iteration of two, as recommended by the software developer. Future studies could compare the performance of multiple algorithms to determine their impact on the results. Regardless, it is important to note that the GUILDify software ranks these genes based on network topology, while the ToppGene software ranks the candidate genes based on functional similarity. These differences may help reduce the potential bias from the algorithm used for GUILDify in this study.

Moreover, this study was limited to GWAS performed using SNP panels, which were more widely available for the Nellore breed. As a result, some significant genomic markers associated with traits of interest may have been overlooked. To address this limitation, future studies should incorporate a broader range of molecular technologies, including whole-genome sequencing, copy number variation analysis, transcriptomics, and epigenomics. Integrating data from these sources could provide a more comprehensive understanding of the genetic architecture of these traits and ensure that all relevant genetic markers are identified.

5. Conclusions

The candidate genes identified in this study enhance our understanding of the genetic architecture and key genetic mechanisms controlling reproduction, carcass, growth, and morphological traits in the Nellore breed. The observed overlap of genes among traits underscores their complexity and indicates that many of these genes have cross-phenotype effects, where a single candidate gene is associated with multiple traits and shares common genetic pathways. The functional candidate genes identified in this study can serve as targets for gene expression and fine mapping studies, which are crucial for validating these findings and further elucidating their roles in trait regulation.