Search

Abstract Background Genome sequence variants affecting complex traits (quantitative trait loci, QTL) are enriched in functional regions of the genome, such as those marked by certain histone modifications. These variants are believed to influence gene expression. However, due to the linkage disequilibrium among nearby variants, pinpointing the precise location of QTL is challenging. We aimed to identify allele-specific binding (ASB) QTL (asbQTL) that cause variation in the level of histone modification, as measured by the height of peaks assayed by ChIP-seq (chromatin immunoprecipitation sequencing). We identified DNA sequences that predict the difference between alleles in ChIP-seq peak height in H3K4me3 and H3K27ac histone modifications in the mammary glands of cows. Results We used a gapped k-mer support vector machine, a novel best linear unbiased prediction model, and a multiple linear regression model that combines the other two approaches to predict variant impacts on peak height. For each method, a subset of 1000 sites with the highest magnitude of predicted ASB was considered as candidate asbQTL. The accuracy of this prediction was measured by the proportion where the predicted direction matched the observed direction. Prediction accuracy ranged between 0.59 and 0.74, suggesting that these 1000 sites are enriched for asbQTL. Using independent data, we investigated functional enrichment in the candidate asbQTL set and three control groups, including non-causal ASB sites, non-ASB variants under a peak, and SNPs (single nucleotide polymorphisms) not under a peak. For H3K4me3, a higher proportion of the candidate asbQTL were confirmed as ASB when compared to the non-causal ASB sites (P 0.05). Conclusions We identified putatively causal sites for asbQTL using the DNA sequence surrounding these sites. Our results suggest that many sites influencing histone modifications may not directly affect gene expression. However, it is important to acknowledge that distinguishing between putative causal ASB sites and other non-causal ASB sites in high linkage disequilibrium with the causal sites regarding their impact on gene expression may be challenging due to limitations in statistical power.

Abstract Background Female fertility is an important trait in dairy cattle. Identifying putative causal variants associated with fertility may help to improve the accuracy of genomic prediction of fertility. Combining expression data (eQTL) of genes, exons, gene splicing and allele specific expression is a promising approach to fine map QTL to get closer to the causal mutations. Another approach is to identify genomic differences between cows selected for high and low fertility and a selection experiment in New Zealand has created exactly this resource. Our objective was to combine multiple types of expression data, fertility traits and allele frequency in high- (POS) and low-fertility (NEG) cows with a genome-wide association study (GWAS) on calving interval in Australian cows to fine-map QTL associated with fertility in both Australia and New Zealand dairy cattle populations. Results Variants that were significantly associated with calving interval (CI) were strongly enriched for variants associated with gene, exon, gene splicing and allele-specific expression, indicating that there is substantial overlap between QTL associated with CI and eQTL. We identified 671 genes with significant differential expression between POS and NEG cows, with the largest fold change detected for the CCDC196 gene on chromosome 10. Our results provide numerous candidate genes associated with female fertility in dairy cattle, including GYS2 and TIGAR on chromosome 5 and SYT3 and HSD17B14 on chromosome 18. Multiple QTL regions were located in regions with large numbers of copy number variants (CNV). To identify the causal mutations for these variants, long read sequencing may be useful. Conclusions Variants that were significantly associated with CI were highly enriched for eQTL. We detected 671 genes that were differentially expressed between POS and NEG cows. Several QTL detected for CI overlapped with eQTL, providing candidate genes for fertility in dairy cattle.

ABSTRACT: The longevity of dairy cattle has economic, animal welfare, and health implications and is influenced by the frequency of mortality on the farm and sale for slaughter. In this study cows removed from the herd due to death or slaughter during the lactation were coded 1 and cows that were not terminated were coded 0. Genetic parameters for mortality rates (MR) and slaughter rates (SR) were estimated for Holstein (H) and Jersey (J) breeds by applying both linear (LM) and threshold (TM) sire models using about 1.2 million H and 286,000 J cows. Estimated breeding values (EBV) for MR and SR were predicted using animal models to assess the opportunity for selection and genetic trends. Cow termination data, recorded between 1990 and 2020 on a voluntary basis by Australian dairy farmers, were analyzed. Cow MR has increased from below 1% in the 1990s to 4.1% and 3.6% in recent years in H and J cows, respectively. Most dead cows (∼36%) left the herd before 120 d of lactation, while cows that were slaughtered left the herd toward the end of the lactation. Using the LM, heritability (h2) estimates for MR were lower (1%) than those for SR (2%–3.5%). When h2 were estimated using a TM, the estimates for both traits varied between 4% and 20%, suggesting that the difference in incidence level is one of the reasons for the difference in the h2 values between MR and SR. Early test-day milk yield (MY) and 305-d MY (305-d MY) have unfavorable genetic correlations (0.32–0.41) with MR in both breeds. The genetic correlations of calving interval with MR were stronger (0.54–0.68) than with SR (0.28–0.45) suggesting that poor fertility can serve as an early indicator of poor cow health that may lead to increased risk of death. High early test-day somatic cell count is genetically associated with increased likelihood of slaughter (0.24–0.46), but not with increased likelihood of death. In H, 305-d protein yield (PY) had the strongest genetic correlation (−0.34 to −0.40) with SR whereas in J, both 305-d PY and fat yield showed high genetic (−0.64 to −0.70) and moderate environmental (−0.35 to −0.37) correlations with SR. The genetic correlation of removal from the herd due to death and slaughter was negative (−0.3) in J and zero in H. Strong selection for improved fertility and survival and less selection emphasis for MY, has led to an improvement in the genetic trend for cow MR in H and the trend in J has stabilized. Although genetic evaluations for cow MR are feasible, the reliabilities of the EBV are low and the level of cow MR in Australia are relatively low compared with similar countries. Therefore, genetic evaluation for survival based on mortality and slaughter data could be sufficient in the current selection circumstances where breeding objectives are broadly defined. Nevertheless, all Australian farmers should be encouraged to continue recording mortality and slaughter data for monitoring of the trends and for future development of genetic evaluations.

Abstract Background Many phenotypes in animal breeding are derived from incomplete measures, especially if they are challenging or expensive to measure precisely. Examples include time-dependent traits such as reproductive status, or lifespan. Incomplete measures for these traits result in phenotypes that are subject to left-, interval- and right-censoring, where phenotypes are only known to fall below an upper bound, between a lower and upper bound, or above a lower bound respectively. Here we compare three methods for deriving phenotypes from incomplete data using age at first elevation (> 1 ng/mL) in blood plasma progesterone (AGEP4), which generally coincides with onset of puberty, as an example trait. Methods We produced AGEP4 phenotypes from three blood samples collected at about 30-day intervals from approximately 5,000 Holstein–Friesian or Holstein–Friesian × Jersey cross-bred dairy heifers managed in 54 seasonal-calving, pasture-based herds in New Zealand. We used these actual data to simulate 7 different visit scenarios, increasing the extent of censoring by disregarding data from one or two of the three visits. Three methods for deriving phenotypes from these data were explored: 1) ordinal categorical variables which were analysed using categorical threshold analysis; 2) continuous variables, with a penalty of 31 d assigned to right-censored phenotypes; and 3) continuous variables, sampled from within a lower and upper bound using a data augmentation approach. Results Credibility intervals for heritability estimations overlapped across all methods and visit scenarios, but estimated heritabilities tended to be higher when left censoring was reduced. For sires with at least 5 daughters, the correlations between estimated breeding values (EBVs) from our three-visit scenario and each reduced data scenario varied by method, ranging from 0.65 to 0.95. The estimated breed effects also varied by method, but breed differences were smaller as phenotype censoring increased. Conclusion Our results indicate that using some methods, phenotypes derived from one observation per offspring for a time-dependent trait such as AGEP4 may provide comparable sire rankings to three observations per offspring. This has implications for the design of large-scale phenotyping initiatives where animal breeders aim to estimate variance parameters and estimated breeding values (EBVs) for phenotypes that are challenging to measure or prohibitively expensive.

Abstract Background The reduction in phenotypic performance of a population due to mating between close relatives is called inbreeding depression. The genetic background of inbreeding depression for semen traits is poorly understood. Thus, the objectives were to estimate the effect of inbreeding and to identify genomic regions underlying inbreeding depression of semen traits including ejaculate volume (EV), sperm concentration (SC), and sperm motility (SM). The dataset comprised ~ 330 K semen records from ~ 1.5 K Holstein bulls genotyped with 50 K single nucleotide polymorphism (SNP) BeadChip. Genomic inbreeding coefficients were estimated using runs of homozygosity (i.e., F ROH > 1 Mb) and excess of SNP homozygosity (F SNP). The effect of inbreeding was estimated by regressing phenotypes of semen traits on inbreeding coefficients. Associated variants with inbreeding depression were also detected by regressing phenotypes on ROH state of the variants. Results Significant inbreeding depression was observed for SC and SM (p

ABSTRACT: Achieving an acceptable level of fertility in herds is difficult for many dairy producers because identifying cows in estrus has become challenging owing to poor estrus expression, increased herd size, and lack of time and skilled labor for estrus detection. As a result, synchronization of estrus is often used to manage reproduction. The aims of this study were (1) to identify artificial inseminations (AI) that were performed following synchronization and (2) to assess the effect of synchronization on genetic parameters and evaluation of fertility traits. This study used breeding data collected between 1995 and 2021 from over 4,600 Australian dairy herds that had at least 30 matings per year. Because breeding methods were not reported, the recording pattern of breeding dates showing a large proportion of the total AI being recorded on a single date of the year served as an indicator of synchronization. First, the proportion of AI recorded on each day of the year was calculated for each herd-year. Subsequently, synchronization was defined when a herd with, for instance, only 30 matings in a year, had at least 0.20 or more AI on the same day. As the number of breedings in a herd-year increased, the threshold for classifying AI was continuously reduced from 0.20 to as low as 0.03 under the assumption that mating of many cows on a single date becomes increasingly difficult without synchronization. From the current data, we deduced that 0.11 of all AI were possibly performed following synchronization (i.e., timed AI, TAI). The proportion of AI classified as TAI increased over time and with herd size. Although the deviation from equal numbers of mating on 7 d of the week was not used for classifying AI, 0.44 of AI being categorized as TAI were performed on just 2 d of the week. When data classified as TAI were used for estimating genetic parameters and breeding values, the interval between calving and first service (CFS) was found to be the most affected trait. The phenotypic and additive genetic variance and heritability, as well as variability and reliability of estimated breeding values of bulls and cows for CFS were lower for TAI than for AI performed following detected estrus (i.e., estrus-detected AI, EAI). For calving interval, first service nonreturn rate (FNRR), and successful calving rate to first service, genetic correlations between the same trait measured in TAI and EAI were close to 1, in contrast to 0.55 for CFS. The lower genetic variances and heritabilities for FNRR and calving interval in TAI than in EAI suggests that synchronization reduces the genetic variability of fertility. In conclusion, TAI makes CFS an ineffective measure of fertility. One approach to minimize this effect on genetic evaluations is to identify TAI (using the method described for example) and then set the CFS of these cows as missing records when running multitrait genetic evaluations of fertility traits that include CFS. In the long term, the most practical and accurate way to reduce the effect of synchronization on genetic evaluations is to record TAI along with mating data.

Abstract Background Sharing individual phenotype and genotype data between countries is complex and fraught with potential errors, while sharing summary statistics of genome-wide association studies (GWAS) is relatively straightforward, and thus would be especially useful for traits that are expensive or difficult-to-measure, such as feed efficiency. Here we examined: (1) the sharing of individual cow data from international partners; and (2) the use of sequence variants selected from GWAS of international cow data to evaluate the accuracy of genomic estimated breeding values (GEBV) for residual feed intake (RFI) in Australian cows. Results GEBV for RFI were estimated using genomic best linear unbiased prediction (GBLUP) with 50k or high-density single nucleotide polymorphisms (SNPs), from a training population of 3797 individuals in univariate to trivariate analyses where the three traits were RFI phenotypes calculated using 584 Australian lactating cows (AUSc), 824 growing heifers (AUSh), and 2526 international lactating cows (OVE). Accuracies of GEBV in AUSc were evaluated by either cohort-by-birth-year or fourfold random cross-validations. GEBV of AUSc were also predicted using only the AUS training population with a weighted genomic relationship matrix constructed with SNPs from the 50k array and sequence variants selected from a meta-GWAS that included only international datasets. The genomic heritabilities estimated using the AUSc, OVE and AUSh datasets were moderate, ranging from 0.20 to 0.36. The genetic correlations (rg) of traits between heifers and cows ranged from 0.30 to 0.95 but were associated with large standard errors. The mean accuracies of GEBV in Australian cows were up to 0.32 and almost doubled when either overseas cows, or both overseas cows and AUS heifers were included in the training population. They also increased when selected sequence variants were combined with 50k SNPs, but with a smaller relative increase. Conclusions The accuracy of RFI GEBV increased when international data were used or when selected sequence variants were combined with 50k SNP array data. This suggests that if direct sharing of data is not feasible, a meta-analysis of summary GWAS statistics could provide selected SNPs for custom panels to use in genomic selection programs. However, since this finding is based on a small cross-validation study, confirmation through a larger study is recommended.

IntroductionThe automated collection of phenotypic measurements in livestock is of interest to both researchers and farmers. Real-time, low-cost, and accurate phenotyping can enhance precision livestock management and could lead to the optimized utilization of pasture and breeding of efficient animals. Wearable sensors provide the tools for researchers to develop novel phenotypes across all production systems, which is especially valuable for grazing conditions. The objectives of this study were to estimate the repeatability and heritability of traits related to grazing and rumination activities and their correlations with other traits.MethodsThis study was conducted on a commercial Merino farm in the west of Victoria, Australia, from 4 May 2020 to 29 May 2020. A total of 160 ActiGraph sensors embedded in halters were attached to the left side of the muzzles of Merino sheep (M = 74, F = 86) aged 10–11 months while the sheep were grazing on pasture. Support vector machine (SVM) algorithms classified the sensor output into the categories of grazing, rumination, walking, idle, and other activities. These activities were further classified into daily grazing time (GT), number of grazing events (NGE), grazing length (GL), rumination time (RT), number of rumination events (NRE), rumination length (RL), walking time (WT), and idle time (IT). The data were analyzed using univariate and bivariate models in ASReml-SA to estimate the repeatability, heritability, and phenotypic correlations among traits.ResultsThe heritability of GT was estimated to be 0.44 ± 0.23, whereas the other traits had heritability estimates close to zero. The estimated repeatability for all traits was moderate to high, with the highest estimate being for GT (0.70 ± 0.03) and the lowest for RT (0.44 ± 0.03). The intraclass correlation or repeatability at a 1-day interval (i.e., 2 consecutive days) was high for all traits, and steadily reduced when the interval between measurements was longer than 1 week.DiscussionThe estimated repeatability for the grazing traits showed that wearable sensors and SVM methods are reliable methods for recording sheep activities on pasture, and have a potential application in the ranking of animals for selective breeding.

BayesR3 samples the polymorphisms affecting complex traits at reduced computational cost to predict the genetic value, breeding value, or individual risk of genotypes.

Selection decisions are generally based on estimated breeding values (EBV) for a combination of traits that are polygenic (e.g. milk production). However, in some cases, there is additional intense selection for a single allele, or SNP, for a perceived benefit, such as selection for polled or A2 milk. Using a dataset where the A2 mutation was imputed using a reference population with whole genome sequence, we tested the hypothesis that intense selection in Australian Holstein cattle for the A2 allele in the β-casein gene may have resulted in increased inbreeding. We also estimated the average difference in performance between animals homozygous for the A1 or A2 allele for a range of traits. Using high-density genotypes we compared differences in genome-wide and regional inbreeding between Holstein cows homozygous for either the A1 or A2 β-casein alleles i.e. A1/A1 or A2/A2. This study shows that between the years 2000 to 2017, the frequency of the A2/A2 genotype increased by 20% in Holstein cows (from 32% to 52%). Our results suggest that selection for homozygosity at the β-casein A2 allele has increased inbreeding both across the genome and on chromosome 6 in A2/A2 Holstein cows. Animals that were A2/A2 were twice as likely to have a run of homozygosity of at least 1Mb long across the β-casein locus compared to animals that were A1/A1. Cows that are homozygous for the A2 allele had an average protein yield EBV advantage of 0.24 genetic standard deviations (SD) compared to A1/A1 homozygous cows. In contrast, A2/A2 homozygous animals were on average 0.2 genetic SD inferior on fertility EBV. As a result, the difference in the overall economic index (that includes traits contributing to profitability) there was only a small advantage of 0.05 SD for A2/A2 cows compared to A1/A1 cows. However, strong selection for the A2 allele has likely led to a higher level of regional and overall inbreeding which in the long term could harm genetic progress for some or all economic traits. Therefore, applying approaches that mitigate rapid inbreeding while selecting for preferred alleles and quantitative traits may be desirable.

Abstract Blood urea nitrogen (BUN) is an indicator trait for urinary nitrogen excretion. Measuring BUN level requires a blood sample, which limits the number of records that can be obtained. Alternatively, BUN can be predicted using mid-infrared (MIR) spectroscopy of a milk sample and thus records become available on many more cows through routine milk recording processes. The genetic correlation between MIR predicted BUN (MBUN) and BUN is 0.90. Hence, genetically, BUN and MBUN can be considered as the same trait. The objective of our study was to perform genome-wide association studies (GWAS) for BUN and MBUN, compare these two GWAS and detect quantitative trait loci (QTL) for both traits, and compare the detected QTL with previously reported QTL for milk urea nitrogen (MUN). The dataset used for our analyses included 2098 and 18,120 phenotypes for BUN and MBUN, respectively, and imputed whole-genome sequence data. The GWAS for MBUN was carried out using either the full dataset, the 2098 cows with records for BUN, or 2000 randomly selected cows, so that the dataset size is comparable to that for BUN. The GWAS results for BUN and MBUN were very different, in spite of the strong genetic correlation between the two traits. We detected 12 QTL for MBUN, on bovine chromosomes 2, 3, 9, 11, 12, 14 and X, and one QTL for BUN on chromosome 13. The QTL detected on chromosomes 11, 14 and X overlapped with QTL detected for MUN. The GWAS results were highly sensitive to the subset of records used. Hence, caution is warranted when interpreting GWAS based on small datasets, such as for BUN. MBUN may provide an attractive alternative to perform a more powerful GWAS to detect QTL for BUN.

Abstract Background Urinary nitrogen leakage is an environmental concern in dairy cattle. Selection for reduced urinary nitrogen leakage may be done using indicator traits such as milk urea nitrogen (MUN). The result of a previous study indicated that the genetic correlation between MUN in Australia (AUS) and MUN in New Zealand (NZL) was only low to moderate (between 0.14 and 0.58). In this context, an alternative is to select sequence variants based on genome-wide association studies (GWAS) with a view to improve genomic prediction accuracies. A GWAS can also be used to detect quantitative trait loci (QTL) associated with MUN. Therefore, our objectives were to perform within-country GWAS and a meta-GWAS for MUN using records from up to 33,873 dairy cows and imputed whole-genome sequence data, to compare QTL detected in the GWAS for MUN in AUS and NZL, and to use sequence variants selected from the meta-GWAS to improve the prediction accuracy for MUN based on a joint AUS-NZL reference set. Results Using the meta-GWAS, we detected 14 QTL for MUN, located on chromosomes 1, 6, 11, 14, 19, 22, 26 and the X chromosome. The three most significant QTL encompassed the casein genes on chromosome 6, PAEP on chromosome 11 and DGAT1 on chromosome 14. We selected 50,000 sequence variants that had the same direction of effect for MUN in AUS and MUN in NZL and that were most significant in the meta-analysis for the GWAS. The selected sequence variants yielded a genetic correlation between MUN in AUS and MUN in NZL of 0.95 and substantially increased prediction accuracy in both countries. Conclusions Our results demonstrate how the sharing of data between two countries can increase the power of a GWAS and increase the accuracy of genomic prediction using a multi-country reference population and sequence variants selected based on a meta-GWAS.

Abstract Background Heat tolerance is a trait of economic importance in the context of warm climates and the effects of global warming on livestock production, reproduction, health, and well-being. This study investigated the improvement in prediction accuracy for heat tolerance when selected sets of sequence variants from a large genome-wide association study (GWAS) were combined with a standard 50k single nucleotide polymorphism (SNP) panel used by the dairy industry. Methods Over 40,000 dairy cattle with genotype and phenotype data were analysed. The phenotypes used to measure an individual’s heat tolerance were defined as the rate of decline in milk production traits with rising temperature and humidity. We used Holstein and Jersey cows to select sequence variants linked to heat tolerance. The prioritised sequence variants were the most significant SNPs passing a GWAS p-value threshold selected based on sliding 100-kb windows along each chromosome. We used a bull reference set to develop the genomic prediction equations, which were then validated in an independent set of Holstein, Jersey, and crossbred cows. Prediction analyses were performed using the BayesR, BayesRC, and GBLUP methods. Results The accuracy of genomic prediction for heat tolerance improved by up to 0.07, 0.05, and 0.10 units in Holstein, Jersey, and crossbred cows, respectively, when sets of selected sequence markers from Holstein cows were added to the 50k SNP panel. However, in some scenarios, the prediction accuracy decreased unexpectedly with the largest drop of − 0.10 units for the heat tolerance fat yield trait observed in Jersey cows when 50k plus pre-selected SNPs from Holstein cows were used. Using pre-selected SNPs discovered on a combined set of Holstein and Jersey cows generally improved the accuracy, especially in the Jersey validation. In addition, combining Holstein and Jersey bulls in the reference set generally improved prediction accuracy in most scenarios compared to using only Holstein bulls as the reference set. Conclusions Informative sequence markers can be prioritised to improve the genomic prediction of heat tolerance in different breeds. In addition to providing biological insight, these variants could also have a direct application for developing customized SNP arrays or can be used via imputation in current industry SNP panels.

Abstract While understanding the genetic basis of heat tolerance is crucial in the context of global warming’s effect on humans, livestock, and wildlife, the specific genetic variants and biological features that confer thermotolerance in animals are still not well characterized. We used dairy cows as a model to study heat tolerance because they are lactating, and therefore often prone to thermal stress. The data comprised almost 0.5 million milk records (milk, fat, and proteins) of 29,107 Australian Holsteins, each having around 15 million imputed sequence variants. Dairy animals often reduce their milk production when temperature and humidity rise; thus, the phenotypes used to measure an individual’s heat tolerance were defined as the rate of milk production decline (slope traits) with a rising temperature–humidity index. With these slope traits, we performed a genome-wide association study (GWAS) using different approaches, including conditional analyses, to correct for the relationship between heat tolerance and level of milk production. The results revealed multiple novel loci for heat tolerance, including 61 potential functional variants at sites highly conserved across 100 vertebrate species. Moreover, it was interesting that specific candidate variants and genes are related to the neuronal system (ITPR1, ITPR2, and GRIA4) and neuroactive ligand–receptor interaction functions for heat tolerance (NPFFR2, CALCR, and GHR), providing a novel insight that can help to develop genetic and management approaches to combat heat stress.

Abstract Background Fourier-transform mid-infrared (FT-MIR) spectroscopy provides a high-throughput and inexpensive method for predicting milk composition and other novel traits from milk samples. While there have been many genome-wide association studies (GWAS) conducted on FT-MIR predicted traits, there have been few GWAS for individual FT-MIR wavenumbers. Using imputed whole-genome sequence for 38,085 mixed-breed New Zealand dairy cattle, we conducted GWAS on 895 individual FT-MIR wavenumber phenotypes, and assessed the value of these direct phenotypes for identifying candidate causal genes and variants, and improving our understanding of the physico-chemical properties of milk. Results Separate GWAS conducted for each of 895 individual FT-MIR wavenumber phenotypes, identified 450 1-Mbp genomic regions with significant FT-MIR wavenumber QTL, compared to 246 1-Mbp genomic regions with QTL identified for FT-MIR predicted milk composition traits. Use of mammary RNA-seq data and gene annotation information identified 38 co-localized and co-segregating expression QTL (eQTL), and 31 protein-sequence mutations for FT-MIR wavenumber phenotypes, the latter including a null mutation in the ABO gene that has a potential role in changing milk oligosaccharide profiles. For the candidate causative genes implicated in these analyses, we examined the strength of association between relevant loci and each wavenumber across the mid-infrared spectrum. This revealed shared association patterns for groups of genomically-distant loci, highlighting clusters of loci linked through their biological roles in lactation and their presumed impacts on the chemical composition of milk. Conclusions This study demonstrates the utility of FT-MIR wavenumber phenotypes for improving our understanding of milk composition, presenting a larger number of QTL and putative causative genes and variants than found from FT-MIR predicted composition traits. Examining patterns of significance across the mid-infrared spectrum for loci of interest further highlighted commonalities of association, which likely reflects the physico-chemical properties of milk constituents.

Time-dependent traits are often subject to censorship, where instead of precise phenotypes, only a lower and/or upper bound can be established for some of the individuals. Censorship reduces the precision of phenotypes but can represent compromise between measurement cost and animal ethics considerations. This compromise is particularly relevant for genetic evaluation because phenotyping initiatives often involve thousands of individuals. This research aimed to: 1) demonstrate a data augmentation approach for analysing censored phenotypes, and 2) quantify the implications of phenotype censorship on estimation of heritabilities and predictions of breeding values. First, we simulated uncensored phenotypes, representing fine-scale “age at puberty” for each individual in a population of some 5,000 animals across 50 herds. Analysis of these uncensored phenotypes provided a gold-standard control. We then produced seven “test” phenotypes by superimposing varying degrees of left, interval, and/or right censorship, as if herds were measured on only one, two or three occasions, with a binary measure categorized for animals at each visit (either pre or post pubertal). We demonstrated that our estimates of heritabilities and predictions of breeding values obtained using a data augmentation approach were remarkably robust to phenotype censorship. Our results have important practical implications for measuring time-dependent traits for genetic evaluation. More specifically, we suggest that data collection can be designed with relatively infrequent repeated measures, thereby reducing costs and increasing feasibility across large numbers of animals.

Heat tolerance is the ability of an animal to maintain production and reproduction levels under hot and humid conditions and is now a trait of economic relevance in dairy systems worldwide because of an escalating warming climate. The Australian dairy population is one of the excellent study models for enhancing our understanding of the biology of heat tolerance because they are predominantly kept outdoors on pastures where they experience direct effects of weather elements (e.g., solar radiation). In this article, we focus on evidence from recent studies in Australia that leveraged large a dataset [∼40,000 animals with phenotypes and 15 million whole-genome sequence variants] to elucidate the genetic basis of thermal stress as a critical part of the strategy to breed cattle adapted to warmer environments. Genotype-by-environment interaction (i.e., G × E) due to temperature and humidity variation is increasing, meaning animals are becoming less adapted (i.e., more sensitive) to changing environments. There are opportunities to reverse this trend and accelerate adaptation to warming climate by 1) selecting robust or heat-resilient animals and 2) including resilience indicators in breeding goals. Candidate causal variants related to the nervous system and metabolic functions are relevant for heat tolerance and, therefore, key for improving this trait. This could include adding these variants in the custom SNP panels used for routine genomic evaluations or as the basis to design specific agonist or antagonist compounds for lowering core body temperature under heat stress conditions. Indeed, it was encouraging to see that adding prioritized functionally relevant variants into the 50k SNP panel (i.e., the industry panel used for genomic evaluation in Australia) increased the prediction accuracy of heat tolerance by up to 10% units. This gain in accuracy is critical because genetic improvement has a linear relationship with prediction accuracy. Overall, while this article used data mainly from Australia, this could benefit other countries that aim to develop breeding values for heat tolerance, considering that the warming climate is becoming a topical issue worldwide.

Previous studies have shown reduced enteric methane emissions (ME) and residual feed intake (RFI) through the application of genomic selection in ruminants. The objective of this study was to evaluate feeding behaviour traits as genetic indicators for ME and RFI in Australian Maternal Composite ewes using data from an automated feed intake facility. The feeding behaviour traits evaluated were the amount of time spent eating per day (eating time; ETD; min/day) and per visit (eating time per event; ETE; min/event), daily number of events (DNE), event feed intake (EFI; g/event) and eating rate (ER; g/min). Genotypes and phenotypes of 445 ewes at three different ages (post-weaning, hogget, and adult) were used to estimate the heritability of ME, RFI, and the feeding behaviour traits using univariate genomic best linear unbiased prediction models. Multivariate models were used to estimate the correlations between these traits and within each trait at different ages. The response to selection was evaluated for ME and RFI with direct selection models and indirect models with ETE as an indicator trait, as this behaviour trait was a promising indicator based on heritability and genetic correlations. Heritabilities were between 0.12 and 0.18 for ME and RFI, and between 0.29 and 0.47 for the eating behaviour traits. In our data, selecting for more efficient animals (low RFI) would lead to higher methane emissions per day and per kg of dry matter intake. Selecting for more ETE also improves feed efficiency but results in more methane per day and per kg dry matter intake. Based on our results, ETE could be evaluated as an indicator trait for ME and RFI under an index approach that allows simultaneous selection for improvement in emissions and feed efficiency. Selecting for ETE may have a tremendous impact on the industry, as it may be easier and cheaper to obtain than feed intake and ME data. As the data were collected using individual feeding units, the findings on this research should be validated under grazing conditions.

This study assessed the accuracy and bias of genomic prediction (GP) in purebred Holstein (H) and Jersey (J) as well as crossbred (H and J) validation cows using different reference sets and prediction strategies. The reference sets were made up of different combinations of 36,695 H and J purebreds and crossbreds. Additionally, the effect of using different sets of marker genotypes on GP was studied (conventional panel: 50k, custom panel enriched with, or close to, causal mutations: XT_50k, and conventional high-density with a limited custom set: pruned HDnGBS). We also compared the use of genomic best linear unbiased prediction (GBLUP) and Bayesian (emBayesR) models, and the traits tested were milk, fat, and protein yields. On average, by including crossbred cows in the reference population, the prediction accuracies increased by 0.01–0.08 and were less biased (regression coefficient closer to 1 by 0.02–0.16), and the benefit was greater for crossbreds compared to purebreds. The accuracy of prediction increased by 0.02 using XT_50k compared to 50k genotypes without affecting the bias. Although using pruned HDnGBS instead of 50k also increased the prediction accuracy by about 0.02, it increased the bias for purebred predictions in emBayesR models. Generally, emBayesR outperformed GBLUP for prediction accuracy when using 50k or pruned HDnGBS genotypes, but the benefits diminished with XT_50k genotypes. Crossbred predictions derived from a joint pure H and J reference were similar in accuracy to crossbred predictions derived from the two separate purebred reference sets and combined proportional to breed composition. However, the latter approach was less biased by 0.13. Most interestingly, using an equalized breed reference instead of an H-dominated reference, on average, reduced the bias of prediction by 0.16–0.19 and increased the accuracy by 0.04 for crossbred and J cows, with a little change in the H accuracy. In conclusion, we observed improved genomic predictions for both crossbreds and purebreds by equalizing breed contributions in a mixed breed reference that included crossbred cows. Furthermore, we demonstrate, that compared to the conventional 50k or high-density panels, our customized set of 50k sequence markers improved or matched the prediction accuracy and reduced bias with both GBLUP and Bayesian models.

Abstract Background The mutations changing the expression level of a gene, or expression quantitative trait loci (eQTL), can be identified by testing the association between genetic variants and gene expression in multiple individuals (eQTL mapping), or by comparing the expression of the alleles in a heterozygous individual (allele specific expression or ASE analysis). The aims of the study were to find and compare ASE and local eQTL in 4 bovine RNA-sequencing (RNA-Seq) datasets, validate them in an independent ASE study and investigate if they are associated with complex trait variation. Results We present a novel method for distinguishing between ASE driven by polymorphisms in cis and parent of origin effects. We found that single nucleotide polymorphisms (SNPs) driving ASE are also often local eQTL and therefore presumably cis eQTL. These SNPs often, but not always, affect gene expression in multiple tissues and, when they do, the allele increasing expression is usually the same. However, there were systematic differences between ASE and local eQTL and between tissues and breeds. We also found that SNPs significantly associated with gene expression (p

Abstract Background Topological association domains (TADs) are chromosomal domains characterised by frequent internal DNA-DNA interactions. The transcription factor CTCF binds to conserved DNA sequence patterns called CTCF binding motifs to either prohibit or facilitate chromosomal interactions. TADs and CTCF binding motifs control gene expression, but they are not yet well defined in the bovine genome. In this paper, we sought to improve the annotation of bovine TADs and CTCF binding motifs, and assess whether the new annotation can reduce the search space for cis-regulatory variants. Results We used genomic synteny to map TADs and CTCF binding motifs from humans, mice, dogs and macaques to the bovine genome. We found that our mapped TADs exhibited the same hallmark properties of those sourced from experimental data, such as housekeeping genes, transfer RNA genes, CTCF binding motifs, short interspersed elements, H3K4me3 and H3K27ac. We showed that runs of genes with the same pattern of allele-specific expression (ASE) (either favouring paternal or maternal allele) were often located in the same TAD or between the same conserved CTCF binding motifs. Analyses of variance showed that when averaged across all bovine tissues tested, TADs explained 14% of ASE variation (standard deviation, SD: 0.056), while CTCF explained 27% (SD: 0.078). Furthermore, we showed that the quantitative trait loci (QTLs) associated with gene expression variation (eQTLs) or ASE variation (aseQTLs), which were identified from mRNA transcripts from 141 lactating cows’ white blood and milk cells, were highly enriched at putative bovine CTCF binding motifs. The linearly-furthermost, and most-significant aseQTL and eQTL for each genic target were located within the same TAD as the gene more often than expected (Chi-Squared test P-value

Abstract Background Using whole genome sequence data might improve genomic prediction accuracy, when compared with high-density SNP arrays, and could lead to identification of casual mutations affecting complex traits. For some traits, the most accurate genomic predictions are achieved with non-linear Bayesian methods. However, as the number of variants and the size of the reference population increase, the computational time required to implement these Bayesian methods (typically with Monte Carlo Markov Chain sampling) becomes unfeasibly long. Results Here, we applied a new method, HyB_BR (for Hybrid BayesR), which implements a mixture model of normal distributions and hybridizes an Expectation-Maximization (EM) algorithm followed by Markov Chain Monte Carlo (MCMC) sampling, to genomic prediction in a large dairy cattle population with imputed whole genome sequence data. The imputed whole genome sequence data included 994,019 variant genotypes of 16,214 Holstein and Jersey bulls and cows. Traits included fat yield, milk volume, protein kg, fat% and protein% in milk, as well as fertility and heat tolerance. HyB_BR achieved genomic prediction accuracies as high as the full MCMC implementation of BayesR, both for predicting a validation set of Holstein and Jersey bulls (multi-breed prediction) and a validation set of Australian Red bulls (across-breed prediction). HyB_BR had a ten fold reduction in compute time, compared with the MCMC implementation of BayesR (48 hours versus 594 hours). We also demonstrate that in many cases HyB_BR identified sequence variants with a high posterior probability of affecting the milk production or fertility traits that were similar to those identified in BayesR. For heat tolerance, both HyB_BR and BayesR found variants in or close to promising candidate genes associated with this trait and not detected by previous studies. Conclusions The results demonstrate that HyB_BR is a feasible method for simultaneous genomic prediction and QTL mapping with whole genome sequence in large reference populations.

Abstract Background Enhancers are non-coding DNA sequences, which when they are bound by specific proteins increase the level of gene transcription. Enhancers activate unique gene expression patterns within cells of different types or under different conditions. Enhancers are key contributors to gene regulation, and causative variants that affect quantitative traits in humans and mice have been located in enhancer regions. However, in the bovine genome, enhancers as well as other regulatory elements are not yet well defined. In this paper, we sought to improve the annotation of bovine enhancer regions by using publicly available mammalian enhancer information. To test if the identified putative bovine enhancer regions are enriched with functional variants that affect milk production traits, we performed genome-wide association studies using imputed whole-genome sequence data followed by meta-analysis and enrichment analysis. Results We produced a library of candidate bovine enhancer regions by using publicly available bovine ChIP-Seq enhancer data in combination with enhancer data that were identified based on sequence homology with human and mouse enhancer databases. We found that imputed whole-genome sequence variants associated with milk production traits in 16,581 dairy cattle were enriched with enhancer regions that were marked by bovine-liver H3K4me3 and H3K27ac histone modifications from both permutation tests and gene set enrichment analysis. Enhancer regions that were identified based on sequence homology with human and mouse enhancer regions were not as strongly enriched with trait-associated sequence variants as the bovine ChIP-Seq candidate enhancer regions. The bovine ChIP-Seq enriched enhancer regions were located near genes and quantitative trait loci that are associated with pregnancy, growth, disease resistance, meat quality and quantity, and milk quality and quantity traits in dairy and beef cattle. Conclusions Our results suggest that sequence variants within enhancer regions that are located in bovine non-coding genomic regions contribute to the variation in complex traits. The level of enrichment was higher in bovine-specific enhancer regions that were identified by detecting histone modifications H3K4me3 and H3K27ac in bovine liver tissues than in enhancer regions identified by sequence homology with human and mouse data. These results highlight the need to use bovine-specific experimental data for the identification of enhancer regions.

We predicted dry matter intake of dairy cows using parity, week of lactation, milk yield, milk mid-infrared (MIR) spectrum, and MIR-based predictions of bodyweight, fat, protein, lactose, and fatty acids content in milk. The dataset comprised 10,711 samples of 534 dairy cows with a geographical diversity (Australia, Canada, Denmark, and Ireland). We set up partial least square (PLS) regressions with different constructs and a one-hidden-layer artificial neural network (ANN) using the highest contribution variables. In the ANN, we replaced the spectra with their projections to the 25 first PLS factors explaining 99% of the spectral variability to reduce the model complexity. Cow-independent 10 × 10-fold cross-validation (CV) achieved the best performance with root mean square errors (RMSECV) of 3.27 ± 0.08 kg for the PLS regression and 3.25 ± 0.13 kg for ANN. Although the available data were significantly different, we also performed a country-independent validation (CIV) to measure the models’ performance fairly. We found RMSECIV varying from 3.73 to 6.03 kg for PLS and 3.69 to 5.08 kg for ANN. Ultimately, based on the country-independent validation, we discussed the developed models’ performance with those achieved by the National Research Council’s equation.

Knowing the body weight (BW) of a cow at a specific moment or measuring its changes through time is of interest for management purposes. The current work aimed to validate the feasibility of predicting BW using the day in milk, parity, milk yield, and milk mid-infrared (MIR) spectrum from a multiple-country dataset and reduce the number of predictors to limit the risk of over-fitting and potentially improve its accuracy. The BW modeling procedure involved feature selections and herd-independent validation in identifying the most interesting subsets of predictors and then external validation of the models. From 1849 records collected in 9 herds from 360 Holstein cows, the best performing models achieved a root mean square error (RMSE) for the herd-independent validation between 52 ± 2.34 kg to 56 ± 3.16 kg, including from 5 to 62 predictors. Among these models, three performed remarkably well in external validation using an independent dataset (N = 4067), resulting in RMSE ranging from 52 to 56 kg. The results suggest that multiple optimal BW predictive models coexist due to the high correlations between adjacent spectral points.

Despite the fact that the sulphur hexafluoride (SF6) tracer technique was developed over 25 years ago to measure methane production from grazing and non-housed animals, no studies have specifically investigated whether ambient wind speed, temperature, relative humidity and rainfall influence the accuracy of the method. The aim of this research was to investigate how these weather factors influence the measurement of enteric methane production by the SF6 technique. Six different cohorts of dairy cows (40 per cohort) were kept outdoors and fed a common diet during spring in 3 consecutive years. Methane production from individual cows was measured daily over the last 5 days of each 32-day period. An automated weather station measured air temperature, wind speed, relative humidity and rainfall every 10 min. Regression analyses were used to relate the average daily wind speed, average daily temperature, average daily relative humidity and total daily rainfall measurements to dry matter intake, average daily methane production and methane yield of each cohort of cows. It was concluded that the modified SF6 technique can be used outdoors during a range of wind speeds, ambient temperatures, relative humidities and rainfall conditions without causing a significant effect on the measurement of methane production or methane yield of dairy cows.

Structural variations (SVs) are large DNA segments of deletions, duplications, copy number variations, inversions and translocations in a re-sequenced genome compared to a reference genome. They have been found to be associated with several complex traits in dairy cattle and could potentially help to improve genomic prediction accuracy of dairy traits. Imputation of SVs was performed in individuals genotyped with single-nucleotide polymorphism (SNP) panels without the expense of sequencing them. In this study, we generated 24,908 high-quality SVs in a total of 478 whole-genome sequenced Holstein and Jersey cattle. We imputed 4489 SVs with R2 > 0.5 into 35,568 Holstein and Jersey dairy cattle with 578,999 SNPs with two pipelines, FImpute and Eagle2.3-Minimac3. Genome-wide association studies for production, fertility and overall type with these 4489 SVs revealed four significant SVs, of which two were highly linked to significant SNP. We also estimated the variance components for SNP and SV models for these traits using genomic best linear unbiased prediction (GBLUP). Furthermore, we assessed the effect on genomic prediction accuracy of adding SVs to GBLUP models. The estimated percentage of genetic variance captured by SVs for production traits was up to 4.57% for milk yield in bulls and 3.53% for protein yield in cows. Finally, no consistent increase in genomic prediction accuracy was observed when including SVs in GBLUP.

Disorders of energy metabolism, which can result from a failure to adapt to the period of negative energy balance immediately after calving, have significant negative effects on the health, welfare and profitability of dairy cows. The most common biomarkers of energy balance in dairy cows are β-hydroxybutyrate (BHBA) and non-esterified fatty acids (NEFA). While elevated concentrations of these biomarkers are associated with similar negative health and production outcomes, the phenotypic and genetic correlations between them are weak. In this study, we used an untargeted 1H NMR metabolomics approach to investigate the serum metabolomic fingerprints of BHBA and NEFA. Serum samples were collected from 298 cows in early lactation (calibration dataset N = 248, validation N = 50). Metabolomic fingerprinting was done by regressing 1H NMR spectra against BHBA and NEFA concentrations (determined using colorimetric assays) using orthogonal partial least squares regression. Prediction accuracies were high for BHBA models, and moderately high for NEFA models (R2 of external validation of 0.88 and 0.75, respectively). We identified 16 metabolites that were significantly (variable importance of projection score > 1) correlated with the concentration of one or both biomarkers. These metabolites were primarily intermediates of energy, phospholipid, and/or methyl donor metabolism. Of the significant metabolites identified; (1) two (acetate and creatine) were positively correlated with BHBA but negatively correlated with NEFA, (2) nine had similar associations with both BHBA and NEFA, (3) two were correlated with only BHBA concentration, and (4) three were only correlated with NEFA concentration. Overall, our results suggest that BHBA and NEFA are indicative of similar metabolic states in clinically healthy animals, but that several significant metabolic differences exist that help to explain the weak correlations between them. We also identified several metabolites that may be useful intermediate phenotypes in genomic selection for improved metabolic health.

Most livestock metabolomic studies involve relatively small, homogenous populations of animals. However, livestock farming systems are non-homogenous, and large and more diverse datasets are required to ensure that biomarkers are robust. The aims of this study were therefore to (1) investigate the feasibility of using a large and diverse dataset for untargeted proton nuclear magnetic resonance (1H NMR) serum metabolomic profiling, and (2) investigate the impact of fixed effects (farm of origin, parity and stage of lactation) on the serum metabolome of early-lactation dairy cows. First, we used multiple linear regression to correct a large spectral dataset (707 cows from 13 farms) for fixed effects prior to multivariate statistical analysis with principal component analysis (PCA). Results showed that farm of origin accounted for up to 57% of overall spectral variation, and nearly 80% of variation for some individual metabolite concentrations. Parity and week of lactation had much smaller effects on both the spectra as a whole and individual metabolites (