QTLs and Candidate Genes Associated with Semen Traits in Merino Sheep

Simple Summary Ram semen traits including volume, gross motility, concentration, and percent post-thaw motility are routinely assessed prior to use in artificial breeding programs, as they have been shown to influence conception outcomes. Semen quality and associated traits are complex but heritable, and as such, identifying genes that underlie variability in these traits may help develop alternative means to improve conception outcomes and therefore reproductive efficiency in sheep. Therefore, the aim of this study was to identify genomic regions and associated genes that may significantly influence semen traits like volume, gross motility, concentration, and percent post-thaw motility in Merino sheep. Assessment of over 20 years’ worth of semen collection data identified 35 genomic regions to be significantly associated with Merino ram semen volume, gross motility, concentration, and percent post-thaw motility. A total of 290 candidate genes were identified within genomic regions found to be significantly associated with Merino ram semen traits. All Merino ram semen traits were also found to be lowly heritable, affirming results from previous studies. Validation of candidate genes identified in the current study could provide novel insights into molecular mechanisms contributing to variability in semen associated traits. Abstract Ram semen traits play a significant role in conception outcomes, which in turn may influence reproductive efficiency and the overall productivity and profitability of sheep enterprises. Since hundreds of ewes may be inseminated from a single ejaculate, it is important to evaluate semen quality prior to use in sheep breeding programs. Given that semen traits have been found to be heritable, genetic variation likely contributes to the variability observed in these traits. Identifying such genetic variants could provide novel insights into the molecular mechanisms underlying variability in semen traits. Therefore, this study aimed to identify quantitative trait loci (QTLs) associated with semen traits in Merino sheep. A genome-wide association study (GWAS) was undertaken using 4506 semen collection records from 246 Merino rams collected between January 2002 and May 2021. The R package RepeatABEL was used to perform a GWAS for semen volume, gross motility, concentration, and percent post-thaw motility. A total of 35 QTLs, located on 16 Ovis aries autosomes (OARs), were significantly associated with either of the four semen traits in this study. A total of 89, 95, 33, and 73 candidate genes were identified, via modified Bonferroni, within the QTLs significantly associated with volume, gross motility, concentration, and percent post-thaw motility, respectively. Among the candidate genes identified, SORD, SH2B1, and NT5E have been previously described to significantly influence spermatogenesis, spermatozoal motility, and high percent post-thaw motility, respectively. Several candidate genes identified could potentially influence ram semen traits based on existing evidence in the literature. As such, validation of these putative candidates may offer the potential to develop future strategies to improve sheep reproductive efficiency. Furthermore, Merino ram semen traits are lowly heritable (0.071–0.139), and thus may be improved by selective breeding.


Introduction
Reproductive efficiency in livestock depends on multiple factors, including insemination, conception, embryonic development, and parturition, and is a key determinant of farm productivity and profitability [1,2]. Reproductive efficiency is influenced by an interplay between the environment and host genetics, which in turn includes both maleand female-side genetic determinants. In this context, significant emphasis has been placed on identifying genetic variants associated with female reproductive efficiency in livestock species in the past [3][4][5], while relatively little emphasis has been placed on identifying male-side genetic variants influencing livestock reproductive efficiency. Furthermore, there is only one past sheep study that has investigated the influence of genetic variants associated with semen traits, which was performed using a sheep breed that is not common across Australia [6]. In livestock species such as sheep, a single ram has the potential to produce hundreds of progeny annually via natural breeding or artificial insemination (AI) [7]. For example, top-ranking Merino sires on different selection indices in MERINOS-ELECT frequently have hundreds of progeny recorded across multiple flocks [8]. Therefore, it is important to assess the breeding soundness of rams prior to natural breeding, or alternatively, assess semen quality prior to artificial breeding programs. Likewise, the sire has a significant influence on the number of lambs born, with the breed of ram also influencing prenatal mortality [9]. Poor conception outcomes, in turn, may result in fewer lambs being born and may also potentially extend the joining and lambing periods, which could have subsequent negative impacts on the productivity and profitability of sheep farming enterprises.
In the context of semen traits, previous genome-wide association studies (GWAS) have identified several QTLs and candidate genes in a range of livestock species, including dairy sheep [6], goats [25], cattle [26][27][28][29][30][31], and pigs [32][33][34][35][36][37]. For example, a GWAS in Assaf sheep reported 20, 23, 76, and 32 single nucleotide polymorphisms (SNPs) on chromosomes 3, 17, 4, and 1 to be associated with semen volume, gross motility, concentration, and number of spermatozoa per ejaculate, respectively [6]. In the same study, SLC9C1 and FUT10 were identified as putative candidate genes for gross motility due to their location within identified QTL regions, and because of previous reports that these genes influenced spermatozoal motility, hyperactivation, and fertilisation, respectively. Similarly, a GWAS involving semen volume and number of sperm collected from Holstein-Friesian bulls identified several candidate genes such as GALC, PHF7, and PRKCD, all of which have been previously identified to play a role in spermatogenesis [38]. Similar investigations in pigs have identified several putative candidates associated with post-thaw motility and membrane integrity, including PLBD1, OXSR1, and EML5 [39], and polymorphisms in such genes have been previously reported to be associated with fertilisation in humans [40], spermatozoal motility in pigs [41], and spermatogenesis in bulls [42], respectively. Despite several past GWAS focused on identifying genomic regions and candidate genes associated with semen traits across a wide range of livestock species, the only previous study in sheep was conducted using Assaf rams [6], and there are no studies involving other common sheep breeds like the Merino. While it is likely that the physiology associated with semen traits is mechanistically conserved across species, the genes that contribute to variability in these traits may vary depending on the genetic background of any given population of animals. In Australia, Merino sheep account for 74% of Australia's ewe breeding flock and 56% of the lamb flock, with 42.4 million Merino ewes bred annually [43]. Therefore, the objective of this study was to perform a genome-wide association study to identify QTL regions and associated candidate genes that significantly influence semen traits in Merino rams.

Phenotype Data
Phenotype data, including semen volume, gross motility, concentration, and percent post-thaw motility, was obtained from a single commercial artificial breeding company (which had two centres over time). The collection and assessment of semen quality, descriptive statistics, and data quality control have been previously described [22]. The dataset for Genome-Wide Association (GWA) analysis consisted of a total of 4506 semen collection records from 246 Merino rams that had been previously genotyped by their flock owners (hence, 246 rams had both semen phenotype data and had been genotyped). Descriptive statistics for semen traits (Table 1) also fall within the normal ranges of semen traits for ram [7]. Note: Volume of ejaculate in millilitres (mL); gross motility (score from 0 to 5); concentration (spermatozoa × 10 6 ); and percent post-thaw motility (%).

Genotype Data and Quality Control
Genotype data for 246 rams included in this study was provided by Sheep Genetics (Meat & Livestock Australia). Animals were genotyped using different Ovine SNP chips ranging from 12k (n = 31), 15k (n = 81), 50k (n = 122), and 700k (n = 12). Low-density chips were subsequently imputed to 50k using a reference population as part of the routine pipeline followed by the Animal Genetics and Breeding Unit (AGBU), which was previously described [44]. Briefly, Beagle software v3.2 [45] was used to impute 12k and 15k to 50k SNPs. PLINK [46] was then used to perform quality control, excluding SNPs based on the following criteria: < 0.01, Hardy-Weinberg equilibrium p-value < 1 × 10 −6 , and if they were located on sex chromosomes. In total, 58,946 SNPs were used for subsequent GWA analysis. The genomic positions of SNPs were taken from the ovine (Ovis aries) reference genome build (Oar_v3.1).

Statistical Model
Preliminary analysis was performed to assess the significance of each fixed effect for individual semen traits via linear regression using the R software package, version 4.1.2 [47]. The fixed effects that were assessed included age of collection (n = 8), collection centre (n = 2), season (n = 4), method of collection (n = 2), and collection number. Once the fixed effects that significantly influenced each semen trait were identified (Supplementary Table S1), a linear mixed model was used to perform genome-wide association analysis using a model that included several fixed effects, a random polygenic effect, and the random effect associated with the stud of breeding and year of collection: where y represents the vector of transformed phenotypes for individual semen traits; u is an intercept term; X SNP and β SNP represent SNP dosage and effects, respectively; Z g and Z p are incidence matrices relating to the individuals and their observed values for random polygenic and random stud-year effects; and e represents the residuals [48].

Identification of SNPs and Associated Candidate Genes
The R package RepeatABEL v1.1 [48] was used for GWA analysis, as it allows for the analysis of datasets with repeated measures. A 'prefitModel' function was initially used in RepeatABEL by fitting a linear mixed effect model that included the genomic relationship matrix (GRM). This accounted for random polygenic effects without including SNP effects. Secondly, the 'rGLS' function was used in RepeatABEL to test for SNP association by accounting for the fixed (Supplementary Table S1) and random effects of fitted models for each semen trait. The p-values for each SNP, as well as the SNP effect of each marker, were determined via Wald statistics within the RepeatABEL software package. Finally, to accurately account for multiple testing levels, the significance threshold was calculated via modified Bonferroni [49] as follows: α √ nSNP where α = 0.05 and nSNP is the total number of SNPs identified across each of the four semen traits. Candidate genes were identified in QTL regions that spanned 0.5 Megabase (Mb) around each significant SNP (p-value < 0.0002). These QTL regions were used to search for candidate genes through the National Centre of Biotechnology Information (NCBI) Genome Data Viewer [50] using the ovine (Ovis aries) genome assembly (Oar_v3.1).
The SNP variance was calculated as the total variance (as previously described [51]) explained by each SNP as a proportion of the additive genetic variance using the following formula: where p and q are the major and minor allele frequencies, respectively; a is the additive SNP effect; and σ 2 a is the additive genetic variance.

Heritability and Variance Component Estimation of Semen Phenotypes
Heritability and variance components were estimated for semen volume, gross motility, concentration, and percent post-thaw motility using the 'prefitModel' in the 'rGLS' function, including the GRM, via the R package RepeatABEL [48]. The same significant fixed effects used for the GWAS for each semen phenotype were fitted in the model, as defined (Supplementary Table S1). Heritability (h 2 ) of semen phenotypes was calculated as follows: where σ a 2 is the estimate for additive genetic variance; σ c 2 is the variance estimate for random permanent environmental genetic effects; σ d 2 is the estimate of random stud-year effect; σ e 2 is the estimate of the residual error variance; and σ pe 2 is the phenotypic variance (being the sum of additive genetic, random permanent environmental, stud-year, and residual variance effects).
Estimates of repeatability (r) were calculated as follows: where σ a 2 , σ c 2 , σ d 2 , and σ pe 2 have been previously defined.

Results
Several SNPs were found to be significantly associated with each semen trait included in this study, i.e., volume, gross motility, concentration, and percent post-thaw motility. Manhattan plots associated with GWA analysis of each of the semen traits are presented in Figure 1.
A total of 35 significant SNPs were found across a total of 16 OARs (Table 2). Of those, nine, seven, eight, and eleven SNPs were associated with volume, gross motility, concentration, and percent post-thaw motility, respectively. Furthermore, SNPs for semen volume, gross motility, concentration, and percent post-thaw motility were found across seven, six, six, and eight unique OARs, respectively. The SNPs with the highest significance identified in each semen trait were 12:48,412,883, 10:16,620,353, 20:42,505,407, and 17:22,669,197 for semen volume, gross motility, concentration, and percent post-thaw motility, respectively.
Overall, 143, 141, 97, and 138 candidate genes were identified within significant genomic regions associated with semen volume, gross motility, concentration, and percent post-thaw motility after searching the QTL regions, which spanned 0.5 Megabase (Mb) around each significant SNP associated with semen traits via the NCBI Genome Data Viewer using the Oar_v3.1 genome assembly. Several genes located within QTLs were not annotated or did not have known orthologs, and as such, these genes had a 'LOC' prefix. A summary of all these genes is presented in Supplementary Table S2. Therefore, a total of 35 QTLs, specifically nine, seven, eight, and eleven regions, were associated with volume, gross motility, concentration, and percent post-thaw motility, respectively. Table 3 presents candidate genes identified within each QTL region associated with semen volume, gross motility, concentration, and percent post-thaw motility. Unique candidate genes found within QTL regions for each of the semen traits include 89 for volume, 95 for gross motility, 33 for concentration, and 73 for percent post-thaw motility. Several candidate genes were identified within the most significant QTLs for each semen trait, including volume PRDM16, ACTRT2, TRNAC-GCA, TTC34, MMEL1, FAM213B, HES5, PANK4, PLCH2, PEX10, RER1, MORN1, SKI, C12H1orf86, PRKCZ, GABRD, CFAP74, TMEM52, CALML6, gross motility; SPERT, SIAH3, ZC3H13, CPB2, LCP1, LRRC63, KIAA0226L, LRCH1, ESD, TRNAG-CCC, and concentration; SIRT5, GFOD1, TBC1D7, PHACTR1, and percent post-thaw motility; no candidate genes identified within QTL.
The heritability estimates for semen volume, gross motility, concentration, and percent post-thaw motility are presented in Table 4. Heritability estimates for volume and gross motility are lowly moderate (0.104 and 0.139, respectively), while concentration and percent post-thaw motility have low heritability estimates (0.071 and 0.092, respectively). The additive genetic, environmental, residual, and random effects variances are higher for concentration in comparison to the other semen traits.     The genes with a superscript A were previously found to be associated with spermatogenesis and male fertility; OAR: Ovis aries autosome; Mb: Megabase.

Discussion
In sheep, semen traits including concentration [19], gross motility [18], and percent post-thaw motility [20] are reported to influence conception outcomes, but the underlying mechanisms are not well understood. Identifying genes and genomic regions associated with such traits can help provide novel insights into these mechanisms. Therefore, the aim of this study was to investigate regions on the ovine genome and candidate genes associated with Merino ram semen traits including volume, gross motility, concentration, and percent post-thaw motility.
This study identified a total of 35 QTLs located on 16 OARs, specifically nine, seven, eight, and eleven QTLs for semen volume, gross motility, concentration, and percent post-thaw motility, respectively ( Table 3). The nine QTLs identified to be associated with semen volume were located on OARs three, five, nine, ten, eleven, twelve, and twenty-five. Eighty-nine unique candidate genes were identified within these QTLs for semen volume. Of these, several candidates were noteworthy because they have previously been reported to be expressed within testes and sperm cells as well as associated with spermatozoal maturation. Examples of such genes include solute carrier family 2 member 8 (SLC2A8), actin-related protein T2 (ACTRT2), ribosomal protein L12 (RPL12), and MIS12 kinetochore complex component (MIS12). SLC2A8, also known as Glucose transporter 8 (GLUT8), is highly expressed in the testes and is located within the acrosome of spermatozoal cells as well as in the plasma membrane [52]. In the same study, SLC2A8 mRNA synthesis was found to coincide with spermatozoal maturation in the mouse testes. In other murine studies, spermatozoal motility was significantly impaired in SLC2A8 null mice due to reduced mitochondrial function [53], while oocytes harvested from SLC2A8 null mice had irregular changes to the endometrium during embryo implantation [54]. As such, SLC2A8 may be a putative candidate for ram spermatogenesis and motility. Another such gene, ACTRT2, is reported to be localised in the post-acrosomal region of mouse spermatozoa [55], as well as the post-acrosomal region and midpiece in males with low spermatozoal motility which resulted from obesity [56]. Therefore, ACTRT2 may play a key role in the normal morphological development of the acrosome and membrane integrity, as well as promoting spermatozoal motility. Proteomic analysis of ram semen identified ACTRT2 as being significantly expressed in fresh semen compared to cryopreserved semen [57]. In a pig study, ACTRT2 was significantly downregulated in fresh semen collected during the nonbreeding season [58], which is associated with reduced semen quality in pigs [59]. As time of year influences semen quality in both pigs [60,61] and sheep [10,11], ACTRT2 may be a key gene modulating semen quality during the breeding and non-breeding seasons. Another gene, RPL12, is known to be significantly downregulated in bulls with low fertility, as assessed via conception rate following artificial insemination [62]. RPL12 has also been found to be upregulated in ram testes [63]. Given that spermatogenesis, maturation of spermatozoa, and spermatozoal motility occur within the testes [64], RPL12 may be a positional candidate gene that plays a key role in promoting ram semen maturation and motility. MIS12 has been previously reported to be enriched within networks associated with chromatin condensation following proteomic analysis of turkey testes [65]. Given that chromatin condensation is correlated with the fertilising ability of a sperm [66,67], MIS12 may play a crucial role in sperm quality and conception success. Furthermore, proteomic analysis of seminal plasma identified MIS12 as being commonly expressed across all fertile males [68]. As semen consists of spermatozoa and seminal plasma, which significantly enhance ram sperm motility and membrane integrity [69], MIS12 may be a putative candidate gene within ram seminal plasma that may modulate sperm motility and viability. Thus, associated positional candidate genes identified in QTLs significant for semen volume may also have functional roles in overall semen quality.
Seven QLTs, located on six different OARs, were significantly associated with gross motility. A total of 95 candidate genes were identified within these seven QTL regions, of which several key genes have been previously associated with spermatogenesis and semen quality, including SH2B adaptor protein 1 (SH2B1), mitogen-activated protein kinase 3 (MAPK3), aldolase, fructose-bisphosphate A (ALDOA), 5 -nucleotidase ecto (NT5E), and spermatid associated protein (SPERT). SH2B1 may be a putative candidate influencing spermatozoal motility since SH2B1 was identified within one of the significant QTLs (18.51-18.87 Mb on Sus scrofa chromosome 3) associated with progressive motility in Duroc boars [33]. Furthermore, reduced mRNA expression of SH2B1 in obese mice resulted in reduced motility and hyperactivation compared to control mice [70]. Another gene, MAPK3, is associated with spermatozoal hyperactivation [71], which is a crucial process required prior to sperm fertilising an ova [72]. Moreover, MAPK3 was also identified as a common differentially expressed gene (DEG) in semen collected from Merino rams when contrasting breed specific differences in the spermatozoal transcriptome [17], as well as a candidate gene associated with sire fertility traits following a GWAS in Holstein cattle [73]. Hence, MAPK3 may be a key putative candidate underlying ram spermatozoal motility. Another candidate, ALDOA, was abundantly expressed in spermatozoa collected from rams with high motility compared to rams with low motility [74]. Likewise, ALDOA was abundantly expressed in rams with higher pregnancy rates following natural mating compared to rams with lower pregnancy rates across three breeding seasons [75]. This suggests that expression of ALDOA may contribute to spermatozoa motility in rams and thus may be a key putative candidate for spermatozoal motility. Proteomic analysis of ram spermatozoa identified NT5E to be abundantly expressed in ram semen with high motility following cryopreservation (and therefore, a high tolerance to cryopreservation) compared to semen with a low cryo-tolerance [76]. Therefore, NT5E may be a key functional candidate for spermatozoal motility and cryo-tolerance in sheep. Additionally, NT5E was found to be enriched in both seminal plasma and seminal vesicle fluid of bulls [77]. Finally, SPERT, also known as chibby family member 2 (CBY2), was significantly expressed in the testes collected from buffalo during puberty compared to testicular tissue from both pre-and post-pubertal buffalo [78]. In the same study, SPERT was found to be localised to the acrosome, neck, and midpiece of buffalo spermatozoa. As such, SPERT may be a key candidate for spermatogenesis and spermatozoal morphology, given that it is significantly expressed during puberty and is localised to the sperm head and midpiece.
Thirty-three candidate genes were identified within eight QTLs, located on six OARs, that were found to be associated with semen concentration. Some candidate genes identified in this study have previously been reported in studies related to spermatogenesis and testicular development, including peptidyl arginine deiminase 2 (PADI2), cytochrome P450, family 19 (CYP19), and actin like 8 (ACTL8). PADI2 may be a putative candidate gene influencing spermatozoal concentration and motility since PADI2 is abundantly expressed in the epididymis, predominantly in the caput epididymis and vas deferens following RNA sequencing of murine testes [79]. In several mammalian species, transit through the epididymis is crucial for spermatozoal maturation and attainment of motility [80][81][82][83], as well as ensuring adequate concentration of spermatozoa in semen through fluid absorption [84]. Hence, PADI2 may be a putative candidate associated with spermatozoal concentration and a positional candidate modulating spermatozoal maturation and motility. In the same context, PADI2 was found within one of ten Bos taurus autosomes (BTA), explaining the highest portion of genetic variance for both total and progressive motility (BTA 2:1349.81-1359.81) [29]. In the current study, PADI2 is located within OAR 2:247.37-248.37, and as such, PADI2 may be a significant putative candidate for spermatozoal motility in Merino sheep. Furthermore, PADI2 was identified within a QTL associated with the number of piglets born dead and alive [85]; hence, PADI2 may play a key role in litter size and post-parturition offspring survival in livestock. Another noteworthy candidate gene, CYP19, is associated with spermatozoal concentration and motility in humans, with polymorphisms of CYP19 resulting in reduced sperm concentration within normozoospermic males [86]. The authors suggested that altered aromatase concentrations within the testes resulted in reduced concentrations of spermatozoa within normozoospermic males. Aromatase promotes oestrogen production in humans and mice [87][88][89] and is found to be expressed in several cell types, e.g., Sertoli and Leydig cells, spermatids in rats [90], as well as semen in humans [91,92]. A past pig study revealed CYP19 to be abundantly enriched in spermatozoa, which had higher embryonic cleavage rates compared to spermatozoa, which had lower cleavage rates [93]. Therefore, CYP19 may be a key candidate for modulating the concentration of spermatozoa in an ejaculate. Expression of ACTL8 has been reported to be isolated to the testes following a study that undertook exome sequencing in males with non-obstructive azoospermia [94], a condition where no sperm is present in an ejaculate due to failures during spermatogenesis [95]. In the same study, mutations in ACTL8 were not present in the control patients. Hence, ACTL8 may significantly influence spermatogenesis and spermatozoal concentration within an ejaculate. Thus, noteworthy genes identified in QTLs significantly associated with concentration of Merino ram semen, such as PADI2, CYP19, and ACTL8, may be putative candidates for semen traits like concentration and motility in sheep.
A total of eleven QTLs on eight different OARs were found to be significantly associated with percent post-thaw motility. A total of 73 candidate genes were identified within these QTL regions. Of these 73 candidates, several have previously been reported to be associated with semen quality and male-side influence on conception outcomes, including MON1 homolog B, secretory trafficking associated (MON1B), RAB3B, member of the RAS oncogene family (RAB3B), ribosomal protein L30 (RPL30), and sorbitol dehydrogenase (SORD). MON1B, located on OAR 14, was identified within genomic regions associated with both cow-and daughter-conception rates following a GWAS in Holstein bulls with known conception rates [96]. In a similar study, MON1B was found to be one of twelve genes associated with embryonic development following a bovine GWAS for cryotolerance and embryonic development in cattle [97]. Thus, MON1B may be an important candidate gene that may modulate conception success. RNA Sequencing identified RAB3B as being significantly expressed in good-quality semen compared to low-quality semen following assessment of semen traits in bulls [98]. Hence, RAB3B may be a key putative candidate for modulating spermatozoal motility in sheep. Moreover, RAB3-peptides, including RAB3B, are key regulators of the acrosome reaction in sheep [99]. Given that the acrosome reaction is a necessary process for fertilisation [100], RAB3B may be an important candidate gene for promoting sperm binding to the ova and promoting conception success. In a bull study, RPL30 was abundantly expressed in bulls with high sperm motility compared to bulls with low motility [101]. Moreover, RPL30 was up-and downregulated in normozoospermic and asthenozoospermic males, respectively, compared to controls [102]. Hence, RPL30 may be a key putative candidate gene modulating spermatozoal motility. Another gene, SORD, was identified within the flagella of sperm extracted from the cauda epididymis and had increased expression during the late stages of spermatogenesis in mice [103]. In the same study, spermatozoa collected from the cauda epididymis incubated in sorbitol exhibited motility after 4 h of incubation, compared to control spermatozoa. In a sheep study, SORD was identified as a key candidate in cryotolerance following proteomic analysis of ram seminal plasma [15]. Hence, SORD may play a role in modulating ram spermatozoal motility and cryotolerance and, therefore, may influence post-thaw motility in sheep.
While there are several published GWAS on semen traits in livestock, only a single study has been published in sheep, which involves identifying QTL regions associated with semen volume, concentration, gross motility, and number of spermatozoa using semen records collected from 429 Assaf rams [6]. There is no overlap between QTLs identified in the present study and QTL regions associated with Assaf ram semen volume, concentration, and gross motility [6]. While the present study identified 35 QTLs, this was based on a relatively lenient genome-wide significance threshold (Modified Bonferroni, p-value < 0.0002). Attempting to use more stringent thresholds, e.g., based on simpleM (p-value < 0.000001) [104], resulted in only a single SNP on OAR 17 (17:22,669,197, genome build Oar_v3.1) crossing genome-wide significance for semen percent post-thaw motility, and none of the other SNPs were found to be associated with any of the other traits at the genome-wide level. This is consistent with the previously published GWAS in sheep, where only one of three semen traits had QTLs that were able to achieve genome-wide significance based on a false discovery rate of 10% [6]. In the same context, a past GWAS focused on semen traits in Chinese Holsteins also used a relatively lenient threshold when the use of a Bonferroni-based threshold did not yield significant QTLs, likely due to high type II errors [105]. The use of a relatively lenient threshold subsequently resulted in the identification of 19 SNPs that were significantly associated with Chinese Holstein bull semen traits, including volume, percent motility, concentration, number of spermatozoa, and number of motile spermatozoa.
In the present study, significant SNPs associated with semen volume, gross motility, concentration, and percent post-thaw motility accounted for 17.08%, 7.59%, 4.69%, and 28.11% of the additive genetic variance (Table 2). Previous GWAS in livestock species have also reported the proportion of genetic variance explained by QTLs found to be significantly associated with semen traits. In cattle, SNPs significantly associated with total motility were found to account for 15.98% of additive genetic variance in crossbred bulls [106], and in another study in dairy cattle, a single SNP significantly associated with semen volume was found to account for 12.22% of the additive genetic variance [107]. In pigs, a GWAS involving almost 100,000 semen collection records from 2022 pigs identified SNPs significantly associated with volume, concentration, and percent motility, accounting for 2.23% and 2.48%, 1.95%, and 6.15% of genetic variance, respectively [108]. Overall, the proportion of genetic variance explained by significant QTLs in any study is likely influenced by a variety of factors, including the genetic background of the resource population and the degree of environmental variance associated with the phenotypes in question. The semen quality data used in this study was collected over a roughly 20 year period and was hence subject to significant environmental variance, which likely contributed to the low heritability associated with these phenotypes. However, this study also involved sheep from multiple environmental locations and management practices, indicating a relatively high degree of genetic variance that likely contributed to individual SNPs explaining a proportion of the additive genetic variance that is largely consistent with previous reports in livestock species. Overall, the QTLs and associated candidate genes identified in the present study could provide insights into the underlying biological mechanisms influencing ram semen traits and possibly be used to develop future strategies to improve the reproductive efficiency of sheep.
It is worth noting that the resource population used in the current study represents a subset of animals that have previously been used to estimate genetic parameters associated with semen traits [22]. The animals included in this study represent the subset of animals for which genotype data was available. Therefore, while heritability of semen traits has been estimated in both studies, they differ in terms of the breed(s) of sheep used, the size of the resource population (n = 11,470 vs. n = 4506), as well as the method used to estimate the relationship between animals in the population (Pedigree vs. GRM). Despite these differences, the heritability estimates computed in both studies are comparable (Table 4) and are also largely in alignment with previous sheep studies [21][22][23][24]. Comparable low to moderate estimates of heritability for semen traits have also been reported in a range of livestock species such as cattle [109,110], goats [111], and pigs [112,113]. Given that the dataset was collected over a roughly 20 year period, environmental variability likely had a significant effect on the variability observed in semen traits. Furthermore, while standardised procedures for assessment of semen traits were used across the collection period, it is likely that subjectivity associated with different technicians that assessed semen during the collection period also contributed to the variability in observed phenotypes. The existence of technician bias, especially for semen traits like percent post-thaw motility, has been previously reported [114]. Overall, this would all contribute to increased environmental variance, which aligns with the low heritability estimates of all semen traits reported in this study. These results indicate that while most of the variability observed in semen traits is environmental in nature, genetic variants do contribute to this variability, and therefore, identifying such genetic variants would afford valuable insights into the physiological drivers underlying these traits, thereby contributing to the development of novel opportunities to improve reproductive efficiency.

Conclusions
Overall, 35 QTL regions located on 16 OARs were significantly associated with Merino ram semen volume, gross motility, concentration, and percent post-thaw motility. Several genes identified in the present study, including SORD, SH2B1, and NT5E, have been previously found to play crucial roles in spermatogenesis, spermatozoal motility, and high motility following cryopreservation; they may be promising candidate genes for Merino ram semen traits in future studies. Several genes like PADI2, RAB3B, and ALDOA, which have been previously associated with promoting spermatozoal maturation, the acrosome reaction, and conception success, were also identified. As such, it would be beneficial to validate such putative candidate genes via molecular analysis to characterise their influence on semen traits and possibly conception outcomes. Overall, results indicate that Merino semen traits are lowly heritable and that several QTLs influencing variability in these traits exist. Further studies are needed to characterise these QTL regions and the genes contained within them.

Supplementary Materials:
The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/ani13142286/s1, Table S1: Significance of fitted effects for each semen volume, gross motility, concentration, and percent post-thaw motility; Table S2: Uncharacterised genes associated with each QTL identified in semen volume, gross motility, concentration, and percent post-thaw motility.

Data Availability Statement:
The phenotype data present in this study is available upon request from the corresponding author, while the genotype data is confidential and covered under a data sharing agreement with Charles Sturt University (CSU) and Meat and Livestock Australia (MLA).