http://orcid.org/0000-0003-3286-2441
Figures
Abstract
The analysis of runs of homozygosity (ROH), using high throughput genomic data, has become a valuable and frequently used methodology to characterize the genomic and inbreeding variation of livestock and wildlife animal populations. However, this methodology has been scarcely used in highly inbred domestic animals. Here, we analyzed and characterized the occurrence of ROH fragments in highly inbred (HI; average pedigree-based inbreeding coefficient FPED = 0.164; 0.103 to 0.306) and outbred Retinta bulls (LI; average FPED = 0.008; 0 to 0.025). We studied the length of the fragments, their abundance, and genome distribution using high-density microarray data. The number of ROH was significantly higher in the HI group, especially for long fragments (>8Mb). In the LI group, the number of ROH continuously decreased with fragment length. Genome-wide distribution of ROH was highly variable between samples. Some chromosomes presented a larger number of fragments (BTA1, BTA19, BTA29), others had longer fragments (BTA4, BTA12, BTA17), while other ones showed an increased ROH accumulation over specific loci (BTA2, BTA7, BTA23, BTA29). Similar differences were observed in the analysis of 12 individuals produced by a similar inbred event (FPED3 = 0.125). The correlation between the fraction of the genome covered by ROH (FROH) and FPED was high (0.79), suggesting that ROH-based estimations are indicative of inbreeding levels. On the other hand, the correlation between FPED and the microsatellite-based inbreeding coefficient (FMIC) was only moderate (r = 0.44), suggesting that STR-based inbreeding estimations should be avoided. Similarly, we found a very low correlation (r = -0.0132) between recombination rate and ROH abundance across the genome. Finally, we performed functional annotation analyses of genome regions with significantly enriched ROH abundance. Results revealed gene clusters related to pregnancy-associated proteins and immune reaction. The same analysis performed for regions enriched with recently formed ROH (> 8 Mb) showed gene clusters related to flagellum assembly. In both cases, the processes were related to male and female reproductive functions, which may partially explain the reduced fertility associated with inbred populations.
Citation: Goszczynski D, Molina A, Terán E, Morales-Durand H, Ross P, Cheng H, et al. (2018) Runs of homozygosity in a selected cattle population with extremely inbred bulls: Descriptive and functional analyses revealed highly variable patterns. PLoS ONE 13(7): e0200069. https://doi.org/10.1371/journal.pone.0200069
Editor: Arda Yildirim, Gaziosmanpasa University, TURKEY
Received: September 18, 2017; Accepted: June 19, 2018; Published: July 9, 2018
This is an open access article, free of all copyright, and may be freely reproduced, distributed, transmitted, modified, built upon, or otherwise used by anyone for any lawful purpose. The work is made available under the Creative Commons CC0 public domain dedication.
Data Availability: Data were obtained from the National Association of Breeders of Select Retinta Cattle (ANCRE) via a signed data access agreement which does not allow for data sharing. Interested, qualified researchers may request these data by contacting ANCRE at razaretinta@retinta.es. The authors confirm that they did not have any special access privileges to these data.
Funding: National Association of Breeders of Select Retinta Cattle kindly supplied all the pedigree data and DNA samples. This work was funded by the Veterinary Genetics Diagnostic Laboratory (LDGV) and the Molecular and Applied Animal Cytogenetic Laboratory (MERAGEM Research Group) at the University of Cordoba (Spain) and by the Agencia Nacional de Promoción Científica y Tecnológica (ANPCyT, Argentina) (grant PICT2016-4832). Ester Terán is supported by a doctoral internal fellowship of CONICET.
Competing interests: The authors have declared that no competing interests exist.
Introduction
Inbreeding depression is the reduced survival and fertility of offspring from related individuals. According to Leroy [1], the decreased fitness could be mediated by three different genetic mechanisms: 1) an increased phenotypic expression of recessive deleterious mutations (dominance hypothesis); 2) the lack of phenotypic advantages provided by heterozygote genotypes at loci that are maintained by balancing selection at intermediate frequencies (overdominance hypothesis); and 3) an epistatic interaction in which the probability of favourable gene combinations for heterozygotes is higher. However, Charlesworth and Willis [2] proposed the main cause of inbreeding depression to be the dominance hypothesis rather than the other two theories, since large contributions of overdominant genes and epistatic interactions have not been clearly identified yet. Inbreeding depression has been observed in several species and taxon [3, 4]. In cattle, inbreeding depression has been associated with reductions in productive traits, longevity and the ability of the individuals to cope with environmental challenges [5–8]. However, high levels of homozygosity have also been recently reported to be compatible with life and livelihood in cattle isolated for many generations (even up to 340 years) [9].
The increase in the average inbreeding coefficient observed in cattle during the last 20 years may have been caused by two main reasons: 1) the high selection intensity applied in the breeding schemes of highly productive breeds [10], and 2) the use of small and well-adapted populations in extensive and isolated production systems characterized by absence of breeding management [11]. The inbreeding coefficient F [12] is the most commonly used parameter in breeding schemes since it can be easily determined from pedigree data. Through its analysis, it has been suggested that the individual productivity begins to be affected when F reaches values higher than 5% [13]. For instance, Sewalem, Kistemaker [5] demonstrated that culling risk increased to near 14% when F values exceeded 6.25%. Therefore, individuals with F values higher than 20% are rarely found in commercial cattle herds [6], which constitutes the main source of genomic data destined to research in this species [14]. To our knowledge, the only existing study analyzing highly inbred livestock animals has been performed on a limited sample of Chillingham cattle (n = 16; [9]).
The F coefficient has been conventionally estimated from pedigree data (FPED) [15] or microsatellite markers (FMIC) [16]. However, the development of high-throughput genomic technologies based on SNP genotyping has allowed the establishment of novel approaches to determining inbreeding levels [17]. Nowadays, it has been demonstrated that SNP-based estimations of inbreeding are substantially more accurate and often less biased than FPED, even when large pedigrees are considered [18]. Additionally, SNP techniques allow the estimation of inbreeding values in individuals for which pedigree data are not available [19, 20]. One of the most common approaches for studying inbreeding consists in determining runs of homozygosity (ROH), which are long segments of the genome where the alleles are identical because both parents inherited them from a common ancestor (usually referred as autozygous) [21]. ROH analysis has been widely proposed and validated as a tool to estimate individual’s inbreeding level [14] in several species including beef [22–24] and dairy [25, 26] cattle. Moreover, the mean length of the ROH has been associated with the number of generations since the common ancestor [27], allowing to infer the demographic history of a given population.
Microarray density and parameter settings have been proven essential for accurate ROH identification [28]. In this aspect, proper settings (determined per Mb of ROH length) are crucial to avoid biased results due to genotyping errors and missing calls [14]. It has been demonstrated that medium-density (MD) genotype data tend to overestimate the number of fragments shorter than 4 Mb, since many heterozygous SNPs located within the fragments are not genotyped [28]. On the contrary, high-density (HD) microarrays provide a higher number of genotypes at the expense of more genotyping errors and missing calls, but these limitations can be accounted for by fine adjustment of software parameters, indicating that HD microarrays are a more reliable source of consistent results. Nevertheless, most studies in cattle have been performed using MD arrays (50,000 to 70,000 markers per individual), which is the way breeders associations normally genotype elite individuals for use in genomic-assisted breeding programs [25].
In humans, the accumulation of ROH in certain genomic positions has been used to analyze the demographic history of populations [29]. This strategy has also been employed to compare and characterize beef [22] and dairy cattle breeds [30]. Since ROH are normally abundant in regions under positive selection [31], their accumulation at specific loci or “hotspots” has been studied to identify genomic regions that reflect directional selection in cattle [22]. However, the same analysis has been recently employed to detect functional variants associated with inbreeding depression [32] and QTL´s [33]and the genetic control of reproductive traits [34, 35] and diseases [36]. On the other hand, regions with an unusually low abundance of ROH, (“cold spots” [37]) are thought to harbor loci with critical functions escaping lethal or damaging recessive variants. In both cases, it is worth mentioning that the inheritance of such autozygous fragments is also subject to the stochasticity of the recombination events across the genome. In this context, the accumulation of ROH along the genome has become the starting point for other techniques that aim at identifying biological factors acting behind the depressed phenotypes as alternatives to conventional GWAS studies [22, 33–36, 38].
Retinta is the second largest Spanish cattle breed characterized by the quality of its meat, its rusticity and its adaptation to marginal pasturelands and extreme weather [39]. All the individuals of this breed have been raised in the south part of the Spanish peninsula under dry and hot environmental conditions with scarce foraging resources. [40]. The Retinta cattle breeding plan was established more than three decades ago and has been focused on the selection of animals with increased growth speed and extended productive life, while maintaining its adaptation to harsh environments [41]. However, since inbreeding control is not often considered by breeders, this breed is a rare case of individuals with extremely variable FPED coefficients within the population. This makes the Retinta breed an interesting model for animal inbreeding research.
The aims of this study were to characterize the occurrence of ROH in a Retinta population with extremely inbred and outbred animals. We analyzed the length, number, and genomic distribution of the ROH, obtained from HD genotyping data, as well as the relation between long ROH and recent inbreeding events retrieved from pedigree and microsatellite records. We also determined the differences among ROH patterns in a group of individuals with the same increase in FPED values over the last 3 generations to evaluate the reliability of FPED as a predictor of inbreeding depression. Finally, we estimated the influence of the recombination rate on the ROH patterns and identified loci putatively affected by inbreeding in highly related individuals through a functional annotation analysis.
Materials and methods
Animal samples
All the samples were obtained from the National Association of Breeders of Select Retinta Cattle (ACRE, by its Spanish acronym). DNA was obtained from frozen sperm straws (one per bull) using the HigherPurity™ Tissue Genomic DNA Purification Kit (Canvax Biotech, Cordoba, Spain). Management and treatment of the animals during sperm extraction were complied the ethical guidelines of the Animal Ethics Committee of the University of Cordoba, Spain. No ethics committee approval was necessary since no experimental procedures were carried out in animals during this study.
In total, we analyzed 54 bulls belonging to the artificial insemination (AI) program of the National Breeders Association of Select Retinta Cattle (ACRE, for its Spanish acronym). Individuals were selected based on the number of equivalent complete generations (ECG) and FPED values. FPED, FPED3 (FPED for last three generations) and ECG values were determined as described by Meuwissen and Luo [15] and Maignel, Boichard [42] using ENDOG software [43]. FMIC was determined using a 17 STR-based genomic test recommended by the International Society for Animal Genetics for parentage testing as described by Caballero and Toro [44] using Molkin software [45]. For the purposes of our analyses, bulls were classified into two groups: inbred (HI, n = 32; average FPED = 0.164; 0.103 to 0.306) and outbred (LI, n = 22; average FPED = 0.008; 0 to 0.025).
Genotyping
Samples were genotyped using the Axiom® Genome-Wide BOS 1 Bovine high-density SNP Array (Affymetrix, Santa Clara, Ca, USA) in a GeneTitan® Multi-Channel platform (Affymetrix). Samples were run at the IGEVET Genomics Core Facility (University of La Plata, La Plata, Argentina). Raw data were processed using Axiom™ Analysis Suite software (Affymetrix) and setting call rates and DQC at 0.97. Commonly used genotype filters such as minor allele frequency (MAF) were not applied since they would have led to an underestimation of ROH [28]. MAF filtering was avoided given the extreme inbreeding values existent in our population and the high number of SNPs with very low MAF. Discarding such a high number of SNPs would have interfered with some of the main purposes of our study Moreover, skipping such filters is common among ROH-related studies [9, 46]). SNPs that were either orphan or assigned to sex chromosomes or mitochondrial DNA were excluded from the analyses. The number of remaining variants was of 624,737.
Detection and classification of ROH
ROH were estimated using cgaTOH software [47]. The minimum number of SNPs needed to constitute a ROH (L) was calculated according to Purfield, Berry [48], as follows:
where ns is the number of SNPs per individual, ni is the number of individuals, α is the percentage of false positive ROH (0.05), and het is the average SNP heterozygosity across all SNPs. ROH were divided into five different categories according to their length: 1–2 Mb, 2–4 Mb, 4–8 Mb, 8–16 Mb and >16 Mb. Short ROH (< 1 Mb) were not included in the analyses since many of them may derive from the inheritance of common allozygous haplotypes [25].
The number of heterozygous (nH) and missing (nM) genotypes allowed in each ROH category for each chromosome was calculated as
where mL is ROH minimum length, dS is the average distance between SNPs in the chromosome, eG is the genotyping error rate (0.25% according to Affymetrix standard procedures), and mG is the average missing genotype rate in the chromosome. The number of missing and heterozygous SNPs allowed per length category and chromosome are shown in S1 Table. Since each minimum length required a new round of ROH identification, fragments were merged in case of overlaps to avoid the underestimation of long ROH by using the reduce function of the GenomicRanges R package [49].
ROH and ROH-based inbreeding coefficient (FROH) characterization
The mean, minimum and maximum length of the ROH and the standard error were determined for each sample group globally and also for each chromosome, individual and length category. Statistical differences were analyzed using a T standard test with a significance threshold of 0.05. FROH values were estimated for each individual and chromosome as the sum of all ROH divided by the genome (or chromosome) length, according to McQuillan, Leutenegger [21]. Results were expressed as Avg. ± standard error of the mean (S.E.M.)
We evaluated the reliability of FROH as inbreeding estimator compared with FPED, FPED3, and FMIC through Spearman correlations, since the distribution of the SNP data among samples was not normal (Kolmogorov–Smirnov test). The correlation between FROH>8Mb and FPED (six generations), as well as between FROH>16Mb and FPED3 (three generations), were particularly estimated to evaluate the relationship between the length of ROH and the number of generations since the common ancestor, as was proposed by Fisher [50].
To determine if ROH were accumulated over loci with low recombination rates, we estimated the correlation between the number of ROH detected at a given position and the recombination rate at the same locus. To this end, recombination rate at each genomic position was obtained from data reported by Ma, O'Connell [51]. This was performed in the whole population and each group separately, both globally and by chromosome, and considering all ROH and ROH>8Mb.
Gene clustering and functional analysis
To identify putative effects on the biological functions of the highly inbred animals, we performed a functional annotation analysis considering genes within loci that were statistically enriched with ROH. These loci were determined by performing a permutation test under the null hypothesis that ROH are equally distributed along the genome. The number of ROH at each position was randomized one million times using newly designed scripts in JULIA [52]. Then, the p-value for each position was calculated as the relative frequency of randomizations that produced a more extreme number of ROH than the one observed. Significantly enriched genomic intervals were defined as consecutive significant SNPs separated by no more than 1 Mb. These intervals were used to retrieve gene lists from ENSEMBL BioMart v89 [53]. The in silico analysis was conducted using the Functional Annotation Clustering tool implemented in DAVID [54] considering the following annotation categories: Cog_Ontology, Up_Keywords, Up_Seq_Feature, Goterm_Bp_Direct (biological processes), Goterm_Cc_Direct (cellular components), Goterm_Mf_Direct (molecular functions), Kegg_Pathway, Interpro, and Smart. The classification stringency was set to medium and groups were defined by enrichment scores greater than 1.301, which is equal to–log10 (0.05). The analysis was performed on the full set of statistically enriched regions (considering all ROH, regardless of their length) and on a reduced set comprising only ROH>8Mb to evaluate the effect of recently formed fragments.
Results
ROH characterization
The number of ROH was statistically different between the HI and LI groups (106.96 ± 37.56 in HI vs 39.63 ± 23.67 in LI; P<0.001). This statistical difference was also observed at the length category level (Table 1). The HI/LI relation in ROH>8Mb is more than four times higher than in short ROH, demonstrating a clearly higher occurrence of long ROH in the inbred group.
Differences between groups were analyzed using a T-test. Results are expressed as mean ± S.E.M.
A similar analysis performed by chromosome showed similar results, with statistical differences in 25 different chromosomes (BTA5, BTA14, BTA25, and BTA27 were not significant; P>0.05; Fig 1), suggesting that some chromosomes are more propense to accumulate ROH than others.
Results are expressed as average number of ROH per chromosome ± S.E.M. N.S.: Not significant; P>0.05; T-test.
The average ROH length was also significantly different between groups (3.51 ± 0.08 Mb in HI vs 2.23 ± 0.10 Mb in LI; P < 0.001) (Fig 2A), but HI showed a higher variation among individuals (Fig 2B). In this sense, maximum length varied from 4.71 (H24) to 66.78 Mb (H4). Samples H7, H8, H10, H14, H24, and H30 showed a lower variability, lacking ROH > 11.5 Mb even though their FPED3 values were higher than 0.10. The LI group showed a lower variation, although samples L2, L8, L9 and L22 seemed slightly more variable in length of ROH (Fig 2B). Samples L2, L8, L9, L16, L17, and L22 showed ROH>16 Mb even though their FPED was of zero.
Results are expressed as average length ± S.E.M. by individual.
Significant differences were also observed in the length of ROH by chromosome but in a lesser extent. Only seven chromosomes showed statistical differences between groups (BTA1, BTA7, BTA10, BTA13, BTA18, BTA22, BTA24; P < 0.05; Fig 3). The number of chromosomes containing ROH longer than 20 Mb was higher in HI (25 chromosomes) than in LI (6 chromosomes). Surprisingly, a very long ROH (56.73 Mb) was detected on BTA11 in the L9 sample. Details on the distribution and length of ROH by individual and chromosome are available in S2 and S3 Tables, respectively.
Results are expressed as average length ± S.E.M. by individual. *: P < 0.05 T-test by group.
Comparison of FPED, FMIC, and FROH as identity by descent (IBD) estimators
The comparison between IBD estimators showed significant differences for all the inbreeding estimators, but for FMIC, the p-value was at least two orders higher (Table 2). Additionally, FMIC showed the highest value in the LI group as well as the lowest HI/LI ratio, suggesting that STR-based inbreeding coefficients may tend to overestimate the value in outbred individuals.
In the HI group, FROH values ranged from 0.06 to 0.28 (0.15 ± 0.09), showing a high correlation with FPED (0.79). On the contrary, correlations between FPED and FMIC (0.44) and FROH and FMIC (0.36) were much lower, which is consistent with comparisons performed between groups. Correlations between FPED and FROH>8Mb, FPED3 and FROH>16Mb, and FPED3>0 (individuals with a common ancestor in the last three generations) and FROH>16Mb were 0.60, 0.54 and 0.42, respectively.
Analysis of FROH estimates by fragment length and chromosome
FROH values were significantly higher in the HI group (P<0.001). The percentage of FROH value explained by each length category is shown in Fig 4. Short ROH (1–2 Mb) accounted for nearly 50% of the FROH value in LI, showing a decrease in contribution towards the 8–16 Mb category. On the other hand, a similar contribution was observed for each length category in HI. The percentage of FROH explained by short fragments (1–2 Mb) was significantly higher in LI, whereas the percentage explained by long fragments (4–8 Mb, 8–16 Mb, and >16 Mb) was significantly higher in HI.
*: P<0.05. T-tests by group.
The analysis by chromosome showed a higher FROH in HI in 26 chromosomes (Fig 5; P<0.05). Only BTA21, BTA26, and BTA27 were not significant, probably due to the variability observed among individuals, rather than by a decreased difference between means. The variability within groups was high, but some chromosomes, such as BTA14, showed a higher abundance of ROH than others (e.g. BTA20) in both groups. It is also noteworthy that some outbred individuals showed unexpectedly high abundance in some of the chromosomes. For example, sample L9 had 57% of chromosome BTA11 covered by a single ROH (data not shown).
Only chromosomes marked as N.S. were nonsignificant (P > 0.05).
Effect of the recombination rate on ROH abundance
The correlation between the occurrence of a ROH at a given position (defined as the number of individuals with a ROH at that position) and the recombination rate at the same locus was -0.0132. Similar results were obtained when the analysis was performed by group (r = -0.0178 in HI and r = -0.02715 in LI) and chromosome (S4 Table). When all samples from both groups were analyzed at the same time, the highest correlation was found at BTA23 (r = 0.11). BTA28 showed the highest correlation in HI (r = 0.11) and BTA26 the highest one in LI (0.08). No strong correlations were either observed in the analysis of putative recent ROH (r = 0.0025; 0.0003 and 0.028 for the full sample set, HI and LI respectively) nor in the case of individual chromosomes (S4 Table).
Genomic localization of ROH´s
The accumulation of ROH across the genome of the highly inbred samples was analyzed by length category: all ROH lengths (Fig 6A), ROH>8Mb (Fig 6B) and ROH>16Mb (Fig 6C). Although the distribution of the ROH was relatively even and accumulation was moderate in general, we found a few outstanding peaks with a high occurrence of ROH. For example, 26 animals presented ROH at the first million bases (1,575–1,606,000) of BTA23, as well as in BTA2 (71,882,569–73,223,813), BTA7 (51,003,562–53,466,651) and BTA29 (38,490,349–39,617,454) to a lower extent. In the analysis of ROH > 8 Mb (hypothetically generated during the last six generations), the highest occurrence (nine individuals) was detected at BTA24 (39,376,427–58,300,844). Lastly, when we analyzed ROH>16Mb (three generations since the common ancestor), the highest occurrence was found at BTA 7, BTA13, BTA14, BTA18, BTA22 and BTA24. Interestingly, two chromosomes (BTA26 and BTA29) showed a complete absence of ROH> 16Mb.
Number of ROH detected at each SNP position of the array considering ROH of different length categories. Footnote: A: Total FROH; B: 8–16 MB FROH; C: >16 Mb FROH.
To further show the high variability within chromosomes, we illustrated the distribution of the ROH on BTA24, which showed the highest average FROH and average length of ROH in the HI group (Fig 7). It is noteworthy that some individuals showed no ROH across all the chromosome (H1, H9, H10, H15, and H26). On the contrary, some individuals harbored ROH in more than 60% of their chromosomal extension (H13 and H32).
Finally, the relation of FPED with IBD was analyzed by plotting the genome-wide distribution of ROH>16Mb in 12 individuals with identical FPED3 (0.125) (Fig 8). The pattern showed that the distribution of IBD fragments hypothetically acquired during the last 3 generations was highly uneven and variable among individuals. Furthermore, no ROH>16Mb were detected in the H10 sample nor in seven chromosomes (BTA2, BTA17, BTA21, BTA24, BTA26, BTA28 and BTA29). These results suggest that the phenotypic effect of inbreeding in individuals with the same FPED value could be highly divergent.
Candidate regions and genes
To analyze putative biological functions affected in zones with high ROH abundance, we retrieved a list of genes from the intervals significantly enriched based on the permutation test. As before, to characterize recently formed fragments, we conducted separate analysis for all ROH fragments (no length restriction) and for ROH>8Mb. The lists of genes obtained were then submitted to DAVID Bioinformatics for functional annotation clustering. When all fragment lengths were considered, the number of retrieved known genes was of 769. These genes were distributed across 22 different chromosomes, but eight of them were lightly represented (less than 10 genes). The most represented chromosomes in terms of genes covered by ROH were BTA7, BTA14, and BTA12 with 108, 94 and 77 genes, respectively. Three significant functional clusters were detected in the annotation clustering analysis (Table 3). The one with the highest enrichment score was related to protein catabolic processes and involved the PAG family of genes. These genes are located on BTA29 and code for precursors of the pregnancy-associated glycoprotein family (PAG1, PAG4, PAG7, PAG15, PAG16, PAG19, PAG20, and PAG21), a group of aspartic proteases secreted by the placenta. The second most enriched cluster was related to 5 β-defensin genes, which codify for antimicrobial peptides found in white blood cells.
When focused on regions ROH>8Mb, that hypothetically formed during the last six generations, we retrieved 362 known genes. These genes were distributed over seven chromosomes. BTA22, BTA24, and BTA14 were the most represented ones, with 87, 68 and 64 genes, respectively. Three significant clusters were identified in the functional annotation analysis (Table 4). The most significant one was related to the assembly of the flagellum and cilia microtubules involving 4 genes (DNAAF1, DNAH1, LRRC6, and ZMYND10). The remaining two clusters were related to calcium binding and metabolism, and hyaluronan metabolism.
Discussion
In this study, we characterized the abundance and distribution of ROH in inbred and outbred cattle. To our knowledge, this is the first genomic report to use HD genotyping to analyze the relation between IBD fragments and known recent inbreeding events in individual cattle with extreme F values. This is also the first analysis of a cattle population that has been bred and selected under standardized conditions and can be clustered into two highly divergent groups in terms of individual inbreeding values.
FROH as inbreeding estimator on highly inbred cattle
Our results showed a higher correlation between FPED and FROH in comparison with previous reports [25], but only when the whole set of samples was considered. Correlations were lower when both groups were analyzed separately, showing values more consistent with those reported previously [14]. Three causes were suggested as source of divergence between FROH and FPED: 1) the persistence of ancestral short ROH through time due to low recombination rates, which are ignored in the estimation of FPED [55], 2) the depth and reliability of pedigree information [18, 56] and 3) the stochastic nature of IBD inheritance [57]. In this study, pedigree errors were minimized using molecular parentage tests for all individuals. Similarly, pedigree data were shown reliable since outbred animals (LI) showed low FROH values, which may indicate that inbreeding events that occurred in the population before the fifth parental generation are scarce. For this reason, our correlations are hypothetically reliable. On the contrary, the STR-based inbreeding estimations (FMIC) were poorly correlated with both FPED and FROH, demonstrating that the use of STR is a poor source of information in highly inbred individuals. This fact was previously described in cattle but using a moderately inbred population, as was shown by Baumung and Solkner [58].
Seven years ago, Howrigan, Simonson [27] modeled the relationship between the length of ROH and the number of generations since the common ancestor by large-scale simulations. The study was based on the concept proposed by Fisher [50] in which, the length of the IBD fragments is associated with the number of generations since the common ancestor. For instance, they showed that ROH>16Mb are likely inherited from parents that shared a common ancestor three generations before. More recently, Marras, Gaspa [24] and Ferencakovic, Hamzic [55] analyzed this concept in commercial populations of several cattle breeds but using individuals with low FPED. Since our population had a well-known history of recent inbreeding, we expected a higher correlation between estimates coming from long ROH (FROH derived from ROH>8 and ROH>16) and FPED values from the corresponding number of generations (FPED for six and FPED3 for three). There was a general agreement between increased ROH lengths and inbreeding events in the inbred individuals, showing correlations between r = 0.6 and r = 0.42 when the long ROH fragments were used to estimate FROH. However, these correlations were lower than the overall correlation obtained using the analysis of the whole population. One possible explanation is the stochastic nature of IBD inheritance [57]. However, Kardos, Luikart [18] demonstrated that FPED can be easily underestimated when the pedigree depth is shorter than 20 generations. In our case, FPED was estimated using only 6 generations on average, thus the influence of previous common ancestors was not accounted for.
ROH length distribution
The length of the ROH was highly variable between and within the two sample groups. For instance, five LI samples showed ROH longer than 20 Mb, when presence of long ROH was not expected. Although the mean length of the ROH should at some extent reflect the number of generations since the common ancestor, we should also consider that ROH formation is a complex process involving dynamic rates of recombination along the genome and the stochastic nature of gamete formation processes [38]. In our study, this hypothesis is also supported by the variation detected in the HI group, where some animals had ROH as long as 50 Mb (produced hypothetically by an inbreeding event that occurred one generation ago), which was inconsistent with the FPED. A similar variability was also described in Brown Swiss cattle and Valle de Belice sheep [28, 38]. However, in both studies pedigree data was not available, thus preventing the comparison between FPED and ROH length.
Significant differences were also observed in the number of ROH detected by length category. Shorter ROH predominated in LI, but the number of ROH decreased as the length increased. This pattern may be attributed to a “foundational inbreeding” produced during the creation of the breed, as suggested in other breeds and species [22, 24, 59]. On the other hand, HI showed a number of ROH [8-16Mb] and ROH>16 Mb that was eight and ten times higher than in LI. To our knowledge, high differences in the number of ROH were not reported previously in livestock. In our case, the experimental design, including two groups of animals which diverged only in terms of inbreeding (same breed, selection scheme, and breeding objectives) could explain the differences observed between groups. In any case, our results were consistent with expectations, since inbreeding events that occurred more recently (6 generations) tended to produce longer ROH fragments. The fact that only 1.5% percent of the genome was covered by short ROH in LI may indicate that the original population from which the Retinta breed was founded 50 generations ago was ample and scarcely selected.
Chromosomal distribution of ROH
In a recent study, Zavarez, Utsunomiya [59] reported a reduced variation in the distribution of the ROH among the chromosomes of Nelore cows, with FROH values ranging from 0.05 to 0.1. On the contrary, Ferencakovic, Solkner [28] and Szmatoła, Gurgul [22] showed that the length of the ROH and their localization in the genome could be extremely variable in cattle. In our case, the degree of variation among chromosomes was high, which agrees with the latter. This variation within and between chromosomes has been also described in dairy [56] and beef [22] cattle populations, but with low inbreeding levels and with a lesser extent. These findings may suggest that some genomic loci are less able to sustain accumulation of IBD than others. In our case, it may be explained by the existence of deleterious alleles located in ROH located in specific regions (For instance in regions with systematically lower genomic constraint [60]). Those alleles will decrease the ROH abundance in those locations since if a variant in a population is lethal in homozygous form, inbreeding will greatly increase the chance of generating a lethal genotype. In this experiment, it is reasonable to think that the hotspots detected might have been produced mostly because of mating close relatives in absence of additional conditions that can differentiate the groups analyzed. But also, the fact that recombination rates, which have been also pointed as the cause of ROH increase in specific regions, seemed to exert no effects on ROH abundance, is supporting the previous idea. Finally, variation was also high when we analyzed the genome distribution of ROH>8Mb and ROH>16Mb. Interestingly, some genomic regions were found to lack long ROH (e.g. BTA26 and BTA29 for ROH>16Mb). These findings may suggest that some genomic areas cannot be affected by inbreeding within a short time lapse.
The uneven genome-wide distribution of the ROH is particularly important since breeders usually associate an increase in the inbreeding coefficient with a proportional detrimental effect [61, 62]. However, the phenotypic effects of the IBD blocks will be determined by their localization on the genome. Forty years ago, Franklin [63] modeled and determined variation expected in terms of homozygous blocks in individuals resultant from a similar inbred event in D. melanogaster. In that study, the author demonstrated that the chance of finding a homozygous genotype at a specific position of the genome after an event of inbreeding depends on the recombination rate of the locus and the length of the chromosome. According to the author, the abundance of inbreeding blocks should be higher in chromosomes longer than 1M, but the variability should be lower in chromosomes shorter than 0.5M. To show that, we analyzed the distribution of recently formed ROH (>8Mb and >16Mb) on 12 individuals with the same increase in F during the last three generations (FPED3 = 0.125). Once again, results showed an uneven genome distribution with a moderate accumulation of ROH over specific loci (up to nine ROH>8Mb on some loci of BTA24). Besides, seven chromosomes showed a total absence of ROH >16Mb. Despite the limited number of samples, our study shows that high inbreeding levels might not always lead to inbreeding depression, as was recently also shown in a study of isolated cattle [9]. In fact, even the mating of unrelated animals might possibly lead to inbreeding depression if, by chance, they happened to present long ROH at some biologically important loci. We showed that although accumulation of fragments over specific loci is not common, such long ROH are easy to find among outbred individuals. Furthermore, we showed that the genome-wide distribution of the ROH was generally uneven, which means that the phenotypic effects of inbreeding depression may vary between animals regardless of their FPED.
Functional annotation clustering
The functional study of genomic regions significantly enriched with ROH constitutes a practical approach to identify metabolic pathways putatively affected by inbreeding. In our case, we performed two separate analyses: one considering all the ROH fragments and another one considering only putatively modern ROH. The purpose of this separation was to distinguish the putative effects of short-term inbreeding from the rest. But also, it was suggested that the analysis of biological functions affected by recent inbreeding (ROH >8 Mb) is a more powerful tool than the use of ROH of shorter lengths to analyze inbreeding depression [18].
The analysis including all the ROH showed a highly enriched annotation cluster including several precursors of the PAG family. These proteins are used as pregnancy-status indicators since their expression varies through pregnancy stages, and they are associated with placental mass, fetus number, and birth weight in cattle [64, 65]. Nowadays, the detection of PAG proteins is being used as an early pregnancy test in beef [66] and dairy cattle [67, 68], and also as an early marker of twin pregnancies [69]. Recently, the same protein family has been associated with milk yield, clinical diseases [70], pregnancy loss [71] and retained membranes [72]. In our study, PAG catabolism is the most significant pathway affected in inbreed individuals, suggesting that the reduced fertility observed in inbred bulls of this breed [73] could be partially mediated by the metabolism of those proteins.
In cattle, inbreeding has been also associated with changes in the immune system [74] and an increased occurrence of respiratory diseases [75]. Our results also identified a metabolic pathway related to B-defensin genes. Those proteins have been associated with embryos that showed a delayed developmental stage at day 7 after insemination [76] and with cows suffering Staphylococci mastitis [77]. Retarded embryos are often low quality and less viable embryos [78, 79], as well as mastitis has been recently associated with an impaired developmental competence in oocytes [80], even in sub-clinical presentations [81]. Both cases could also partially explain the mechanisms involved in the reduction of fertility observed in inbred individuals [73].
The analysis of the genomic positions affected by recently formed ROH showed three enriched clusters, but with a lower significance compared with the previous analysis. The most affected cluster was associated with microtubule structures, flagellum assembly, inositol polyphosphate kinase activity, chemotaxis and ATP-binding features. This cluster included the DNAAF1, DNAH1, LRRC6 and ZMYND10 genes. All those processes have been related to an impaired axonal assembly of the dynein arms that produced an abnormal movement of flagellum and cilia [82, 83]. Recently, mutations in DNAH1 have been heavily associated with dysplasia of the sperm fibrous sheet [84], several flagellar defects and asthenozoospermia [85] in humans. These findings agree with a previous study which demonstrated that highly inbred Retinta bulls present increased hyperactivation-like motility and a reduced reproductive performance in field conditions [73]. We hypothesize that the accumulation of deleterious homozygote variants through excessive inbreeding events may have triggered this biological condition.
Conclusion
In this study, we analyzed the use of FROH as IBD predictor in extremely inbred Retinta bulls. We showed that although it is correlated to the pedigree-based inbreeding coefficient (FPED), ROH seemed to distribute in an uneven fashion among individuals with similar FROH and FPED coefficients. Furthermore, we found that individuals with a recent history of high inbreeding showed an increased number of long ROH (>8Mb). In this population, specific regions of the genome showed high accumulation of ROH, which was not associated to the recombination rate. ROH harbored genes related to pregnancy-associated proteins, cell motility, skeletal reorganization and immune system. These results are in line with a previous study where we showed that some of these bulls presented an increased hyperactivation-like motility pattern. Additionally, we demonstrated that animals with the same inbreeding coefficient may present different phenotypes, which may not always lead to a detriment in the production trait of interest.
Supporting information
S1 Table. ROH identification parameters at each length size per chromosome.
Footnote: L: minimum number of SNPs; nH: number of heterozygous genotypes allowed per ROH; nM: number of missing genotypes allowed per ROH.
https://doi.org/10.1371/journal.pone.0200069.s001
(DOCX)
S2 Table. Descriptive statistics of ROH occurrence per chromosome and group of individuals.
Footnote: Statistical differences were estimated between groups per chromosome using a T test.
https://doi.org/10.1371/journal.pone.0200069.s002
(DOCX)
S3 Table. Descriptive statistics of ROH length per chromosome and group of individuals.
Footnote: Statistical differences were estimated between groups per chromosome using a T test.
https://doi.org/10.1371/journal.pone.0200069.s003
(DOCX)
S4 Table. Correlation between recombination rate and ROH abundance per chromosome.
Footnote: Data was analyzed by chromosome in the whole population (WP) as well in HI and LI separately, using ROH abundance of any length or ROH>8Mb. Spearman correlation. Correlations were estimated per SNP within chromosome.
https://doi.org/10.1371/journal.pone.0200069.s004
(DOCX)
Acknowledgments
We thank the National Association of Breeders of Select Retinta Cattle for kindly supplying the pedigree data and sperm samples. We also thank Adriana Di Maggio for language and manuscript edition.
References
- 1. Leroy G. Inbreeding depression in livestock species: review and meta-analysis. Anim Genet. 2014;45(5):618–28. Epub 2014/07/01. pmid:24975026.
- 2. Charlesworth D, Willis JH. The genetics of inbreeding depression. Nat Rev Genet. 2009;10(11):783–96. Epub 2009/10/17. pmid:19834483.
- 3. Charlesworth D, Charlesworth B. Inbreeding Depression and Its Evolutionary Consequences. Annu Rev Ecol Syst. 1987;18:237–68. PubMed PMID: WOS:A1987K958800011.
- 4. Keller LF, Waller DM. Inbreeding effects in wild populations. Trends Ecol Evol. 2002;17(5):230–41. PubMed PMID: WOS:000175024300013.
- 5. Sewalem A, Kistemaker GJ, Miglior F, Van Doormaal BJ. Analysis of inbreeding and its relationship with functional longevity in Canadian dairy cattle. J Dairy Sci. 2006;89(6):2210–6. Epub 2006/05/17. pmid:16702287.
- 6. Gonzalez-Recio O, Lopez de Maturana E, Gutierrez JP. Inbreeding depression on female fertility and calving ease in Spanish dairy cattle. J Dairy Sci. 2007;90(12):5744–52. Epub 2007/11/21. pmid:18024768.
- 7. Ercanbrack SK, Knight AD. Effects of inbreeding on reproduction and wool production of Rambouillet, Targhee, and Columbia ewes. J Anim Sci. 1991;69(12):4734–44. Epub 1991/12/01. pmid:1808170.
- 8. Saura M, Fernandez A, Varona L, Fernandez AI, de Cara MA, Barragan C, et al. Detecting inbreeding depression for reproductive traits in Iberian pigs using genome-wide data. Genet Sel Evol. 2015;47(1):1. Epub 2015/01/18. pmid:25595431; PubMed Central PMCID: PMCPMC4297446.
- 9. Williams JL, Hall SJ, Del Corvo M, Ballingall KT, Colli L, Ajmone Marsan P, et al. Inbreeding and purging at the genomic Level: the Chillingham cattle reveal extensive, non-random SNP heterozygosity. Anim Genet. 2016;47(1):19–27. Epub 2015/11/13. pmid:26559490.
- 10. Stachowicz K, Sargolzaei M, Miglior F, Schenkel FS. Rates of inbreeding and genetic diversity in Canadian Holstein and Jersey cattle. J Dairy Sci. 2011;94(10):5160–75. Epub 2011/09/29. pmid:21943766.
- 11. Cecchi F, Ciampolini R, Castellana E, Ciani E. Genetic diversity within and among endangered local cattle breeds from Tuscany (Italy). Large Animal Review. 2012;18(2):79–85. PubMed PMID: WOS:000303323600006.
- 12. Wright S. Coefficients of inbreeding and relationship. Am Nat. 1922;56(645):330–8. PubMed PMID: WOS:000200856800029.
- 13. Hansen LB. Consequences of selection for milk yield from a geneticist's viewpoint. J Dairy Sci. 2000;83(5):1145–50. Epub 2000/05/23. pmid:10821591.
- 14. Peripolli E, Munari DP, Silva M, Lima ALF, Irgang R, Baldi F. Runs of homozygosity: current knowledge and applications in livestock. Anim Genet. 2017;48(3):255–71. Epub 2016/12/03. pmid:27910110.
- 15. Meuwissen THE, Luo Z. Computing Inbreeding Coefficients in Large Populations. Genet Sel Evol. 1992;24(4):305–13. PubMed PMID: WOS:A1992JW09400002.
- 16. Fernandez J, Toro MA. A new method to estimate relatedness from molecular markers. Mol Ecol. 2006;15(6):1657–67. Epub 2006/04/25. pmid:16629818.
- 17. Kristensen TN, Pedersen KS, Vermeulen CJ, Loeschcke V. Research on inbreeding in the 'omic' era. Trends Ecol Evol. 2010;25(1):44–52. Epub 2009/09/08. pmid:19733933.
- 18. Kardos M, Luikart G, Allendorf FW. Measuring individual inbreeding in the age of genomics: marker-based measures are better than pedigrees. Heredity (Edinb). 2015;115(1):63–72. Epub 2015/06/11. pmid:26059970; PubMed Central PMCID: PMCPMC4815495.
- 19. de Roos AP, Schrooten C, Veerkamp RF, van Arendonk JA. Effects of genomic selection on genetic improvement, inbreeding, and merit of young versus proven bulls. J Dairy Sci. 2011;94(3):1559–67. Epub 2011/02/23. pmid:21338821.
- 20. Villa-Angulo R, Matukumalli LK, Gill CA, Choi J, Van Tassell CP, Grefenstette JJ. High-resolution haplotype block structure in the cattle genome. BMC Genet. 2009;10(11):19. Epub 2009/04/28. pmid:19393054; PubMed Central PMCID: PMCPMC2684545.
- 21. McQuillan R, Leutenegger AL, Abdel-Rahman R, Franklin CS, Pericic M, Barac-Lauc L, et al. Runs of homozygosity in European populations. Am J Hum Genet. 2008;83(3):359–72. Epub 2008/09/02. pmid:18760389; PubMed Central PMCID: PMCPMC2556426.
- 22. Szmatoła T, Gurgul A, Ropka-Molik K, Jasielczuk I, Ząbek T, Bugno-Poniewierska M. Characteristics of runs of homozygosity in selected cattle breeds maintained in Poland. Livest Sci. 2016;188:72–80.
- 23. Ferencakovic M, Hamzic E, Gredler B, Curik I, Sölkner J. Runs of homozygosity reveal genomewide autozygosity in the Austrian Fleckvieh cattle. Agric Conspec Sci. 2011;76(4):325–8.
- 24. Marras G, Gaspa G, Sorbolini S, Dimauro C, Ajmone-Marsan P, Valentini A, et al. Analysis of runs of homozygosity and their relationship with inbreeding in five cattle breeds farmed in Italy. Anim Genet. 2015;46(2):110–21. Epub 2014/12/23. pmid:25530322.
- 25. Kim ES, Sonstegard TS, Van Tassell CP, Wiggans G, Rothschild MF. The Relationship between Runs of Homozygosity and Inbreeding in Jersey Cattle under Selection. PLoS One. 2015;10(7):e0129967. Epub 2015/07/15. pmid:26154171; PubMed Central PMCID: PMCPMC4496098.
- 26. Gurgul A, Szmatola T, Topolski P, Jasielczuk I, Zukowski K, Bugno-Poniewierska M. The use of runs of homozygosity for estimation of recent inbreeding in Holstein cattle. J Appl Genet. 2016;57(4):527–30. Epub 2016/01/25. pmid:26803654.
- 27. Howrigan DP, Simonson MA, Keller MC. Detecting autozygosity through runs of homozygosity: a comparison of three autozygosity detection algorithms. BMC Genomics. 2011;12(1):460. Epub 2011/09/29. pmid:21943305; PubMed Central PMCID: PMCPMC3188534.
- 28. Ferencakovic M, Solkner J, Curik I. Estimating autozygosity from high-throughput information: effects of SNP density and genotyping errors. Genet Sel Evol. 2013;45(1):42. Epub 2013/10/31. pmid:24168655; PubMed Central PMCID: PMCPMC4176748.
- 29. Kirin M, McQuillan R, Franklin CS, Campbell H, McKeigue PM, Wilson JF. Genomic runs of homozygosity record population history and consanguinity. PLoS One. 2010;5(11):e13996. Epub 2010/11/19. pmid:21085596; PubMed Central PMCID: PMCPMC2981575.
- 30. Mastrangelo S, Tolone M, Di Gerlando R, Fontanesi L, Sardina MT, Portolano B. Genomic inbreeding estimation in small populations: evaluation of runs of homozygosity in three local dairy cattle breeds. Animal. 2016;10(5):746–54. Epub 2016/04/15. pmid:27076405.
- 31. Kardos M, Qvarnstrom A, Ellegren H. Inferring Individual Inbreeding and Demographic History from Segments of Identity by Descent in Ficedula Flycatcher Genome Sequences. Genetics. 2017;205(3):1319–34. Epub 2017/01/20. pmid:28100590; PubMed Central PMCID: PMCPMC5340341.
- 32. Zhang Q, Guldbrandtsen B, Bosse M, Lund MS, Sahana G. Runs of homozygosity and distribution of functional variants in the cattle genome. BMC Genomics. 2015;16(1):542. Epub 2015/07/23. pmid:26198692; PubMed Central PMCID: PMCPMC4508970.
- 33. Biscarini F, Biffani S, Stella A. Más allá del GWAS: alternativas para localizar QTLs (In Spanish). ArXiv. 2015. doi: arXiv:1504.03802v1 [q-bio.GN].
- 34. Ferencakovic M, Solkner J, Kaps M, Curik I. Genome-wide mapping and estimation of inbreeding depression of semen quality traits in a cattle population. J Dairy Sci. 2017;100(6):4721–30. Epub 2017/04/25. pmid:28434751.
- 35. Metzger J, Karwath M, Tonda R, Beltran S, Agueda L, Gut M, et al. Runs of homozygosity reveal signatures of positive selection for reproduction traits in breed and non-breed horses. BMC Genomics. 2015;16(1):764. Epub 2015/10/11. pmid:26452642; PubMed Central PMCID: PMCPMC4600213.
- 36. Biscarini F, Biffani S, Morandi N, Nicolazzi EL, Stella A. Using runs of homozygosity to detect genomic regions associated with susceptibility to infectious and metabolic diseases in dairy cows under intensive farming conditions. arXiv. 2016. doi: arXiv:1601.07062v1.
- 37. Pemberton TJ, Absher D, Feldman MW, Myers RM, Rosenberg NA, Li JZ. Genomic patterns of homozygosity in worldwide human populations. Am J Hum Genet. 2012;91(2):275–92. Epub 2012/08/14. pmid:22883143; PubMed Central PMCID: PMCPMC3415543.
- 38. Mastrangelo S, Tolone M, Sardina MT, Sottile G, Sutera AM, Di Gerlando R, et al. Genome-wide scan for runs of homozygosity identifies potential candidate genes associated with local adaptation in Valle del Belice sheep. Genet Sel Evol. 2017;49(1):84. Epub 2017/11/16. pmid:29137622; PubMed Central PMCID: PMCPMC5684758.
- 39. Morales R, Menendez-Buxadera A, Aviles C, Molina A. Direct and maternal genetic effects for preweaning growth in Retinta cattle estimated by a longitudinal approach throughout the calving trajectory of the cow. J Anim Breed Genet. 2013;130(6):425–34. Epub 2013/11/19. pmid:24236605.
- 40.
Jiménez JM, Criado M, Molina A. Las razas bovinas andaluzas de fomento: Retinto, Lidia. (In Spanish). Patrimonio Ganadero Andaluz—Book II. II. Sevilla, España: Junta de Andalucía; 2007. p. 9–52.
- 41. Morales R, Phocas F, Sole M, Demyda-Peyras S, Menendez-Buxadera A, Molina A. Breeding beef cattle for an extended productive life: Evaluation of selection criteria in the Retinta breed. Livest Sci. 2017;204:115–21. PubMed PMID: WOS:000413127500018.
- 42. Maignel L, Boichard D, Verrier E. Genetic variability of French dairy breeds estimated from pedigree information. Interbull Bulletin. 1996;14:49–54.
- 43. Gutierrez JP, Goyache F. A note on ENDOG: a computer program for analysing pedigree information. J Anim Breed Genet. 2005;122(3):172–6. Epub 2005/09/01. pmid:16130468.
- 44. Caballero A, Toro MA. Analysis of genetic diversity for the management of conserved subdivided populations. Conserv Genet. 2002;3(3):289–99. PubMed PMID: WOS:000177718900007.
- 45. Gutierrez JP, Royo LJ, Alvarez I, Goyache F. MolKin v2.0: a computer program for genetic analysis of populations using molecular coancestry information. J Hered. 2005;96(6):718–21. Epub 2005/10/28. pmid:16251515.
- 46. Kukuckova V, Moravcikova N, Ferencakovic M, Simcic M, Meszaros G, Solkner J, et al. Genomic characterization of Pinzgau cattle: genetic conservation and breeding perspectives. Conserv Genet. 2017;18(4):893–910. PubMed PMID: WOS:000405212200013.
- 47. Zhang L, Orloff MS, Reber S, Li S, Zhao Y, Eng C. cgaTOH: extended approach for identifying tracts of homozygosity. PLoS One. 2013;8(3):e57772. Epub 2013/03/08. pmid:23469237; PubMed Central PMCID: PMCPMC3585782.
- 48. Purfield DC, Berry DP, McParland S, Bradley DG. Runs of homozygosity and population history in cattle. BMC Genet. 2012;13:70. Epub 2012/08/15. pmid:22888858; PubMed Central PMCID: PMCPMC3502433.
- 49. Lawrence M, Huber W, Pages H, Aboyoun P, Carlson M, Gentleman R, et al. Software for computing and annotating genomic ranges. PLoS Comput Biol. 2013;9(8):e1003118. Epub 2013/08/21. pmid:23950696; PubMed Central PMCID: PMCPMC3738458.
- 50. Fisher RA. A Fuller Theory of Junctions in Inbreeding. Heredity. 1954;8(2):187–97. PubMed PMID: WOS:A1954YC17800002.
- 51. Ma L, O'Connell JR, VanRaden PM, Shen B, Padhi A, Sun C, et al. Cattle Sex-Specific Recombination and Genetic Control from a Large Pedigree Analysis. PLoS Genet. 2015;11(11):e1005387. Epub 2015/11/06. pmid:26540184; PubMed Central PMCID: PMCPMC4634960.
- 52. Bezanson J, Edelman A, Karpinski S, Shah VB. Julia: A Fresh Approach to Numerical Computing. SIAM Rev. 2017;59(1):65–98. PubMed PMID: WOS:000396692500004.
- 53. Aken BL, Ayling S, Barrell D, Clarke L, Curwen V, Fairley S, et al. The Ensembl gene annotation system. Database: the journal of biological databases and curation. 2016;2016. Epub 2016/06/25. pmid:27337980; PubMed Central PMCID: PMCPMC4919035.
- 54. Huang da W, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nature protocols. 2009;4(1):44–57. Epub 2009/01/10. pmid:19131956.
- 55. Ferencakovic M, Hamzic E, Gredler B, Solberg TR, Klemetsdal G, Curik I, et al. Estimates of autozygosity derived from runs of homozygosity: empirical evidence from selected cattle populations. J Anim Breed Genet. 2013;130(4):286–93. Epub 2013/07/17. pmid:23855630.
- 56. Zhang Q, Calus MP, Guldbrandtsen B, Lund MS, Sahana G. Estimation of inbreeding using pedigree, 50k SNP chip genotypes and full sequence data in three cattle breeds. BMC Genet. 2015;16(1):88. Epub 2015/07/22. pmid:26195126; PubMed Central PMCID: PMCPMC4509611.
- 57. Curik I, Ferencakovic M, Solkner J. Genomic dissection of inbreeding depression: a gate to new opportunities. Revista Brasileira De Zootecnia-Brazilian Journal of Animal Science. 2017;46(9):773–82. PubMed PMID: WOS:000415388200010.
- 58. Baumung R, Solkner J. Pedigree and marker information requirements to monitor genetic variability. Genet Sel Evol. 2003;35(4):369–83. Epub 2003/08/21. pmid:12927072; PubMed Central PMCID: PMCPMC2697993.
- 59. Zavarez LB, Utsunomiya YT, Carmo AS, Neves HH, Carvalheiro R, Ferencakovic M, et al. Assessment of autozygosity in Nellore cows (Bos indicus) through high-density SNP genotypes. Front Genet. 2015;6(JAN):5. Epub 2015/02/18. pmid:25688258; PubMed Central PMCID: PMCPMC4310349.
- 60. Szpiech ZA, Xu J, Pemberton TJ, Peng W, Zollner S, Rosenberg NA, et al. Long runs of homozygosity are enriched for deleterious variation. Am J Hum Genet. 2013;93(1):90–102. Epub 2013/06/12. pmid:23746547; PubMed Central PMCID: PMCPMC3710769.
- 61. Bjelland DW, Weigel KA, Vukasinovic N, Nkrumah JD. Evaluation of inbreeding depression in Holstein cattle using whole-genome SNP markers and alternative measures of genomic inbreeding. J Dairy Sci. 2013;96(7):4697–706. PubMed PMID: WOS:000320745200064. pmid:23684028
- 62. Bjelland DW, Weigel KA, Coburn AD, Wilson RD. Using a family-based structure to detect the effects of genomic inbreeding on embryo viability in Holstein cattle. J Dairy Sci. 2015;98(7):4934–44. Epub 2015/05/11. pmid:25958282.
- 63. Franklin IR. The distribution of the proportion of the genome which is homozygous by descent in inbred individuals. Theor Popul Biol. 1977;11(1):60–80. pmid:404725
- 64. Zoli AP, Beckers JF, Wouters-Ballman P, Closset J, Falmagne P, Ectors F. Purification and characterization of a bovine pregnancy-associated glycoprotein. Biol Reprod. 1991;45(1):1–10. Epub 1991/07/01. pmid:1908709.
- 65. Butler JE, Hamilton WC, Sasser RG, Ruder CA, Hass GM, Williams RJ. Detection and partial characterization of two bovine pregnancy-specific proteins. Biol Reprod. 1982;26(5):925–33. Epub 1982/06/01. pmid:6807365.
- 66. Roberts JN, Byrem TM, Grooms DL. Application of an ELISA Milk Pregnancy Test in Beef Cows. Reprod Domest Anim. 2015;50(4):651–8. Epub 2015/06/11. pmid:26058919.
- 67. Mayo LM, Moore SG, Poock SE, Silvia WJ, Lucy MC. Technical note: Validation of a chemical pregnancy test in dairy cows that uses whole blood, shortened incubation times, and visual readout. J Dairy Sci. 2016;99(9):7634–41. Epub 2016/07/11. pmid:27394956; PubMed Central PMCID: PMCPMC5772973.
- 68. Commun L, Velek K, Barbry JB, Pun S, Rice A, Mestek A, et al. Detection of pregnancy-associated glycoproteins in milk and blood as a test for early pregnancy in dairy cows. J Vet Diagn Invest. 2016;28(3):207–13. Epub 2016/03/06. pmid:26944877.
- 69. Garcia-Ispierto I, Rossello-Visa MA, Serrano-Perez B, Mur-Novales R, de Sousa NM, Beckers JF, et al. Plasma concentrations of pregnancy-associated glycoproteins I and II and progesterone on day 28 post-AI as markers of twin pregnancy in dairy cattle. Livest Sci. 2016;192:44–7. PubMed PMID: WOS:000385605600007.
- 70. Mercadante PM, Ribeiro ES, Risco C, Ealy AD. Associations between pregnancy-associated glycoproteins and pregnancy outcomes, milk yield, parity, and clinical diseases in high-producing dairy cows. J Dairy Sci. 2016;99(4):3031–40. Epub 2016/02/08. pmid:26851856.
- 71. Gatea AO, Smith MF, Pohler KG, Egen T, Pereira MHC, Vasconselos JLM, et al. The ability to predict pregnancy loss in cattle with ELISAs that detect pregnancy associated glycoproteins is antibody dependent. Theriogenology. 2018;108:269–76. Epub 2017/12/25. pmid:29275034.
- 72. Hooshmandabbasi R, Zerbe H, Bauersachs S, de Sousa NM, Boos A, Klisch K. Pregnancy-associated glycoproteins in cows with retained fetal membranes. Theriogenology. 2018;105:158–63. Epub 2017/10/06. pmid:28982025.
- 73. Dorado J, Cid RM, Molina A, Hidalgo M, Ariza J, Moreno-Millan M, et al. Effect of inbreeding depression on bull sperm quality and field fertility. Reprod Fertil Dev. 2017;29(4):712–20. Epub 2015/12/19. pmid:26678044.
- 74. Macedo AA, Bittar JF, Bassi PB, Ronda JB, Bittar ER, Panetto JC, et al. Influence of endogamy and mitochondrial DNA on immunological parameters in cattle. BMC Vet Res. 2014;10:79. Epub 2014/04/04. pmid:24694022; PubMed Central PMCID: PMCPMC4230805.
- 75. Snowder GD, Van Vleck LD, Cundiff LV, Bennett GL. Influence of breed, heterozygosity, and disease incidence on estimates of variance components of respiratory disease in preweaned beef calves. J Anim Sci. 2005;83(6):1247–61. Epub 2005/05/14. pmid:15890802.
- 76. Beltman ME, Forde N, Lonergan P, Crowe MA. Altered endometrial immune gene expression in beef heifers with retarded embryos. Reprod Fertil Dev. 2013;25(6):966–70. Epub 2012/10/05. pmid:23034328.
- 77. Kosciuczuk EM, Lisowski P, Jarczak J, Krzyzewski J, Zwierzchowski L, Bagnicka E. Expression patterns of beta-defensin and cathelicidin genes in parenchyma of bovine mammary gland infected with coagulase-positive or coagulase-negative Staphylococci. BMC Vet Res. 2014;10(1):246. Epub 2014/10/08. pmid:25286984; PubMed Central PMCID: PMCPMC4194403.
- 78. Beltman ME, Forde N, Furney P, Carter F, Roche JF, Lonergan P, et al. Characterisation of endometrial gene expression and metabolic parameters in beef heifers yielding viable or non-viable embryos on Day 7 after insemination. Reprod Fertil Dev. 2010;22(6):987–99. Epub 2010/07/02. pmid:20591333.
- 79. Hockett ME, Rohrbach NR, Schrick FN. Alterations in embryo development in progestogen-supplemented cows administered prostaglandin F2α. Prostaglandins Other Lipid Mediat. 2004;73(3–4):227–36. pmid:15287154
- 80. Roth Z, Dvir A, Kalo D, Lavon Y, Krifucks O, Wolfenson D, et al. Naturally occurring mastitis disrupts developmental competence of bovine oocytes. J Dairy Sci. 2013;96(10):6499–505. Epub 2013/08/21. pmid:23957998.
- 81. Roth Z, Asaf S, Furman O, Lavon Y, Kalo D, Wolfenson D, et al. Subclinical mastitis disrupts oocyte cytoplasmic maturation in association with reduced developmental competence and impaired gene expression in preimplantation bovine embryos. Reprod Fertil Dev. 2015;28(11):1653–62. Epub 2015/04/22. pmid:25891636.
- 82. Moore DJ, Onoufriadis A, Shoemark A, Simpson MA, zur Lage PI, de Castro SC, et al. Mutations in ZMYND10, a gene essential for proper axonemal assembly of inner and outer dynein arms in humans and flies, cause primary ciliary dyskinesia. Am J Hum Genet. 2013;93(2):346–56. Epub 2013/07/31. pmid:23891471; PubMed Central PMCID: PMCPMC3738835.
- 83. Liu L, Luo H. Whole-Exome Sequencing Identified a Novel Compound Heterozygous Mutation of LRRC6 in a Chinese Primary Ciliary Dyskinesia Patient. Biomed Res Int. 2018;2018:1854269. Epub 2018/03/08. pmid:29511670; PubMed Central PMCID: PMCPMC5817365.
- 84. Sha Y, Yang X, Mei L, Ji Z, Wang X, Ding L, et al. DNAH1 gene mutations and their potential association with dysplasia of the sperm fibrous sheath and infertility in the Han Chinese population. Fertil Steril. 2017;107(6):1312–8 e2. Epub 2017/06/05. pmid:28577616.
- 85. Coutton C, Vargas AS, Amiri-Yekta A, Kherraf ZE, Ben Mustapha SF, Le Tanno P, et al. Mutations in CFAP43 and CFAP44 cause male infertility and flagellum defects in Trypanosoma and human. Nat Commun. 2018;9(1):686. Epub 2018/02/17. pmid:29449551; PubMed Central PMCID: PMCPMC5814398.