Effect of genotyping density on the detection of runs of homozygosity and heterozygosity in cattle

Abstract Runs of homozygosity (ROHom) are contiguous stretches of homozygous regions of the genome. In contrast, runs of heterozygosity (ROHet) are heterozygosity-rich regions. The detection of these two types of genomic regions (ROHom and ROHet) is influenced by the parameters involved in their identification and the number of available single-nucleotide polymorphisms (SNPs). The present study aimed to test the effect of chip density in detecting ROHom and ROHet in the Italian Simmental cattle breed. A sample of 897 animals were genotyped at low density (50k SNP; 397 individuals), medium density (140k SNP; 348 individuals), or high density (800k SNP; 152 individuals). The number of ROHom and ROHet per animal (nROHom and nROHet, respectively) and their average length were calculated. ROHom or ROHet shared by more than one animal and the number of times a particular SNP was inside a run were also computed (SNPROHom and SNPROHet). As the chip density increased, the nROHom increased, whereas their average length decreased. In contrast, the nROHet decreased and the average length increased as the chip density increased. The most repeated ROHom harbored no genes, whereas in the most repeated ROHet four genes (SNRPN, SNURF, UBE3A, and ATP10A) previously associated with reproductive traits were found. Across the 3 datasets, 31 SNP, located on Bos taurus autosome (BTA) 6, and 37 SNP (located on BTA21) exceeded the 99th percentile in the distribution of the SNPROHom and SNPROHet, respectively. The genomic region on BTA6 mapped the SLIT2, PACRGL, and KCNIP4 genes, whereas 19 and 18 genes were mapped on BTA16 and BTA21, respectively. Interestingly, most of genes found through the ROHet analysis were previously reported to be related to health, reproduction, and fitness traits. The results of the present study confirm that the detection of ROHom is more reliable when the chip density increases, whereas the ROHet trend seems to be the opposite. Genes and quantitative trait loci (QTL) mapped in the highlighted regions confirm that ROHet can be due to balancing selection, thus related to fitness traits, health, and reproduction, whereas ROHom are mainly involved in production traits. The results of the present study strengthened the usefulness of these parameters in analyzing the genomes of livestock and their biological meaning.


Introduction
Runs of homozygosity (ROHom) are contiguous stretches of homozygous segments within genomes, first recognized by Broman and Weber (1999) in human populations.They reflect autozygosity (McQuillan et al., 2008) because the homozy-gous segments may have been inherited from common ancestors, i.e., these segments are identical by descent (Purfield et al., 2012) and can be used as a predictor of inbreeding (Ferenčaković et al., 2013a).Their length is a temporal indicator of the inbreeding occurrence (Kirin et al., 2010): recent inbreeding results in longer ROHom, whereas ancient inbreeding is associated with shorter ROHom because recombination events break the segments over each generation.The ROHom can be due to natural and artificial selection because of the fixation of homozygous favorable alleles at selected loci; linkage disequilibrium can extend the variation of allele frequency to neighboring loci, resulting in an increase of homozygosity and the proportion of the genome covered by ROHom (Macciotta et al., 2021).For this reason, ROHom shared within a population can assist in identifying breed-specific regions potentially under selection (Mastrangelo et al., 2018a).Runs of heterozygosity (ROHet) cannot be defined as true runs but rather as heterozygosity-rich regions (Marras et al., 2018) because they are not as dense and contiguous as ROHom are.The ROHet were first introduced by Williams et al. (2016) in a study on Chillingham white cattle, which had only 9.1% polymorphic loci compared to the 62% to 90% found in commercial cattle breeds.Interestingly, these few polymorphic loci clustered in specific chromosomal regions, named ROHet blocks, that contain genes resilient to genetic drift with possible effects on fitness (Williams et al., 2016).These heterozygous regions could be associated with survival rate, fertility, and other fitness traits (McParland et al., 2009).ROHet islands can be investigated to search for balancing selection, defined as the natural selection that maintains genetic diversity via heterozygote advantage (Fijarczyk and Babik, 2015).Since the definition of ROHet is a relatively recent concept, the literature about this metric is not as abundant as for ROHom.After the first study by Williams et al. (2016), four other studies about ROHet in cattle can be found in the literature (Ferenčaković et al., 2016;Biscarini et al., 2020;Hidalgo et al., 2021;Mulim et al., 2022).Moreover, few studies have been conducted on other livestock species, e.g., in sheep (Selli et al., 2021), goats (Li et al., 2022;Chessari et al., 2024), turkey (Marras et al., 2018), horses (Santos et al., 2021(Santos et al., , 2023)), and pigs (Chen et al., 2022;Ruan et al., 2022;Bordonaro et al., 2023), with a lack of consensus in establishing the criteria to define ROHet.Therefore, one main challenge is defining the parameters (i.e., minimum run length, number of consecutive heterozygote markers, and missing single-nucleotide polymorphisms, SNPs) for their correct identification.Moreover, the density of the SNP chip used is another factor affecting autozygosity and heterozygosity estimates.Together with the lack of golden standards parameters for ROHom and ROHet, the genotyping density can introduce bias in their detection.
The ROHom tend to be more abundant in inbred and strongly selected populations (Kim et al., 2018), whereas in theory, ROHet may be more common in outbred or less selected populations (Chessari et al., 2024).The latter can be identified in local and dual-purpose populations in which the rate of artificial selection is usually lower than cosmopolitan breeds.
The Italian Simmental cattle breed represents a good livestock model to investigate both ROHom and ROHet.This breed is farmed mostly in small herds in the mountainous areas of Northeastern Italy (Cesarani et al., 2020) and it is the third largest Italian cattle breed (www.vetinfo.it, Sistema Informativo Veterinario 2022).The aim of this study was to investigate the differences in the detection of ROHom and ROHet using three different BeadChip genotyping arrays and the signals of selection highlighted by the two approaches.

Materials and Methods
Animal Care and Use Committee approval was not needed as data were obtained from preexisting databases.

Animals and genotypic data
Three datasets were used and consisted of a total of 897 Italian Simmental cattle: 397 genotypes were obtained using the Illumina BovineSNP50k array and denoted as low density (50k SNP, LD), 348 genotypes were obtained using the customized Geneseek genomic profiler GGP-HDv3 and denoted as medium-density (140k SNP, MD), and 152 genotypes were obtained using the Illumina BovineHD and denoted as high density (800k SNP, HD).After quality control carried out using PLINK v. 1.9 (Purcell et al., 2007;Chang et al., 2015), SNP were retained for the analysis if: there was no statistical deviation from the Hardy-Weinberg equilibrium (P < 1e −6 ), minor allele frequency (MAF) was > 0.01, individual animal and SNP call rates were > 95%.In addition, SNP mapped on allosomes or unmapped according to the ARS-UCD1.3were discarded.After quality control, 43,431 SNP were retained for the LD, 113,042 for the MD, and 583,637 for the HD datasets.

Detection of runs
ROHom and ROHet were identified using the "consecutive" algorithm implemented in the "detectRUNS" R package (Biscarini et al., 2019).To minimize the number of falsepositive ROHom, the minimum number of SNP to define an ROHom was computed using the following formula from Purfield et al. (2012): where n s is the number of SNP per individual, n i is the number of individuals, α is the percentage of false positive (0.05), and het is the average heterozygosity.The same equation (1) was used to compute the minimum number of SNP that constituted an ROHet, considering homozygosity instead of heterozygosity in the denominator.In order to account for the lower number of heterozygote genotypes in the genome, the number of opposite SNP allowed in a ROHet (i.e., homozygotes inside the run) was computed as where nSNP ROHom is the minimum number of SNP in an ROHom and nSNP ROHet is the minimum number of SNP in an ROHet, both computed using equation (1).From this, the following parameters were adopted: (i) ROHom = 50 homozygote SNP and no heterozygote or missing allowed within the run; (ii) ROHet = 18 heterozygote SNP, 3 opposite (i.e., homozygote), and 0 missing.
For both ROHom and ROHet, the minimum length and the maximum gap between adjacent SNP were set to 1 Mb.The number of ROHom and ROHet per animal (nROHom and nROHet, respectively) and their average length were identified.ROHom and ROHet were grouped according to their length in five different classes: 1 to 2 Mb, 2 to 4 Mb, 4 to 8 Mb, 8 to 16 Mb, and > 16 Mb.
The number of regions (nROHom and nROHet) and their average length across densities were compared using ANOVA.
The ROHom-or ROHet-based coefficients (i.e., F ROHom and D ROHet ) were computed as the ratio between the total sum of ROHom or ROHet length per animal divided by the genome length covered by SNP.As reported by Bordonaro et al. (2023), the D ROHet can be used as coefficient of diversity.The F ROHom and D ROHet values computed in the three different datasets (LD, MD, and HD) were compared using ANOVA.

Runs of homozygosity and runs of heterozygosity detected on the same animals
To avoid a possible sampling bias, for the animals genotyped at HD, the SNP in common with the LD (43,431 markers) and MD (113,042 markers) datasets were also used to compute both ROHom and ROHet.Moreover, the level of linkage disequilibrium was computed in the three densities to evaluate its potential effect on the runs detection.The squared correlation coefficient of allele frequencies at pairs of loci (r2) was estimated for all pairwise combinations of SNP between 0 and 1,000 kb apart (McKay et al., 2007) using Haploview (Barrett et al., 2005).

Repeated regions and islands
As proposed by Cesarani et al. (2018) and Macciotta et al. (2021), the identified genomic regions starting and ending at the same position found in more than one animal were regarded as repeated (ROHom REP or ROHet REP ) and those in the top 0.1% of the distribution of animals sharing the repeated region were of interest.Finally, the number of animals with an SNP in a run was computed as SNP ROHom and SNP ROHet for the homozygous and heterozygous runs, respectively.The markers exceeding the 99 th percentile of these distributions (i.e., SNP ROHom and SNP ROHet values) were considered as ROHom or ROHet islands.These two sets of values were compared using ANOVA.

Gene and quantitative trait loci enrichment
Using the NCBI online database (National Center for Biotechnology Information, www.ncbi.nlm.nih.gov), the genes mapped in or close (± 250 kb; Manca et al., 2020) to the most repeated ROHom and ROHet and to the highlighted islands were identified and described according to the available literature.Moreover, for the same regions, the quantitative trait loci (QTL) were identified using the GALLO R package (Fonseca et al., 2020a), which was also used to carry out an enrichment analysis using Bonferroni correction for the P-value (0.05).

Runs of homozygosity
Results of ROHom, in terms of number of regions and their average length, are reported in Table 1.The number of ROHom per animal showed a non-linear trend: the highest value (47.61 ± 19.80) was observed in the MD and the lowest (16.00 ± 6.52) in LD.On the contrary, the average length decreased by about four times (from LD to HD), as the density of genotyping increased.The same was true for the distribution of ROHom in the different length classes (Table 1): in the  Different superscript letters within column within ROHom or ROHet indicate significant differences for P < 0.001.

Runs of heterozygosity
In contrast to ROHom detection, the largest nROHet was found in the LD dataset (3.71 ± 1.89), whereas the lowest average nROHet was found in MD (1.53 ± 0.72).The average ROHet length increased (from 1.29 to 1.45 Mb) moving from LD to MD.In the HD dataset, no relevant ROHet were found (Table 1).The ROHet were found only in the first two classes of length (1 to 2 and 2 to 4 Mb), with only 3% (LD) and 18% (MD) of the identified ROHet found in the 2 to 4 Mb length class.As reported in Table 1, the BeadChip density had an impact on the average length of the ROHet, which significantly (P < 0.001) increased moving from LD (1.29 ± 0.88 Mb) to MD (1.45 ± 0.48 Mb).
The D ROHet values were very low (all below 1%): the maximum D ROHet observed in LD was 0.005, whereas the maximum value (0.003) was estimated in MD.As shown in Figure 1, the values computed for the LD dataset were significantly larger (P < 0.001) than those computed in the MD dataset.

Runs of homozygosity and runs of heterozygosity detected on the same animals
Results of ROHom and ROHet identified in the 152 animals genotyped with the three densities were reported in Table 2.These results agree, in terms of both the average number of regions per animal and their length, with those observed in the other two groups of animals (Table 1).Moreover, to test a possible effect of the differences in the linkage disequilibrium among the three datasets (LD, MD, and HD), the values of r 2 were compared among different densities (Supplementary Figure S1): no differences were observed in the linkage disequilibrium decay.

Repeated regions
In this study, a total of 21 regions exceeding the 0.1% of the ROHom REP distribution and two regions exceeding the 0.1% ROHet REP distribution were found (Table 3).Among the top 21 ROHom REP , four were detected in the LD, 13 in the MD, and four in the HD dataset, respectively.The most shared ROHom was found in 27 different animals (~8% of the animals genotyped with the MD BeadChip), and it was located at 38.43 to 39.46 Mb on BTA6 (Table 3), in which no genes were mapped.However, this genomic region overlapped with 414 QTL (Supplementary Table S1).According to the enrichment analysis, a total of 27 different terms were highlighted, of which 13 were significant (Supplementary Table S1).The significant terms were associated with exterior (22 QTL), meat and carcass (114), production (209), and reproduction (4) traits.
The most repeated ROHet was found in 41 animals (~10% of the animals in the LD dataset) and it was located at 2.15 to 3.25 Mb on BTA21 (Table 4), where four genes (SNRPN, SNURF, UBE3A, and ATP10A) were mapped (± 250 kb downstream and upstream from the repeated ROHet).The same region overlapped with 762 QTL (Supplementary Table S1), which were grouped in 14 terms by the enrichment analysis.Among them, only one, calving ease, was significant (738 QTL)., respectively.Across datasets, the average probability of having a particular SNP inside a ROHom (3.90 ± 2.56%) was significantly higher (P < 0.001) compared to ROHet (0.09 ± 0.58%).Thus, an SNP had a larger likelihood to be within an ROHom rather than ROHet.The SNP ROHom values computed in the three datasets were moderately and positively correlated each other: 0.51 (LD vs MD), 0.42 (LD vs HD), and 0.53 (MD vs HD).The correlation between SNP RO- Het computed in the LD and MD datasets was lower (0.33).

ROHom and ROHet islands
A total of 419, 1,102, and 5,626 SNP exceeded the 1% of the SNP ROHom distribution in the LD, MD, and HD datasets, respectively: 31 SNP were in common among the three different datasets.These markers mapped on BTA6 between 38.20 and 40.63 Mb and could be considered as one large ROHom island (Table 4).It is interesting to note that this island included the most repeated ROHom (38.43 to 39.46 Mb, Table 3).The genomic region identified by the 31 SNP exceeding the 99th percentile of the SNP ROHom mapped three genes (SLIT2, PACRGL, and KCNIP4) and overlapped with 1,750 QTL.The latter were divided into 61 terms, of which 15 were significant.The significantly enriched terms were related to exterior (22 QTL), meat and carcass (461), production (1,113), and reproduction (4).
A total of 409 (in LD dataset) and 642 SNP (in MD dataset) exceeded the SNP ROHet distribution.There were 37 SNP ROHet in common among the two datasets, of which 19 were mapped on BTA16 and 18 on BTA21 (Table 4).On BTA16, four different ROHet islands were identified; in these regions, a total of 39 genes were present (Supplementary Table S2).
On BTA21, two different islands were highlighted.The longest ROHet island (17 SNP) completely overlapped with the most shared ROHet mentioned above (Table 3).The second peak, highlighted by a single SNP, had one gene, the gamma-aminobutyric acid type A receptor subunit beta3 (GABRB3).The ROHet islands found in the two chromosomes (BTA16 and 21) partially or totally overlapped with 819 QTL (Supplementary Table S1) that were divided in 24 enriched terms.Among them, only 2 were significant: calving ease (738 QTL) and interval to first estrus after calving (20), both in the reproduction group.

Discussion
The availability of high-throughput sequencing or genotyping SNPs data has opened the possibility of characterizing individual segments of the genome in terms of their homozygosity (regions of homozygosity (ROH)) or heterozygosity  Different superscript letters within column within ROHom or ROHet indicate significant differences for P < 0.001.
(heterozygosity-rich regions) (Biscarini et al., 2020).These genomic regions are both a consequence of the selection that shapes the genomic structure of livestock populations (Bordonaro et al., 2023).The occurrence of ROHet avoids the deleterious effects of continuous homozygous genotypes and favors the heterozygote advantage in immune-related genes as well as in productive and reproductive traits (Sanglard et al., 2021;Chen et al., 2022;Chessari et al., 2024).
On the contrary, the study of ROH distribution has been a helpful tool to detect regions potentially under selection (e.g., Mastrangelo et al., 2017;Cesarani et al., 2021;Macciotta et al., 2021).Despite such advancements, both ROHom and ROHet detection, are sensitive to various parameters, such as the genotyping density, and so far, only few studies addressed this issue in ROHom (Purfield et al., 2012;Ferenčaković et al., 2013b;Hillestad et al., 2018) or ROHet (Mulim et al., 2022).Therefore, in this study, the effect of genotyping density on ROHom and ROHet features was tested, and the signals highlighted across datasets were further investigated.Genes and associations with QTL were searched within ± 250 kb downstream and upstream of the reported positions.
2 Details and references listed in Supplementary Tables S1 and S2.

Runs of homozygosity
The nROHom decreased and their average length increased as the genotyping density increased.This is because a higher SNP density improves the discovery resolution, by reducing the detection of false long ROHom.However, for the same reason, the use of the HD panel could lead to an underestimation of number of segments longer than 8 Mb (Ferenčaković et al., 2013b).This was confirmed by the distribution of ROHom in the different length classes (Table 1).In fact, when a denser SNP panel is used, a larger number of opposite markers (i.e., heterozygote for ROHom and homozygote for ROHet, respectively) can break a long region to shorter ones, thus reducing the number of long regions (Hillestad et al., 2018).On the contrary, Purfield et al. (2012) reported that using an LD BeadChip could lead to an overestimation of regions shorter than 4 Mb.Also, Ferenčaković et al. (2013b) showed that the 50k panel revealed an abundance of small segments and overestimated the numbers of segments 1 to 4 Mb long, suggesting that it is not sensitive enough for the precise determination of small segments.Thus, the LD and HD BeadChips are not able to precisely identify small or long segments, respectively.To evaluate the detection power of the three investigated datasets, the ROHom statistics were compared with the values reported in the literature for the same breed, to avoid bias due to different breeding management, selection pressure, or evolution history since ROHom patterns can be associated with these phenomena (Purfield et al., 2012;Gaspa et al., 2014;Forutan et al., 2018;Macciotta et al., 2021).Most of the available literature on ROHom detection in this breed involved the use of a 50k SNP, which correspond to the LD dataset in the present study.Most of these studies used an arbitrary fixed number of markers (i.e., 15 SNP) and 1 Mb as minimum values to define a ROHom.For example, Cesarani et al. (2021) and Marras et al. (2015) identified the ROHom in Italian Simmental, and they reported the nROHom as 77.47 ± 15.23 and 94.30 ± 12.20, with an average length of 2.45 ± 3.24 and 2.2 Mb, respectively.Ferenčaković et al. (2011) reported similar values (i.e., nROHom = 96.79± 13.37, with an average length of 2.4 Mb) in Austrian Simmental bulls.Finally, Szmatola et al. (2016) found nROHom to be 81.5 ± 11.8 with an average length of 2.49 Mb in Polish Simmental.Zhao et al. (2021) used a HD BeadChip array to investigate ROHom in Chinese Simmental.These authors defined a ROH with at least 100 SNP covering 500 kb and they allowed two homozygous and one missing genotypes.The average nROHom was 99.03, with an average length of 1.18 Mb.All these values, especially regarding nROHom, are greater than those reported here, probably because of the stricter parameters adopted in the present study for ROHom identification (e.g., 50 as minimum number of SNP, and no heterozygotes or missing markers).The MD dataset showed the greatest nRO-Hom (with a standard deviation lower than that reported for HD) and an average size (2.72 ± 2.61 Mb) similar to the values reported in the literature.
The F ROHom values computed in the present study were compared with the estimates reported in literature to investigate which density led to more consistent results.In literature, different F ROHom values are reported for the Simmental breed: 0.08 ± 0.04 (Szmatola et al., 2016), 0.09 ± 0.02 (Ferenčaković et al., 2011), 0.07 ± 0.03 (Cesarani et al., 2021), and 0.08 (Marras et al., 2015).However, values comparable to those estimated in this study were reported in Italian Simmental by Mastrangelo et al. (2018b), with a mean F ROHom of 0.03 ± 0.02.In the present study, the values closer to the values found in the literature were computed using the MD dataset.The lower values here computed could be ascribed to the old age of genotyped animals, especially in the HD.In animal populations under genetic improvement, old animals usually show inbreeding coeffi-cients lower than young animals (Makanjuola et al., 2020;Guinan et al., 2023).Moreover, the number and length of ROHom tend to increase over time: Forutan et al. (2018) analyzed ROHom in North American Holstein cattle and found that the rate of increase of ROHom longer than 1 Mb in the last 5 yr was almost double that of the previous 5-yr period.
Since the LD and HD datasets could lead to bias in the number of detected ROHom and the ROHom statistics (in terms of both average length and inbreeding coefficient) obtained in the MD are closer to the literature, the latter density could represent the best option to detect ROHom in the Simmental breed, also because MD BeadChip is cheaper than HD.

Runs of heterozygosity
The ROHet identified in the present study were fewer and shorter compared to ROHom.In particular, only ROHet shorter than 4 Mb were found.A similar result was reported in a study on horses, where only ROHet shorter than 2 Mb were identified (Santos et al., 2021).Recently, Chessari et al. (2024) in a study on goats showed an average ROHet length < 1Mb.Moreover, Biscarini et al. (2020) found just two ROHet longer than 2 Mb in Maremmana semi-feral cattle.Ruan et al. (2022) analyzed the distribution of ROHet in two Duroc pig populations and they found only about 3%-5% of regions in the length class > 4 Mb.As already observed for ROHom, the BeadChip density had an impact on both nROHet and their average length, which increased from LD to MD.This result is in agreement with Mulim et al. (2022), who also found shorter ROHet using lower density.Since no ROHet were detected in the HD dataset and the use of LD could overestimate long runs, the MD could represent the optimum array for the ROHet detection, as already pointed out for ROHom.ROHet statistics were not compared with values from the literature, because no studies about ROHet on Simmental cattle were found.The low D ROHet values computed in the present study are due to the low number of heterozygous regions and their short length.Very low ROHet-based coefficients agree with the few reports available in literature.For example, Bordonaro et al. (2023) analyzed diversity indices estimated from ROHet in pigs and they found average values ranging from 0.0001 to 0.0047.A recent study (Chessari et al., 2024) on ROHet in Italian goat populations reported low values of a similar magnitude (0.0024 ± 0.0003).

Runs of homozygosity and runs of heterozygosity detected on the same animals
Since the results on ROHom and ROHet were obtained on different datasets of animals, results reported in the present study for different SNP densities may be affected by a sampling bias.To check this hypothesis, analyses were repeated using only the 152 animals genotyped with 800k SNP: LD and MD densities were then mimicked by retaining only the markers included in the lower densities.Results on this subset of animals, in terms of the number of regions and length, were similar to those obtained using the three different datasets.As expected, no differences were observed in the linkage disequilibrium among densities since all animals belonged to the same breed.

Repeated regions
A genomic region shared among different animals of the same breed could be associated with a selection pressure on portions of the genome that control economically important traits as well as other important animal characteristics such as disease resistance or general immune competence.In particular, if the shared genomic region is characterized by a high level of homozygosity, it could be due to directional selection, both artificial or natural (Kim et al., 2013;Gorssen et al., 2021), whereas a high level of heterozygosity could be associated to balancing selection (Fijarczyk and Babik, 2015).
In the present study, the most shared ROHom was located on BTA6, where no genes were mapped.For the same genomic region, several associations with QTL were found in the enrichment analysis (Supplementary Table S1).The discovery of QTL associated with production and meat and carcass traits was expected due to the breeding goals of this breed.
In the most repeated ROHet, on BTA21, four genes (SNRPN, SNURF, UBE3A, and ATP10A) were mapped.The ubiquitin protein ligase E3A (UBE3A) gene was previously found to be associated with stillbirth and calving ease in cattle (Mészáros et al., 2016).The ATPase phospholipid transporting 10A, putative (ATP10A) has been associated with calving ease (Frischknecht et al., 2017) and milking speed (Marete et al., 2018).UBE3A and SNRPN (small nuclear ribonucleoprotein polypeptide N) have been associated with cattle temperament by Costilla et al. (2020) and, together with the SNURF (SNRPN upstream open reading frame) gene, with the age at first calving (Alves et al., 2022).As Suzuki et al. (2009) pointed out, the SNRPN and SNURF constituted a bicistronic gene (SNRPN-SNURF), which has been extensively studied in mice and humans and has been associated with neurodevelopmental disorders.Finally, all four genes were found related to the occurrence of early pregnancy in Nellore cattle (Irano et al., 2016).Also the QTL overlapping with this genomic region were mainly associated with reproduction traits and, in particular, with calving ease.Heterozygosity and, thus, ROHet could be mainly associated with balancing selection rather than directional selection.Indeed, genes and QTL found to be associated with this metric were mostly related to functional phenotypes.

ROHom and ROHet islands
As shown in Figures 2 and 3, an SNP is more likely to be in an ROHom than an ROHet.This was expected because of the larger number of homozygous genotypes along the genome, which was reflected in the results reported in Table 1.The ROHom island identified by the SNPROHom exceeding the 99th percentile of the distribution in all the three datasets was located on BTA6 and contained three genes.The slit guidance ligand 2 (SLIT2) gene has been previously reported to be involved in several weight traits: in particular, internal organ (especially spleen) weight in Simmental cattle (An et al., 2018), bone weight in beef cattle (Niu et al., 2021), birth, yearling, and weaning weights in US Red Angus cattle (Smith et al., 2022), and birth weight in US Gelbvieh cattle (Smith et al., 2019).Moreover, the same gene was also associated with the infection of tropical theileriosis parasite (Larcombe et al., 2022) and with female fertility in Nordic Red cattle (Höglund et al., 2015).The parkin coregulated like (PACRGL) has been associated with height and stature of cattle (Doyle et al., 2020).The potassium voltage-gated channel interacting protein 4 (KCNIP4) was reported to be related to milk fat percentage (Pedrosa et al., 2021), fertility (Tarekegn et al., 2021), birth weight and yearling weight (Smith et al., 2022), and to backfat thickness and carcass weight (Srikanth et al., 2020).Moreover, all the three genes were reported to be associated with clinical or subclinical ketosis by Soares et al. (2021).The majority of QTL flagged by this genomic region (Supplementary Table S1) was significantly enriched in two main categories: production-average daily gain, body weight and body weight gain, dry matter intake, metabolic body weight-and meat and carcass-biceps brachii weight, bone weight, carcass weight, lean meat yield, liver weight, longissimus muscle area, and subcutaneous fat thickness.As expected, the ROHom islands harbored genes and QTLs with similar functions than those found in the most repeated ROHom regions.The relationship with meat traits was expected because they have a weight of 24% in the breeding program of Italian Simmental (Cesarani et al., 2020).At the same time, the non-significance of QTLs associated with milk production traits (N = 104, Supplementary Table S1) was quite surprising since the latter has a larger weight (44%) in the breeding program of this breed.

FIGURE 3 .
FIGURE 3. Manhattan plot of the number of times a particular SNP was inside a run of heterozygosity (SNP ROHet ) detected in the three datasets.Green dots represent SNP with a SNP ROHet value exceeding the 99th percentile of its distribution.

TABLE 1 .
Runs of homozygosity (ROHom) and runs of heterozygosity (ROHet) detected with three densities

TABLE 2 .
Runs of homozygosity (ROHom) and runs of heterozygosity (ROHet) detected with three densities on the same 152 animals

TABLE 4 .
SNP exceeding the 99th percentile of the number of times a particular SNP was inside a run of homozygosity (SNP ROHom ) or heterozygosity (SNP ROHet ) distribution in each density