Runs of Homozygosity and Quantitative Trait Locus/Association for Semen Parameters in Selected Chinese and South African Beef Cattle

Simple Summary Genome-wide runs of homozygosity (ROH) are excellent in understanding population history, estimating genomic inbreeding, and deciphering genetic architecture of complex traits and diseases, as well as identifying genes of agro-economic traits. ROH are defined as continuous region on the chromosome where an individual is homozygous across the genome. This study identified distribution of ROH in the six selected beef cattle breeds, Chinese Simentaler (CSI), Belgian Blue (BEL), and South African Angus (ANG), Nguni (NGU), Bonsmara (BON), and Simentaler (SIM), using Bovine BeadChip markers. Furthermore, nine candidate genes, CDF9, MARCH1, WDR19, SLOICI, ST7, DOP1B, CFAF9, INHBA, and ADAMTS1, were suggested to be associated with semen QTL traits and reported moderate inbreeding in some breeds with high to low correlation inbreeding between breeds. The study findings will allow proper guidelines for breeder’s societies. Abstract In this study, runs of homozygosity (ROH) and quantitative trait locus/association (QTL) for semen parameters in selected Chinese and South African beef cattle breed were estimated. The computed results showed 7516 ROH were observed between classes 0–5 Mb with no ROH observed in classes >40 Mb. Distribution of ROH showed high level of genomic coverage for ANG, NGU, CSI, and BEL breeds. Approximately 13 genomic regions with QTL were controlling sperm motility, sperm concentration, semen volume, sperm count, sperm head abnormalities, sperm tail abnormalities, sperm integrity, and percentage of abnormal sperm traits. Nine candidate genes, CDF9, MARCH1, WDR19, SLOICI, ST7, DOP1B, CFAF9, INHBA, and ADAMTS1, were suggested to be associated with above mentioned QTL traits. The results for inbreeding coefficient showed moderate correlation between FROH vs FHOM at 0.603 and high correlation between FROH 0–5 Mb 0.929, and lowest correlation for 0–>40 Mb 0.400. This study suggested recent inbreeding in CSI, BEL, ANG, BON, SIM, and NGU breeds. Furthermore, it highlighted varied inbreeding levels and identified QTL for semen traits and genes of association. These results can assist in implementation of genetic improvement strategies for bulls and provide awareness and proper guidelines in developing breeding programs.


Introduction
Breeding programs for cattle have been known to be widely implemented in a lot of countries in the world. Several countries have regarded breeding as an important aspect in beef cattle, while some countries have already implemented those programs. These programs include breeds like South African Bonsmara, an indigenous cross between European Shorthorn, Hereford, and Afrikaner [1], and the indigenous Southern African (South Africa, Eswatini, Namibia, Zimbabwe, Botswana, and Angola) Nguni breed a cross between Bos taurus and Bos indicus. The Simmentaler, a cross between German cattle and Swiss cattle is still one of the most utilized breeds in breeding programs. The Scottish Angus [2] has been in South Africa over a hundred years, even part of the South African Studbook. The Belgian blue originally from Belgium and is the first known Belgian breeds [3]. These programs lead to fast genetic progress, but they also lead to the accumulation of inbreeding via heavy impact of a few selected individuals [4]. Therefore, using genome-wide runs of homozygosity (ROH) is excellent in understanding population history and estimating genomic inbreeding, solving genetic architecture of complex traits and diseases, as well as identifying genes of agro-economic traits [5]. ROH can be defined as continuous regions on the chromosome where an individual is homozygous across the genome [5]. They occur when both haplotypes transmitted from parents are identical and inherited from a common ancestry [5]. An increase in homozygous loci and regions of homozygosity may be an indication of loss in genetic diversity [6], or population going through bottleneck.
A few studies have reported the negative impact of high homozygosity on fertility traits, including bull semen quality [7], calving rate [8], stillbirths, and dystocia [9]. A previous study [10] showed distinct differences in the length, quantity, and frequency of ROH between 11 Polish breeds as well as highlighting level of genomic inbreeding within breeds. The authors of [7] reported that longer ROH (>100 Mb) indicate recent inbreeding that has occurred within a population; however, shorter ROH may indicate ancient inbreeding that has occurred in population. Long ROH may persist in out bred individuals, due to unusual mutation, linkage disequilibrium (LD), and recombination rates at certain genomic locations [8]. The investigation of ROH in farm animals also suggested their importance contributions to complex traits [6]. Runs of homozygosity in American Holstein cattle revealed candidate genes associated with reproduction traits affecting fertility [11]. There are several factors which affects the quality of ROH calling rate; these include the marker density, their distribution across the genome, the quality of the genotype calling rates, and minor allele frequency [12]. Single nucleotide polymorphism provide information about both past and more recent demographic variations of a population [12], allowing a comparison of the degree of homozygosity among populations with varying degrees of isolation and inbreeding [5]. The authors [13] stated different approaches to estimating inbreeding at an individual and population level using genetic information. The objective of the study is to use distribution of ROH to identify inbred individuals in both South African and Chinese beef bulls and identify QTL for semen traits and their associated genes. This will assist in understanding population history and reveal genes linked to with semen traits. It will further show levels of inbreeding in commercial breeds, especially beef breeds [10].

Estimation of Runs of Homozygosity (ROH)
Runs of homozygosity were estimated using PLINK v1.07 [14] following the criteria: (1) the minimum length was 1000 kb; (2) the proportion of homozygous overlap window was 0.05; (3) the minimum number of SNPs included in a ROH was 100; (4) the minimum SNP density was set to 50 kb/SNP; (5) the maximum gap between continuous homozygous SNPs was 1000 kb; and (6) a maximum of one SNPs with missing genotypes and up to one heterozygous genotype were allowed in a ROH. RStudio software was used, package "detectRUNS" version 0.9.6 was used to summarize ROH into five classes: (0-5 Mb), (5-10 Mb), (10)(11)(12)(13)(14)(15)(16)(17)(18)(19)(20), , and >40 for all six breeds. Runs of homozygosity were calculated per breed and ROH coverage in each chromosome was estimated as the sum of the total length of the chromosome covered by ROH of all individuals in a population. Furthermore, the percentage of SNPs present in ROH were calculated based on the frequency of a SNP in a ROH across individuals.

QTL and Genomic Regions in ROH
PLINK v1.07 [14] (Cambridge, MA, USA) was used to estimate the consensus regions across individuals, which represent ROH pools of overlapping and potentially matching segments. The software Bovine Quantitative Trait Locus (QTL) Analysis (QTL) [15] https: //www.animalgenome.org.cgi-bin/qtldb/btsearch (accessed on 28 January 2022) [16] was used for information, gene information, and animal information. Then, National Center for Biotechnology Information (NCBI) https://www.ncbi.nlm.nih.gov (accessed on 27 March 2022) was used to confirm ID of the genes found from published QTL/association and functions. We then used KEGG [17] https://www.genome.jp/kegg-bin/showorganism?org= bta (accessed on 28 January March 2022) to identify the position and length of the genes.

Classification of Runs and Inbreeding Co-Efficient
Measuring homozygosity per individual was calculated following the method of [18].
LROH is the total length of ROH and LAUTO is the length of the autosomal genome [19]. The ROH length categories were as follows: >0 Mb, >5 Mb, >10 Mb, >20 Mb, and >40 Mb. Furthermore, inbreeding based on homozygous SNPs was determined using PLINK v1.07 software [14]. The inbreeding coefficient for an individual (FHOM) was computed. Correlations of the inbreeding coefficient for two methods were estimated using the Pearson correlation from Minitab ™ version 17.1.0 [20] (State College, PA, USA).

The proportion of SNPs in ROH
A Manhattan plot was used to compute the number of significant SNP in a ROH (Figure 4). The results revealed SNP in ROH across all individual breeds ( Figure 4A) with BEL showing significant SNP in ROH on BTA 3, 4, and 14. (Figure 4B), CSI on BTA 6 and 14, ( Figure 4C); for BON on BTA 6 and 14, (Figure 4F) ANG on BTA 1, 3, 5, 7, 13, 14, and 20, and ( Figure 4E) for SIM, significant SNP were detected on BTA 6 and 14, (Figure 4D), with no significant SNP for NGU breeds. The results further highlighted most significant SNP were higher in all breeds on BTA 14, excluding NGU and SIM. The South African ANG had the most significant SNP in ROH based on BTA coverage. Further results supporting Manhattan plot are shown on (Table 1).    Runs of homozygosity length segments were used to detect QTL for semen traits using cattle QTL online database. QTL were identified for sperm motility, sperm concentration, semen volume, sperm count, sperm head abnormalities, sperm tail abnormalities, sperm integrity, and percentage of abnormal sperm, with the majority QTL suggested for sperm motility accounting for 80% of the results (Table 2). On BTA 6, the study revealed suggestive QTL for sperm motility, sperm head abnormalities, and sperm tail abnormalities located on the WDR 19 gene. MARCH1 gene also shared an association on BTA 6 to suggestive QTL for sperm count, semen volume, and sperm motility, with the majority of QTL on BTA 6 reported for Chinese CSI and South African SIM breeds covering 68 Mb of the genome position. Secondly on BTA 1 suggestive QTL for sperm concentration revealed an association to CFAP9 gene and QTL for sperm motility and sperm count revealed an association to DOP1B gene, with suggestive QTL for sperm motility revealing an associated with ADAMTS1 and CRYZL1 gene. Other suggestive QTL for ANG breed were observed on BTA 7 for sperm acrosome integrity rate, sperm motility and sperm count which revealed an association to GDF9 gene. However, some of the suggestive QTL for percentage of normal sperm on BTA 13 and sperm motility on BTA 3 showed no association to any gene for ANG breeds. The South African BON on BTA 5 revealed QTLs for sperm motility and it's associated with SLCO1C1 gene and BEL on BTA 5 suggested QTL for sperm motility and revealed an association with SLCO1C1 and WDR 19 genes, also INHBA and ST7 genes on BTA 4 found an associated to sperm count. The results reported 13 genomic regions had identified QTL and their association; additionally, BTA 14 had the greatest coverage of SNP in ROH in CSI, BEL, BON, and ANG, but there were no detection of QTL or association on those regions.

Inbreeding Coefficient FROH vs. FHOM and ROH Regions
Lengths of ROH were further used to estimate inbreeding coefficient (FROH), shown on ( Figure 5); inbreeding coefficient (FROH) was highest in NGU (0.16), followed by ANG (0.14), BEL (0.13), CSI (0.12), BON (0.11), and, finally, SIM (0.08). The Pearson correlation method was used to estimate the linear correlation coefficient between five classes and at a genomic level shown on (Table 3). Correlation estimated for FROH vs FHOM was 0.603, showing a significant moderate correlation. However, the highest correlation was observed between Froh _5 vs. Froh _0 at (0.92), and the lowest correlation was observed between Froh _40 vs. Froh _5 (0.40). This shows that the most correlation was observed the first classes (0-5 Mb).

Inbreeding Coefficient FROH vs. FHOM and ROH Regions
Lengths of ROH were further used to estimate inbreeding coefficient (FROH), shown on ( Figure 5); inbreeding coefficient (FROH) was highest in NGU (0.16), followed by ANG (0.14), BEL (0.13), CSI (0.12), BON (0.11), and, finally, SIM (0.08). The Pearson correlation method was used to estimate the linear correlation coefficient between five classes and at a genomic level shown on (Table 3). Correlation estimated for FROH vs FHOM was 0.603, showing a significant moderate correlation. However, the highest correlation was observed between Froh _5 vs. Froh _0 at (0.92), and the lowest correlation was observed between Froh _40 vs. Froh _5 (0.40). This shows that the most correlation was observed the first classes (0-5 Mb).

Discussion
Inbreeding is known to negatively affects the reproductive performance of male animals [21], and when expressed at high levels can cause poor semen quality [22], [23], [24]- [25]. However, not all inbreeding reported is harmful; the authors of [26] reported that recent inbreeding is more harmful than ancient inbreeding due to selection decreasing the frequency of deleterious alleles over time. Many studies have explored the genome wide distribution of ROH and inbreeding depression in cattle populations [9] [27], [28]using high density Illumina BovineHD BeadChip microarrays. This includes studies such as [29], whose findings showed that ROH are frequent across all breeds (Angus, Belgian Blue, Charolais, Friesian, Hereford, Holstein, Holstein-Friesian, Limousin, and Simmental);

Discussion
Inbreeding is known to negatively affects the reproductive performance of male animals [21], and when expressed at high levels can cause poor semen quality [22][23][24][25]. However, not all inbreeding reported is harmful; the authors of [26] reported that recent inbreeding is more harmful than ancient inbreeding due to selection decreasing the frequency of deleterious alleles over time. Many studies have explored the genome wide distribution of ROH and inbreeding depression in cattle populations [9,27,28] using high density Illumina BovineHD BeadChip microarrays. This includes studies such as [29], whose findings showed that ROH are frequent across all breeds (Angus, Belgian Blue, Charolais, Friesian, Hereford, Holstein, Holstein-Friesian, Limousin, and Simmental); however, difference in patterns of ROH and variation is based on breed origin and recent management. The authors of [10] reported the distinct differences in length, quantity, and frequency of ROH between breeds and levels of genomic inbreeding. The study highlighted higher levels of inbreeding observed in commercial breeds, especially beef breeds and identified a number of genes confirmed to influence production traits [10]. Homozygous regions of the genome [8] have been reported in different species, not only in cattle [19,30], but these regions have been used to quantify individual inbreeding in humans [18], goats [31], buffalo [32], sheep [33], and pigs [12]; based on results of ROH in these studies, it appears to be more accurate than traditional pedigree-based estimates [13]. In this study, more than 1000 Mb ROH were estimated in six selected Chinese and South African beef cattle breeds at a genomic length of 2.6 Gb, with the highest ROH observed in the NGU, ANG, BEL, CSI, SIM, and BON breeds. Furthermore, the highest ROH were observed between classes 0-5 Mb compared to the rest of the classes (5-10 Mb, [10][11][12][13][14][15][16][17][18][19][20] with no ROH observed in >40 class. Similar results were reported by [29], who compared European cattle breeds to the European bison, with average ROH observed for Angus and Hereford, and the most ROH were observed for the European bison breeds between classes 1-5 Mb, suggesting limited recent inbreeding for Angus and Hereford, but reporting high levels of inbreeding for European bison. These results confirm earlier studies that reported population history involving a severe bottleneck [29]. Runs of homozygosity can be affected by demographic events [34], e.g., age, geography, breed history and origin, area of distribution, climate, introduction of breeds, environment, production system, etc.; the fact that breeds from this study are distributed in two geographically different environments might explain the inbreeding. The authors of [35] reported that inbred individuals are particularly sensitive to environment changes this can explaining most of the breeds in the study are inbred due to moving from one captive environment to another captive environment. Analysis at a genomic level revealed significant distribution of ROH within the selected six Chinese and South African beef cattle breeds. Chinese CSI and BEL had significant ROH observed on chromosomes 5, 6, 11, and 14 for both breeds. South Africa BON showed the highest ROH on chromosomes 5, 7, and 14, and SIM displayed the most ROH on chromosomes 5, 6, and 14. The highest genome coverage in ROH observed for ANG on chromosomes 1-6, 11, and 14, and NGU on chromosomes 1-3, 5-7, 11, and 14. The study further used length and number of ROH to identify SNPs in those ROH. The authors of [30] reported that greatest number of ROH per chromosomes was observed for chromosome 1 across all 867 animals, and ROH per chromosome tended to decrease with chromosome length. By identifying the number of SNPs in ROH, the study identified overlapping ROH that were above the threshold of 0.05. Chromosome 14 had the most SNPs above the threshold; still, it did not report any QTL or association to any candidate gene in all breeds. Additionally, no SNP in ROH were identified for NGU; possibly markers were standardized for commercial breeds, e.g., (Simentaler) however Nguni is an indigenous Southern African breed. In total, we reported 13 genomic regions that were associated to nine candidate genes (CDF9, MARCH1, WDR19, SLOICI, ST7, DOP1B, CFAF9, INHBA, and ADAMTS1). The authors of [36] reported that genomic regions located on BTA 7, BTA 14, BTA 16, and BTA 18 were characterized by a high frequency of ROH occurrence and included important genes related to immune traits, muscularity, and ease of calving. The authors of [10] reported a similar region, BTA 6, was identified for Polish Red, Limousin, and Simmental breeds. MARCH1 [19] revealed an association with QTL for sperm motility and semen volume on BTA 6 in SIM. The MARCH 1 gene was also reported by [37] to be significantly associated with semen production traits (semen volume per ejaculate, number of sperm per ejaculate, and number of motile sperm per ejaculate). Another gene associated with QTL for sperm motility revealed on BTA 5 is SLCO1C1 [17] for BON and BEL. BTA 7 revealed the growth differentiation GDF9 gene was associated with QTL for sperm motility, sperm count, and sperm integrity. This was also reported by [38] who revealed significant association of GDF9 with sperm quality traits in Holstein bulls. The results also stated that GDF9 is involved in the initiation or maintenance of spermatogenesis; however, further verification is needed. WD repeat domain 19 (WDR 19) gene was detected on BTA 5 and BTA 6 for CSI showed an association to QTL sperm motility, sperm head abnormalities, and sperm tail abnormalities; this gene was reported to be associated to ejaculate volume, sperm concentration, sperm motility, and sperm head and tail anomalies [39]. On BTA 4, ST7 gene was significant to QTLs for sperm count and sperm motility, this gene was also reported by [40] as the nearest genes for poor sperm motility. QTL associated with semen traits were identified and validated by previously published literature. This study highlighted genes within QTL regions of semen for beef bulls. Identification of the QTL regions associated with these traits provides the knowledge necessary to enrich these regions [35]. The authors of [41] reported on three QTL found to be related with abnormal sperm frequencies at a significant p < 0.01. The authors of [40] reported several candidate genes associated with sperm concentration, sperm motility and sperm volume in Holstein-Friesian populations.
The study further showed high to moderate correlation between classes Froh_0 −5 Mb and correlation at Froh_ > 40. Moderate correlation was observed between Froh vs Fhom selected breeds. Several authors reported on strong correction between inbreeding coefficient [9]. It should be underlined that the occurrence of ROH in an individual may be the result of inbreeding events, but they may also be present in outbreed populations as result of other phenomena.

Conclusions
The study identified distribution of ROH in the six selected Chinese and south African beef cattle breeds (CSI, BEL, ANG, NGU, BON, and SIM) using Bovine HD Bead Chip makers. The study also showed that ANG, NGU, and CSI showed the highest frequency and length of long ROH (0-5 Mb), indicating higher recent inbreeding in all breeds, including SIM. Several QTL and genes were related to semen traits (CDF9, MARCH1, WDR19, SLOICI, ST7, DOP1B, CFAF9, INHBA, and ADAMTS1). These genes can be used as target genes for future marker-assisted selection.