PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of gseBioMed CentralBiomed Central Web Sitesearchsubmit a manuscriptregisterthis articleGenetics, Selection, Evolution : GSEJournal Front Page
 
Genet Sel Evol. 2017; 49: 44.
Published online 2017 May 12. doi:  10.1186/s12711-017-0319-0
PMCID: PMC5427631

Exploring the genetic architecture and improving genomic prediction accuracy for mastitis and milk production traits in dairy cattle by mapping variants to hepatic transcriptomic regions responsive to intra-mammary infection

Abstract

Background

A better understanding of the genetic architecture of complex traits can contribute to improve genomic prediction. We hypothesized that genomic variants associated with mastitis and milk production traits in dairy cattle are enriched in hepatic transcriptomic regions that are responsive to intra-mammary infection (IMI). Genomic markers [e.g. single nucleotide polymorphisms (SNPs)] from those regions, if included, may improve the predictive ability of a genomic model.

Results

We applied a genomic feature best linear unbiased prediction model (GFBLUP) to implement the above strategy by considering the hepatic transcriptomic regions responsive to IMI as genomic features. GFBLUP, an extension of GBLUP, includes a separate genomic effect of SNPs within a genomic feature, and allows differential weighting of the individual marker relationships in the prediction equation. Since GFBLUP is computationally intensive, we investigated whether a SNP set test could be a computationally fast way to preselect predictive genomic features. The SNP set test assesses the association between a genomic feature and a trait based on single-SNP genome-wide association studies. We applied these two approaches to mastitis and milk production traits (milk, fat and protein yield) in Holstein (HOL, n = 5056) and Jersey (JER, n = 1231) cattle. We observed that a majority of genomic features were enriched in genomic variants that were associated with mastitis and milk production traits. Compared to GBLUP, the accuracy of genomic prediction with GFBLUP was marginally improved (3.2 to 3.9%) in within-breed prediction. The highest increase (164.4%) in prediction accuracy was observed in across-breed prediction. The significance of genomic features based on the SNP set test were correlated with changes in prediction accuracy of GFBLUP (P < 0.05).

Conclusions

GFBLUP provides a framework for integrating multiple layers of biological knowledge to provide novel insights into the biological basis of complex traits, and to improve the accuracy of genomic prediction. The SNP set test might be used as a first-step to improve GFBLUP models. Approaches like GFBLUP and SNP set test will become increasingly useful, as the functional annotations of genomes keep accumulating for a range of species and traits.

Electronic supplementary material

The online version of this article (doi:10.1186/s12711-017-0319-0) contains supplementary material, which is available to authorized users.

Background

In general, genetic variation in complex or quantitative traits is considered to be governed by a large number of loci with small to moderate effects, which are individually undetectable by genome-wide association studies (GWAS) with stringent significance thresholds [15]. A better understanding of the genetic architecture that underlies complex traits (e.g. the distribution of causal variants and their effects) could improve the predictive ability of models [4, 69]. This would be beneficial for genomic prediction of disease risk in humans and for estimating genetic values in livestock and plant species of agricultural importance [4, 69].

The genomic best linear unbiased prediction (GBLUP) assumes that all genomic markers contribute equally to variability of a trait [10] and ignores any prior biological knowledge on genetic architecture of the trait. However, genomic markers that are associated with a complex trait may not be uniformly and randomly distributed over the genome, but rather be clustered in genes that are part of interconnected biological pathways and networks [2, 11, 12]. The genomic regions that are likely to be enriched in variants affecting a trait are defined as genomic features. Based on different biological hypotheses, genomic features can be defined from various sources of biological knowledge, such as genes, gene ontologies, biological pathways, or other types of external evidence. Incorporating this biological information may improve the predictive abilities of models. We extended the GBLUP model to implement this strategy by including a separate random effect for the joint action of single nucleotide polymorphisms (SNPs) within a genomic feature [8], which we call a genomic feature BLUP (GFBLUP) model. As a result, individual SNP relationships can be weighted differentially in GFBLUP according to the variance explained by SNPs within and outside the genomic feature [8]. The GFBLUP model has been applied to three complex traits (i.e. chill coma recovery, starvation resistance and startle response) in the unrelated inbred lines of Drosophila melanogaster populations [8]. Compared to GBLUP, the prediction accuracy with GFBLUP was substantially improved when incorporating several gene ontology (GO) categories as genomic features [8]. A possible increase in prediction accuracy with GFBLUP would depend on whether the genomic feature is enriched in causal mutations.

The GFBLUP model is computationally intensive for evaluating many genomic features [8]. Therefore, it is important to develop a computationally fast approach. The SNP set test based on GWAS-derived single-SNP test statistics could be one such approach. It would be of interest to investigate the relationship between the significance of a genomic feature based on the SNP set test and the predictive ability of the GFBLUP model.

To date, there are many genes that are yet neither functionally characterized nor mapped to any biological databases [1316], in particular in livestock populations. For example, in cattle only ~20% of the genes are annotated in Kyoto Encyclopaedia of Genes and Genomes (KEGG) pathways [17]. However, transcriptomics studies have been conducted on small-scale experimental populations to investigate the dynamic state of the transcriptome in particular tissues, revealing thousands of genomic features (e.g. genes and pathways) that are engaged in the biological processes of complex traits [1820]. Such transcriptomics studies provide tissue-specific genomic features that are likely to be enriched in genomic variants affecting specific traits.

Mastitis, an inflammatory condition of the mammary gland, is often caused by invading pathogens. It is the most costly disease in the dairy industry due to treatment cost, reduction in milk production and milk quality, and in some cases culling of the affected cows [21]. Gram-negative Escherichia coli (E. coli) is a common mastitis-causing bacteria [22], and the lipopolysaccharides (LPS) released by E. coli induce acute inflammatory responses [23]. Genes with expression levels that are significantly affected during the early stage of infection have also been suggested to be involved in overall metabolism [19, 2326]. Moreover, it is well established that mastitis is unfavorably correlated with milk production traits [25]. Since liver plays key roles in innate immune response and metabolic regulation [27], we hypothesized that hepatic transcriptomic regions that are responsive to intra-mammary infection (IMI) may be enriched in genomic variants that impact mastitis and milk production traits. Using these regions as genomic features might provide more predictive GFBLUP models compared to the GBLUP model. In addition, since gene expression patterns and molecular interaction networks are consistent across breeds [28], we further hypothesized that the use of transcriptomic data obtained on one breed may contribute to improve genomic prediction in other breeds.

In the current study, mastitis and three milk production traits (i.e. milk, fat and protein yield) from Nordic Holstein (HOL, n = 5056) and Jersey (JER, n = 1231) cattle were analyzed using imputed sequence genotype data (~15 million SNPs) and hepatic transcriptome data from an IMI study. Our main objectives were to apply the GFBLUP model and SNP set test: (1) to investigate the genomic variance explained by transcriptomic regions that are responsive to IMI; (2) to improve the accuracy of within-breed and across-breed genomic prediction using GFBLUP compared to GBLUP; and (3) to investigate the relationship between the predictive ability of GFBLUP and the significance of genomic features based on the SNP set test.

Methods

Intra-mammary infection (IMI) study

The IMI experimental design and collection of liver biopsies were reported previously [23, 29]. In brief, eight healthy HOL dairy cows in their first lactation (9 to 12 weeks after calving) were selected for the experiment. The udder quarters of all studied cows were free from mastitis pathogens based on bacteriological examinations. Milk somatic cell count (SCC) for each studied quarter was <100,000. The right front quarter was infected with 200 µg of E. coli LPS (0111:B4) (Sigma-Aldrich, Brøndby, Denmark) dissolved in 10 mL of a 0.9% NaCl solution, while the left front quarter was used as a control and challenged with 10 mL of 0.9% NaCl solution only. Clinical signs, data on production traits together with milk and blood parameters associated with LPS infection were recorded throughout the trial and confirmed that mastitis inflammation was induced. Liver biopsies collected 22 h before and 3, 6, 9, 12 and 48 h after LPS infection in three cows were used for RNA extraction. Sampling procedures for liver biopsies were described previously [30]. Finally, 18 RNA-Seq libraries (at each time point with three biological replicates) were sequenced using 100-bp paired-end sequencing in Illumina Hiseq2000 sequencing technology.

Statistical analysis of RNA-Seq data

Statistical approaches used for analysing RNA-Seq data were described previously [31]. Briefly, sequence reads of each sample were aligned to the bovine reference genome assembly (UMD 3.1), using a sensitive and efficient mapping program based on the seed-and-vote algorithm implemented in the Rsubread package in R/Bioconductor [32] _ENREF_65. The number of reads that were mapped to 24,616 Ensemble genes **(ftp://ftp.ensembl.org/pub/release-86/gtf/bos_taurus) was counted using the function Feature-Counts in the Rsubread package with default settings. The average mapping rate across all samples was approximately 68%. Analysis of differential gene expression was conducted using edgeR [33]. A small number of highly expressed genes in a sample can cause an RNA composition effect, i.e. a substantial proportion of the total library size could be consumed by these highly expressed genes, which results in the remaining genes to be under-sampled [33]. Therefore, the most recommended weighted trimmed means of M-values (TMM) were used to normalize the total count data (i.e. the total library size) between each pair of samples, in order to adjust for RNA composition effect [33]. After normalization of the total library size, a negative binomial generalized linear model (GLM) was applied for each gene, because the count data of genes follow non-normal distributions, which commonly exhibit a quadratic mean–variance relationship [33]. The relevant factors in the experimental design were also adjusted by the GLM, and gene differential expression was determined using a likelihood ratio test [33]. In the GLM model, where the number of reads mapped to gene g in sample i is denoted as ygi and the total number of mapped reads is denoted as Ni, it is assumed that ygi ∼  NB(μgiϕg), where μgi and ϕg are the location and the dispersion parameters of the negative binomial distribution, respectively. To ensure stable inference for each gene, an empirical Bayes method was used to compress gene-wise dispersions towards a common dispersion for all genes [33]. Statistical tests for each analysis were adjusted for multiple-testing using the FDR method as implemented in R (version 3.2.4).

Defining genomic features

The differentially-expressed genes (DEG) (i.e. the hepatic transcriptome regions responsive to IMI) that were obtained from the above RNA-Seq analyses were used to define genomic features. First, 30 genomic features were defined using six false discovery rate (FDR) cut-off values (i.e. ≤5×10−2, 10−2, 10−3, 10−6, 10−8, and 10−10) in each of the five experimental comparisons (i.e. 3 vs. −22 h, 6 vs. −22 h, 9 vs. −22 h, 12 vs. −22 h and 48 vs. −22 h), respectively. In addition, since the biological functions of up-regulated and down-regulated genes can be quite different, each of these 30 genomic features was further divided into four subsets based on four log2(fold-change)s cut-off values (i.e. ≤−2, ≤−1, ≤1, and >2). Therefore, another 115 genomic features were built, because five conditions were without DEG. In total, 145 genomic features were defined. The number of DEG in each genomic feature is summarized in Table S1 (see Additional file 1: Table S1).

Phenotypic data

The phenotypes were de-regressed breeding values (DRP) from routine genetic evaluations by the Nordic Cattle Genetic Evaluation (NAV, http://www.nordicebv.info/), and were available for 5056 HOL and 1231 JER cattle. Detailed information of these phenotypes was previously described in [34, 35]. Heritabilities for milk, fat and protein yields and mastitis were equal to 0.39, 0.39, 0.39 and 0.04, respectively in HOL, and very similar in JER [34, 35]. The average reliabilities of the DRP for milk, fat and protein yields and mastitis were equal to 0.95, 0.95, 0.95 and 0.83, respectively in HOL; and 0.92, 0.92, 0.92, and 0.76, respectively in JER.

Genotypic data

Imputation from Illumina BovineSNP50 BeadChip (50 K) to Illumina BovineHD BeadChip (high-density, HD) genotypes for these individuals and further to whole-genome sequence variants was described previously [36, 37]. Briefly, genotypes from the 50 K SNP chip for each individual were first imputed to HD genotypes using a multi-breed reference of 3383 animals (1222 HOL, 1326 Nordic Red, and 835 JER). A total of 648,219 SNPs were obtained after imputation to the HD chip. These imputed HD genotypes were then imputed to the whole-genome sequence level using a multi-breed reference population of 1228 individuals from Run4 of the 1000 Bull Genomes Project [38] and additional whole-genome sequences from Aarhus University including 368 HOL, 86 Nordic Red, and 88 JER individuals [39]. Genotype imputation was done using Minimac2 [40]. In total, 22,751,039 biallelic variants (SNPs and Indel) were included in the imputed sequence genotypic data. The accuracy of imputation was above 0.85 for the across-breed imputation of 19,498,365 SNPs. Detailed information about imputation accuracy was previously reported in [37]. For each breed, SNPs with a large deviation from Hardy–Weinberg proportions (P < 10−6) or with minor allele frequency (MAF) <0.01 were further excluded. A total of 15,355,382 and 13,403,916 SNPs remained for the HOL and JER datasets, respectively. The SNP locations were based on the UMD3.1 reference genome (http://www.ensembl.org/Bos_taurus/Info/Index). A SNP was considered to be linked with a genomic feature if its chromosome position was within the open reading frame of DEG in the particular genomic feature.

Training and validation populations

For within-breed prediction, each of the datasets (i.e. HOL and JER) was divided into training and validation sets based on birth-year of the animal to access prediction accuracy. The birth-year cut-off was 2006 for HOL and 2004 for JER, and the younger animals were assigned to the validation dataset (Table 1). We chose this validation strategy considering routine animal breeding practice where the young bulls breeding values are predicted using a training population of older animals. For across-breed prediction, the complete HOL population (n = 5056) was used as training data to predict breeding values for all JER bulls (n = 1231). Both GBLUP and GFBLUP models were fitted to compare prediction accuracies.

Table 1
Overview of training and validation population sizes for genomic predictions

Genomic models

For each genomic feature as defined before, SNPs were partitioned into two sets (i.e. within and outside the genomic feature), followed by the GFBLUP model analysis:

y = 1μgfg-fe

where y is the vector of phenotypic observations, 1 is a vector of 1s, μ is the overall mean, gf is the vector of genomic values captured by the SNPs within a genomic feature, g-f is the vector of genomic values captured by SNPs outside the genomic feature (i.e. the rest of genome), and e is the vector of residuals. Assumptions for all random effects are given by:

gfg-feN000,Gfσf2000G-fσ-f2000Dσe2,

where Gf and G-f are genomic relationship matrices that are built using the SNPs within and outside the genomic feature, respectively, which were calculated using the second method described in [41]. Briefly, let M be the marker matrix that specifies which alleles the individual inherits, and P be the matrix that contains the frequencies of the second allele at locus (pi) expressed as a difference from the 0.5 value and multiplied by 2, that is, the column i of P is 2(pi - 0.5). Matrix Z was obtained as MP, which allows mean values of the allele effects to be equal to 0. Then, GZTZ, where T is a diagonal matrix with Tii=1m2pi1-pi. D is a diagonal matrix with diagonal elements equal to 1-r2r2, where r2 is the reliability of DRP, σf2,σ-f2 and σe2 are the variance components accounted for by the SNPs within and outside the genomic feature, and by the residuals, respectively.

The standard GBLUP model includes only one random genomic effect:

y = 1μge

with the same notation as above except for g, which is the vector of genomic values captured by all genomic SNPs. The random genomic values and the residuals were assumed to be independently distributed: gN0,Gσg2 and eN0,Dσe2.

Estimation of genomic parameters

The variance components, σf2, σ-f2, σg2 and σe2, were estimated using an average information restricted maximum-likelihood (AI-REML) procedure [42] implemented in DMU [43]. The proportion of genomic variance explained by a genomic feature in the GFBLUP model: Hf2=σf2σf2+σ-f2. The proportion of phenotypic variance explained by all SNPs: hGFBLUP2=σf2+σ-f2σf2+σ-f2+σe2 for GFBLUP, and hGBLUP2=σg2σg2+σe2 for GBLUP.

Validation of genomic prediction

Genomic breeding values (GEBV) were predicted using both GFBLUP and GBLUP models. In the GFBLUP and GBLUP models, GEBV is g^total=g^f+g^-f and g^total=g^, respectively. Accuracy of predicted genomic breeding values (r) is calculated as the correlation between GEBV and DRP in the validation population. The bias of the genomic predictions with both GFBLUP and GBLUP was evaluated by the regression of DRP on the GEBV, i.e. bias=covDRP,GEBV/σGEBV2.

Single-marker GWAS

Single-marker GWAS analyses for four traits were only conducted in the HOL training population, followed by SNP set test analyses for testing the associations between genomic features and traits. Single-marker GWAS was performed using a two-step variance component-based method, to account for population stratification, as implemented in EMMAX [44]. In the first step, the polygenic and residual variances were estimated using the following model:

y = 1μae

where y is a vector of phenotypes; 1 is a vector of 1s; μ is the overall mean; a is a vector of breeding values, where aN0,Gσa2, and G is the genome relationship matrix estimated using EMMAX based on HD SNP genotypes, but excluding the SNPs on the chromosome that harbours the SNP the effect of which is being estimated; and e is the vector of residuals, where eN0,Iσe2 and I is an identity matrix. In the second step, the individual effects of SNPs were obtained using a linear regression model:

y = 1μxbη

where y, 1 and μ are as defined above; x is a vector of imputed genotype dosages (ranging from 0 to 2), b is the vector of allele substitution effects (b), and η is a vector of random residual deviates with (co)variance structure Gσa2+Iσe2.

SNP set test

Summary statistic for a genomic feature

The summary statistic of a genomic feature was calculated as the sum of the test statistics (i.e. t2) of all SNPs within DEG (i.e. open reading frame) that belonged to the genomic feature:

Tsum=i=1mftm2,

where mf is the number of SNPs located in a genomic feature, and tm2 is the square of the t-statistics for each SNP in the genomic feature. The t-statistics was calculated as the estimate of the SNP effect (i.e. b) from single-marker GWAS divided by its standard error. This summary statistic is more powerful compared to count-based summary statistics, particularly in situations where genomic features harbor many SNPs each having a small to moderate effect [9, 45].

Testing for association between a genomic feature and a trait

Under the null hypothesis, all SNPs in a genome feature have the same joint effect as those in the randomly selected genomic features. To ensure a null hypothesis is competitive to the alternative hypothesis, the random genomic features must contain the same number of SNPs as the genomic feature being analysed, and the linkage disequilibrium (LD) structure among SNPs should be retained. An empirical distribution of the summary statistics of a genomic feature was therefore obtained by using the following cyclical permutation procedure as described previously [9, 46]. Briefly, the test statistics of SNPs (i.e. t2) were first ordered based on the chromosome position of the SNPs. A test statistic was randomly selected from this vector. All test statistics were then shifted to new positions, where the selected SNP became the first one, and the other SNPs shifted to new positions, but retained their original order. This uncouples any associations between SNPs and the genomic feature, while retaining the LD structure among SNPs. A new summary statistic was then calculated according to the original position of the genomic feature. The permutation was repeated 1000 times for each genomic feature, and an empirical P value was then calculated based on one-tailed tests of the proportion of randomly sampled summary statistics that were larger than that observed.

Biological function enrichment analysis

In order to investigate the biological function of a genomic feature, functional enrichment analysis of DEG in the particular genomic feature was conducted using a web-based tool, KOBAS2.0 (http://kobas.cbi.pku.edu.cn/home.do) [47], where a hypergeometric gene set enrichment test, based on a gene ontology (GO) database, was applied. The FDR method [48] was used for adjusting multiple tests.

Results

The results for RNA-Seq analyses at different time-point comparisons (i.e. 3 vs. −22 h, 6 vs. −22 h, 9 vs. −22 h, 12 vs. −22 h and 48 vs. −22 h) are summarized in Table S2 (see Additional file 2: Table S2). The −log10(P) values of imputed sequence-level SNPs from single-marker GWAS for mastitis and milk production traits on the HOL training population are shown in the Manhattan plots of Figure S1 (see Additional file 3: Figure S1). The GFBLUP and GBLUP models were compared for all four traits in within-breed (i.e. HOL and JER) genomic prediction, followed by across-breed prediction (i.e. HOL as the training population and JER as the validation population). The degree of enrichment (i.e. −log10(P values)) of genomic features based on the SNP set test in the HOL training population was compared with the changes in prediction accuracy of GFBLUP within- and across-breed predictions, respectively.

GBLUP, GFBLUP and SNP set test analyses for Holstein population

Genomic parameters

As shown in Fig. 1a, 128, 106, 99, and 90 of the 145 genomic features explained larger proportions of the total genomic variance (Hf2) compared to their SNP-proportion over the whole genome for mastitis, protein, milk and fat yield, respectively. Detailed information is summarized in Tables S3, S4, S5 and S6 (see Additional file 4: Tables S3, S4, S5 and S6). These results demonstrated that the genomic variance of the traits studied is not uniformly distributed along the genome, but appears to be enriched in a subset of hepatic transcriptomic regions that are responsive to IMI. Therefore, the assumption of the GBLUP approach that a priori all markers contribute equally to trait variability does not hold good.

Fig. 1
Proportion of genomic variance explained by the genomic features. Each point represents one of the 145 genomic features. a is for Holstein; b is for Jersey; the x axis represents the proportion of SNPs over the whole genome that are located in genomic ...

Prediction accuracy

Prediction accuracy of GBLUP was equal to 0.504 (bias = 0.864) for mastitis, 0.602 (bias = 0.775) for protein yield, 0.635 (bias = 0.862) for milk yield, and 0.607 (bias = 0.808) for fat yield. Compared to the GBLUP model, 27, 44, 17 and 13 of the 145 genomic features resulted in higher prediction accuracies with GFBLUP (Δr ≥ 0.01) for mastitis, protein, milk and fat yield, respectively (see Additional file 4: Tables S3, S4, S5 and S6). Among these, we found 8 (9) up- (down-) regulated genomic features for mastitis, 26 (4) for protein yield, 2 (9) for milk yield, and 4 (9) for fat yield (Fig. 2). These results indicate that down-regulated genes could be more often associated with milk and fat yield than up-regulated genes during IMI. The regression coefficient of DRP on GEBV (bias) for all GFBLUP analyses ranged from 0.862 to 0.873 for mastitis, from 0.772 to 0.783 for protein yield, from 0.857 to 0.866 for milk yield, and from 0.778 to 0.821 for fat yield (see Additional file 4: Tables S3, S4, S5 and S6). The absolute value of (1-bias) tended to be negatively correlated with the change in genomic prediction accuracy with GFBLUP across four traits (see Additional file 5: Figure S2), which indicates that more predictive genomic features lead to less biased predictions. The top five predictive genomic features for each of the four traits are presented in Table 2. The average increase in prediction accuracy with the best-performing genomic feature across the four traits was 0.018, which corresponds to an increase of 3.2% relative to GBLUP.

Fig. 2
Number of up- (down-) regulated genomic features that result in higher prediction accuracy (Δr > 0.01) with GFBLUP in Holstein population. Up represents up-regulated genomic features; down represents down-regulated genomic ...
Table 2
Top five predictive genomic features for mastitis, protein, milk and fat yield in Holstein cattle

Comparisons between degrees of enrichment based on the SNP set test and changes in prediction accuracy of GFBLUP

The results of SNP set tests for all 145 genomic features across four traits in the HOL training population are summarized in Tables S3, S4, S5 and S6 (see Additional file 4: Tables S3, S4, S5 and S6). The changes in prediction accuracy of GFBLUP (Δr) were significantly (P < 0.05) positively correlated with –log10(P) of genomic features based on the SNP set test across all four traits (Fig. 3). Correlations of 0.69 (P < 2.2 × 10−16), 0.46 (P = 4.4 × 10−9), 0.46 (P = 4.4 × 10−9) and 0.44 (P = 3.6 × 10−8) were found between changes in accuracy and −log10(P value) for mastitis, protein yield, milk yield, and fat yield, respectively. These results demonstrated that the SNP set test could be used as a computationally simple way to develop more predictive GFBLUP models.

Fig. 3
Comparisons between degree of enrichment from the SNP set test in the Holstein (HOL) training (reference) population and changes in prediction accuracy with GFBLUP in the HOL validation population. Each point represents one of the 145 genomic features ...

GBLUP and GFBLUP analyses for the Jersey population

Genomic parameters

As in the analyses for the HOL population (Fig. 1b), we observed that 125, 115, 99, and 83 of the 145 genomic features for the JER population explained a larger proportion of the total genomic variance relative to their SNP-proportion over the whole genome for mastitis, protein yield, milk yield, and fat yield, respectively. Detailed information is in Tables S7, S8, S9 and S10 (see Additional file 6: Tables S7, S8, S9 and S10). It should be noted that all genomic features were defined based on gene expression data that were obtained in HOL cattle. These results imply that a subset of hepatic transcriptomic regions responsive to IMI found for HOL were also enriched in genomic variants for mastitis, protein, milk and fat yield in JER.

Prediction accuracy

Prediction accuracy of the GBLUP model was equal to 0.549 (bias = 0.916) for mastitis, 0.530 (bias = 0.760) for protein yield, 0.597 (bias = 0.796) for milk yield, and 0.433 (bias = 0.669) for fat yield. Compared to the GBLUP model, 21, 14 and 2 genomic features resulted in higher prediction accuracy (Δr ≥ 0.01) with GFBLUP for mastitis, protein, and milk yield, respectively (see Additional file 6: Tables S7, S8, S9 and S10), among which 7, 13 and 0 were in common with those found for HOL, respectively. No genomic features resulted in an increase >0.005 in prediction accuracy for fat yield in JER. The regression coefficient of DRP on GEBV (i.e. bias) for all the GFBLUP analyses ranged from 0.891 to 0.930 for mastitis, from 0.727 to 0.807 for protein yield, from 0.760 to 0.809 for milk yield, and from 0.599 to 0.677 for fat yield. As observed in HOL, the absolute value of (1-bias) was negatively correlated with the change in prediction accuracy for all four traits in JER (see Additional file 7: Figure S3). The top five predictive genomic features for each of the four traits are summarized in Table 3. The average increase in prediction accuracy (Δr) with the best-performing genomic feature across the four traits was 0.020, which corresponds to a 3.9% increase compared to GBLUP. These results indicate that the use of gene expression data obtained from one breed may improve marginally the genomic prediction accuracy in other breeds. It should be noted that, for JER, the increase in prediction accuracy with GFBLUP for milk and fat yield was very small (Table 3).

Table 3
Top five predictive genomic features for mastitis, protein, milk and fat yield in Jersey cattle

Comparisons between degree of enrichment from the SNP set test and changes in prediction accuracy of GFBLUP

The changes in prediction accuracy with GFBLUP on the JER validation population were also significantly positively correlated with −log10(P) based on the SNP set test on the HOL training population for mastitis and protein yield (Fig. 4). Correlations of 0.59 (P = 3.0 × 10−15), 0.52 (P = 3.1 × 10−11), 0.19 (P = 0.02) and 0.06 (P = 0.5) were found between changes in accuracy and −log10(P) for mastitis, protein yield, milk yield, and fat yield, respectively.

Fig. 4
Comparisons between degree of enrichment from the SNP set test in the Holstein (HOL) training (reference) population and changes in prediction accuracy with GFBLUP in the Jersey (JER) validation population. Each point represents one of the 145 genomic ...

GBLUP and GFBLUP for across-breed genomic prediction

When the complete HOL population was considered as training population to predict the genomic values of individuals in the JER population, prediction accuracy of GBLUP was very low, i.e. prediction accuracies were equal to −0.058 (bias = −0.343) for mastitis, 0.098 (bias = 0.622) for protein yield, 0.160 (bias = 0.762) for milk yield, and 0.070 (bias = 0.482) for fat yield. Compared to the GBLUP model, 60, 68, 71 and 44 of the 145 genomic features resulted in higher prediction accuracy with GFBLUP (Δr ≥ 0.01) for mastitis, protein, milk and fat yield, respectively (see Additional file 8: Tables S11, S12, S13 and S14). The regression coefficient (i.e. bias) of DRP on GEBV for all GFBLUP analyses ranged from −0.463 to 0.277 for mastitis, from 0.151 to 1.265 for protein yield, from 0.413 to 0.826 for milk yield, and from 0.002 to 0.577 for fat yield. It should be noted that more predictive genomic features lead to less biased predictions across the four traits (see Additional file 9: Figure S4). In addition, for mastitis, protein and milk yield, the changes in accuracy with GFBLUP in across-breed prediction were significantly correlated with the −log10(P) of SNP set test in the HOL training population (Fig. 5). The top five predictive genomic features for each of the four traits are summarized in Table 4. The absolute average increase in prediction accuracy (Δr) with the best-performing genomic feature across four traits was 0.111, which corresponds to a 164.4% increase relative to GBLUP. Compared to within-breed prediction, the relative improvement in genomic prediction accuracy seems to be clearer in across-breed prediction.

Fig. 5
Comparisons between degree of enrichment from the SNP set test in the Holstein (HOL) training (reference) population and changes in prediction accuracy with GFBLUP in the across-breed prediction. Each point represents one of the 145 genomic features
Table 4
Top five predictive genomic features for mastitis, protein, milk and fat yield in across-breed prediction

Discovery of gene sets associated with protein yield

Genomic features can be ranked based on the predictive ability of GFBLUP. Therefore, our GFBLUP can also be used to map gene sets that are associated with complex traits. For instance, a highly up-regulated genomic feature with 34 DEG (FDR < 10−6; log2(fold-change) > 2) that were detected in the 48 vs. −22 h comparison resulted in an increase of 0.204, 0.020 and 0.041 in prediction accuracy for protein yield among across-breed, and within HOL and JER predictions, respectively (see Additional file 10: Table S15). These 34 DEG, which include <0.01% of the total number of SNPs, explained 1.84 and 4.59% of the genomic variance for protein yield in HOL and JER, respectively. In addition, they explained 0.44 and 0.50% of the genomic variance for mastitis in HOL and JER, respectively, but did not improve genomic predictions for mastitis. Detailed information of GFBLUP analyses for these 34 DEG across three prediction scenarios is in Table 5. The P values based on the SNP set test were 0.021 and 0.18 for protein yield and mastitis, respectively, on the HOL training population. The functional enrichment analysis of these 34 DEG revealed that they were significantly (FDR < 0.05) enriched in innate immune response and negative regulation of endopeptidase activity and protein metabolism (Fig. 6).

Table 5
GFBLUP analyses of 34 genes detected in the comparison 48 h vs. −22 h (FDR < 10−6; log2(fold-change) > 2) for mastitis, protein, milk and fat yield
Fig. 6
Significantly enriched (FDR < 0.05) biological processes (BP) for the 34 genes detected in the comparison 48 versus −22 h (FDR < 10−6; log2(fold-change) >2). The significance ...

Discussion

In the current study, we demonstrated that a subset of the hepatic transcriptomic regions responsive to IMI was enriched in genomic variants associated with mastitis and milk production traits. When using these regions as genomic features, the genomic prediction accuracy with GFBLUP was improved marginally compared to GBLUP. In theory, both the GFBLUP model and SNP set test can easily be extended to incorporate other types of biological information as genomic features, such as sequence annotation, biological pathways and eQTL.

Dissection of the genetic architecture and improvement of prediction accuracy for mastitis and milk production traits in dairy cattle

It has been suggested that milk production and disease resistance traits are controlled by several hundred up to several thousand loci in cattle, most of which have a very small effect [4, 49, 50]. Multiple studies, using different strategies, have been conducted to investigate the genetic architecture that underlies such complex phenotypes, and to improve genomic prediction accuracy within and across breeds [6, 17, 49, 51, 52].

Genetic architecture and biological interpretation

The approaches that partition genomic variance based on adjacent genomic regions (e.g. 50-SNP genomic segments) or single chromosomes may not provide enough biological insights into the genetic architecture of a trait [6, 51, 53]. Our results provide evidence that results from gene expression experiments can give additional information about the biological and genetic basis of complex traits. In the current study, we used RNA-Seq data from an IMI experiment as an example to study the genetic and biological basis of mastitis and milk production traits. We found that a subset of hepatic transcriptomic regions responsive to IMI is enriched in genomic variants associated with these traits. We also found that down-regulated genes are more often associated with milk and fat yield, which together with the fact that the liver is a crucial organ for host immune responses and metabolism, including lipogenesis, gluconeogenesis, and cholesterol metabolism [54, 55], implies that the immune responses in the liver during mastitis impair milk production. This is in agreement with a recent study that demonstrated that immune relevant pathways (e.g. leukocyte trans endothelial migration and chemokine signalling pathways) are strongly associated with milk and fat yield in HOL [17].

Within-breed prediction

In populations with a high degree of linkage disequilibrium (LD), such as highly selected dairy cattle breeds, the genomic relationship based on genome-wide markers provides accurate information about the genomic variation of the traits [56], although it does not use any prior biological information. In addition, the LD structure makes it more difficult to partition genomic variance based on genomic features. Therefore, the increase in prediction accuracy with GFBLUP is small compared to GBLUP, i.e. we observed average increases of 0.018 and 0.022 across four traits within HOL and JER, respectively. This is consistent with a recent study [52] that applied a Bayesian genomic feature model (i.e. BayesRC) to milk production traits. Incorporating 790 candidate genes associated with milk production traits as a genomic feature, they found that the increases in within-breed prediction accuracy with BayesRC were quite small (<0.01) compared to BayesR, which ignores any prior biological information [52].

Across-breed prediction

Across-breed genomic prediction accuracies for milk production traits were close to zero, when HOL was used as training population to predict genomic values for JER using the GBLUP approach. This is in agreement with observations in [50, 56]. When validation and training populations are distantly related (i.e. the LD structure becomes weak), genomic feature modelling approaches such as GFBLUP and BayesRC are expected to perform better than models that ignore prior biological information such as GBLUP and BayesR, provided that the genomic feature is enriched in the genomic variants of the traits across breeds [8, 52]. Therefore, shifting the focus from the complete set of genomic markers to those that are more likely to have functional effects might contribute to improve across-breed genomic predictions [7], as observed in our study. However, breed differences in the segregation of quantitative trait loci (QTL), minor allele frequencies and breed-specific SNP effects could add to the complexity in across-breed prediction.

GFBLUP and alternatives

Factors that influence the performance of GFBLUP

The assumption made in the GBLUP model (i.e. the genomic variance is evenly distributed along the whole genome) does not match the real genetic architecture that underlies the traits. It puts equal weights to the elements in the genomic relationship, whereas the GFBLUP allows putting different weights to the individual genomic relationships in the prediction equation according to the estimated genomic parameters [8]. Prediction accuracy of GFBLUP is influenced both by the genomic variance explained by the genomic features and by the number of non-causal SNPs in the feature [8, 9]. The GFBLUP model performs better as the genomic feature contains more causal variants (i.e. explaining more genomic variance) and less non-causal markers [8, 9]. However, if the estimated genomic parameters deviate from the true values, it will lead to reduced prediction accuracy, as shown in the current study (Figs. 3, ,4,4, ,5),5), because too much weight is put on the “wrong” genomic relationships in the prediction equations. Our GFBLUP has two components for genomic effects (i.e. f and - f ), but in theory it is possible to include multiple genomic feature effects [57, 58], which might improve genomic predictions more compared to the current GFBLUP. However, when the correlations among multiple genomic relationship matrices are high, the variance components are not reliably estimated and thus there is no improvement in prediction accuracy [8, 57]. Therefore, further work is needed to investigate the performance of the GFBLUP model with multiple genomic features, in particular in livestock populations with large LD structures.

Bayesian mixture model and Bayesian GF mixture model

Bayesian mixture models, such as BayesR [50], which ignore prior genomic feature information, are considered to be relevant alternative methods. Both GFBLUP and Bayesian mixture models allow assigning markers to different distributions. GFBLUP assigns a marker set (i.e. genomic feature) to a certain distribution [i.e. fN0,Gfσf2 or -fN0,G-fσ-f2] using prior biological knowledge, whereas Bayesian mixture models attempt to assign markers to predefined distributions based on the data themselves. Previous studies demonstrated that an externally informed genomic feature is necessary for a successful partitioning of genomic variance, while the data themselves may not necessarily suggest which marker should have the greatest weight [8, 50]. The external biological information can also be incorporated into Bayesian mixture models, such as BayesRC [52]. All genomic feature models including GFBLUP and BayesRC are computationally intensive, and they do not necessarily perform better than standard models (i.e. GBLUP and BayesR) when genomic features are less enriched in causal variants [8, 59].

SNP set test

The SNP set test based on single-marker test statistics derived from GWAS is a computationally fast way to evaluate a large number of genomic features [60]. The results of the SNP set test could be used to develop more predictive GFBLUP and similar models. The current SNP set test method assesses the association between a genomic feature and a trait based on the sum of t 2 of SNPs within the genomic feature. Another commonly used approach for the SNP set test is based on counting associations exceeding a pre-defined significance threshold within the genomic feature [6163]. One important limitation of this count-based approach is the dichotomization of association signals into significant and non-significant sets, based on a pre-specified significance level, which ignores information regarding the strength of association. Since the genomic variance of mastitis and milk production traits is typically governed by very many markers, each with a small effect [4, 49, 50], the current SNP set test is more likely to match the genetic basis of complex phenotypes, and is more powerful than the count-based approach [9, 45, 46].

Appropriate genomic features facilitate improved biological interpretation

In order to test different biological hypotheses, many genomic features can be constructed using different sources of prior information, such as prior QTL regions, chromosomes, sequence, biological pathways, and other types of external evidence. The gain in biological knowledge of complex traits relies highly on the genomic feature classification strategies. Since associated genomic markers are not evenly, or necessarily physically, clustered along the genome [2, 51], partitioning genomic variance based on adjacent genomic regions (e.g. haplotypes and chromosomes) is not an ideal way to facilitate the interpretation of biological mechanisms underlying the traits. Biological interpretation may be better served by the use of pathways and gene ontologies as genomic features; however, the quantity and quality of the genes that are functionally annotated in current pathway databases are limited [15], particularly for livestock and plant genomes. Here, we used information from gene expression data to define genomic features, providing novel insights into the genetic and biological basis of mastitis and milk production traits and improving genomic prediction accuracy with GFBLUP.

Since mastitis can be caused by various pathogens, the current RNA-Seq data that originate only from E. coli mastitis may be limited to detect all the genes that are functionally relevant with mastitis. Thus, more RNA-Seq data from infections with other types of pathogens could help the detection of genomic features that are associated with mastitis and milk production. In addition, since gene expression patterns depend highly on time, cell types, and tissues, some trait-associated genes might not show differential expression in certain cell types and tissues at a certain physiological stage. Therefore, incorporating more molecular biological information from more tissues (e.g. mammary gland, blood and adipose tissue) and more physiological stages could be important to define the appropriate genomic features that are highly enriched in causal variants.

Conclusions

Compared to GBLUP, GFBLUP models increased the accuracy of genomic prediction for mastitis and milk production traits in dairy cattle by incorporating biological information from gene expression data, and thus provide novel biological insights into the genetic basis of such complex traits. Compared to within-breed prediction, the increase in prediction accuracy seems to be more apparent in across-breed prediction. In addition, the SNP set test can be used as a computationally fast way to develop more predictive GFBLUP or similar models. The current genomic feature modelling approaches provide a general framework for incorporating biological knowledge from independent functional genomics studies to study the genetic architecture and to improve genomic prediction for complex traits. Approaches such as GFBLUP and SNP set test will be increasingly useful as the biological knowledge of functional genomic regions keep accumulating for a range of traits and species.

Additional files

Authors’ contributions

LF conceived the study, designed, performed and evaluated the experiments, analyzed the data and drafted the manuscript. PM contributed to the study design and analyzed the data. GSa, GSu, YY, SZ and ML contributed to the study design and drafted the manuscript. PS conceived and designed the study, and drafted the manuscript. All authors read and approved the final manuscript.

Competing interests

The authors declare that they have no competing interests.

Availability of data and materials

All genomic annotation data defining gene regions is available for download (ftp://ftp.ensembl.org/pub/release-84/gtf/bos_taurus). The whole-genome sequencing data from the 1000 Bull Genomes Project are publicly available, as variations in dbSNP (http://www.ncbi.nlm.nih.gov/projects/SNP/) and sequence data at NBCI with SRA no. SRP039339 (http://www.ncbi.nlm.nih.gov/bioproject/PRJNA238491). Whole-genome sequences from Aarhus University and individual SNP genotype data are available only upon agreement with the commercial breeding organization (http://www.vikinggenetics.com/) and should be requested directly from the authors or the breeding organization.

Ethics approval and consent to participate

All experimental procedures involving animals were approved by the Danish Animal Experiments Inspectorate and complied with Danish Ministry of Justice Laws concerning animal experimentation and care of experimental animals. Animal experiments were conducted in strict accordance with regulations and guidelines established by these committees. Inspection was carried out by members of these committees during the animal infection experiments.

Funding

This study was funded by the Danish Strategic Research Council (GenSAP: Centre for Genomic Selection in Animals and Plants, Contract 12-132452).

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Footnotes

Electronic supplementary material

The online version of this article (doi:10.1186/s12711-017-0319-0) contains supplementary material, which is available to authorized users.

Contributor Information

Lingzhao Fang, kd.ua.gbm@gnaf.oahzgnil.

Goutam Sahana, kd.ua.gbm@anahas.matuog.

Peipei Ma, nc.ude.utjs@am.iepiep.

Guosheng Su, kd.ua.gbm@us.gnehsoug.

Ying Yu, nc.ude.uac@gniyuy.

Shengli Zhang, nc.ude.uac@uaclsgnahz.

Mogens Sandø Lund, kd.ua.gbm@dnul.snegom.

Peter Sørensen, kd.ua.gbm@osp.

References

1. Wang WY, Barratt BJ, Clayton DG, Todd JA. Genome-wide association studies: theoretical and practical concerns. Nat Rev Genet. 2005;6:109–118. doi: 10.1038/nrg1522. [PubMed] [Cross Ref]
2. Lango Allen H, Estrada K, Lettre G, Berndt SI, Weedon MN, Rivadeneira F, et al. Hundreds of variants clustered in genomic loci and biological pathways affect human height. Nature. 2010;467:832–838. doi: 10.1038/nature09410. [PMC free article] [PubMed] [Cross Ref]
3. Yang J, Manolio TA, Pasquale LR, Boerwinkle E, Caporaso N, Cunningham JM, et al. Genome partitioning of genetic variation for complex traits using common SNPs. Nat Genet. 2011;43:519–525. doi: 10.1038/ng.823. [PMC free article] [PubMed] [Cross Ref]
4. Kemper KE, Goddard ME. Understanding and predicting complex traits: knowledge from cattle. Hum Mol Genet. 2012;21:R45–R51. doi: 10.1093/hmg/dds332. [PubMed] [Cross Ref]
5. Visscher PM, Brown MA, McCarthy MI, Yang J. Five years of GWAS discovery. Am J Hum Genet. 2012;90:7–24. doi: 10.1016/j.ajhg.2011.11.029. [PubMed] [Cross Ref]
6. Hayes BJ, Pryce J, Chamberlain AJ, Bowman PJ, Goddard ME. Genetic architecture of complex traits and accuracy of genomic prediction: coat colour, milk-fat percentage, and type in Holstein cattle as contrasting model traits. PLoS Genet. 2010;6:e1001139. doi: 10.1371/journal.pgen.1001139. [PMC free article] [PubMed] [Cross Ref]
7. Snelling WM, Cushman RA, Keele JW, Maltecca C, Thomas MG, Fortes MR, et al. Breeding and genetics symposium: networks and pathways to guide genomic selection. J Anim Sci. 2013;91:537–552. doi: 10.2527/jas.2012-5784. [PubMed] [Cross Ref]
8. Edwards SM, Sørensen IF, Sarup P, Mackay TF, Sørensen P. Genomic prediction for quantitative traits is improved by mapping variants to gene ontology categories in Drosophila melanogaster. Genetics. 2016;203:1871–1883. doi: 10.1534/genetics.116.187161. [PubMed] [Cross Ref]
9. Sarup P, Jensen J, Ostersen T, Henryon M, Sørensen P. Increased prediction accuracy using a genomic feature model including prior information on quantitative trait locus regions in purebred Danish Duroc pigs. BMC Genet. 2016;17:11. doi: 10.1186/s12863-015-0322-9. [PMC free article] [PubMed] [Cross Ref]
10. Meuwissen THE, Hayes BJ, Goddard ME. Prediction of total genetic value using genome-wide dense marker maps. Genetics. 2001;157:1819–1829. [PubMed]
11. Maurano MT, Humbert R, Rynes E, Thurman RE, Haugen E, Wang H, et al. Systematic localization of common disease-associated variation in regulatory DNA. Science. 2012;337:1190–1195. doi: 10.1126/science.1222794. [PMC free article] [PubMed] [Cross Ref]
12. O’Roak BJ, Vives L, Girirajan S, Karakoc E, Krumm N, Coe BP, et al. Sporadic autism exomes reveal a highly interconnected protein network of de novo mutations. Nature. 2012;485:246–250. doi: 10.1038/nature10989. [PMC free article] [PubMed] [Cross Ref]
13. Wang K, Li M, Hakonarson H. Analysing biological pathways in genome-wide association studies. Nat Rev Genet. 2010;11:843–854. doi: 10.1038/nrg2884. [PubMed] [Cross Ref]
14. Khatri P, Sirota M, Butte AJ. Ten years of pathway analysis: current approaches and outstanding challenges. PLoS Comput Biol. 2012;8:e1002375. doi: 10.1371/journal.pcbi.1002375. [PMC free article] [PubMed] [Cross Ref]
15. Ramanan VK, Shen L, Moore JH, Saykin AJ. Pathway analysis of genomic data: concepts, methods, and prospects for future development. Trends Genet. 2012;28:323–332. doi: 10.1016/j.tig.2012.03.004. [PMC free article] [PubMed] [Cross Ref]
16. Sedeño-Cortés AE, Pavlidis P. Pitfalls in the application of gene-set analysis to genetics studies. Trends Genet. 2014;30:513–514. doi: 10.1016/j.tig.2014.10.001. [PMC free article] [PubMed] [Cross Ref]
17. Edwards SM, Thomsen B, Madsen P, Sorensen P. Partitioning of genomic variance reveals biological pathways associated with udder health and milk production traits in dairy cattle. Genet Sel Evol. 2015;47:60. doi: 10.1186/s12711-015-0132-6. [PMC free article] [PubMed] [Cross Ref]
18. Mitterhuemer S, Petzl W, Krebs S, Mehne D, Klanner A, Wolf E, et al. Escherichia coli infection induces distinct local and systemic transcriptome responses in the mammary gland. BMC Genomics. 2010;11:138. doi: 10.1186/1471-2164-11-138. [PMC free article] [PubMed] [Cross Ref]
19. Buitenhuis B, Rontved CM, Edwards SM, Ingvartsen KL, Sorensen P. In depth analysis of genes and pathways of the mammary gland involved in the pathogenesis of bovine Escherichia coli-mastitis. BMC Genomics. 2011;12:130. doi: 10.1186/1471-2164-12-130. [PMC free article] [PubMed] [Cross Ref]
20. Costa V, Aprile M, Esposito R, Ciccodicola A. RNA-Seq and human complex diseases: recent accomplishments and future perspectives. Eur J Hum Genet. 2013;21:134–142. doi: 10.1038/ejhg.2012.129. [PMC free article] [PubMed] [Cross Ref]
21. Aitken SL, Corl CM, Sordillo LM. Immunopathology of mastitis: insights into disease recognition and resolution. J Mammary Gland Biol Neoplasia. 2011;16:291–304. doi: 10.1007/s10911-011-9230-4. [PubMed] [Cross Ref]
22. Wellnitz O, Bruckmaier RM. The innate immune response of the bovine mammary gland to bacterial infection. Vet J. 2012;192:148–152. doi: 10.1016/j.tvjl.2011.09.013. [PubMed] [Cross Ref]
23. Jiang L, Sørensen P, Røntved C, Vels L, Ingvartsen KL. Gene expression profiling of liver from dairy cows treated intra-mammary with lipopolysaccharide. BMC Genomics. 2008;9:443. doi: 10.1186/1471-2164-9-443. [PMC free article] [PubMed] [Cross Ref]
24. Pisoni G, Moroni P, Genini S, Stella A, Boettcher PJ, Cremonesi P, et al. Differentially expressed genes associated with Staphylococcus aureus mastitis in dairy goats. Vet Immunol Immunopathol. 2010;135:208–217. doi: 10.1016/j.vetimm.2009.11.016. [PubMed] [Cross Ref]
25. Loor JJ, Moyes KM, Bionaz M. Functional adaptations of the transcriptome to mastitis-causing pathogens: the mammary gland and beyond. J Mammary Gland Biol Neoplasia. 2011;16:305–322. doi: 10.1007/s10911-011-9232-2. [PubMed] [Cross Ref]
26. Jiang L, Sorensen P, Thomsen B, Edwards SM, Skarman A, Rontved CM, et al. Gene prioritization for livestock diseases by data integration. Physiol Genomics. 2012;44:305–317. doi: 10.1152/physiolgenomics.00047.2011. [PubMed] [Cross Ref]
27. Hotamisligil GS. Inflammation and metabolic disorders. Nature. 2006;444:860–867. doi: 10.1038/nature05485. [PubMed] [Cross Ref]
28. Huang W, Richards S, Carbone MA, Zhu D, Anholt RR, Ayroles JF, et al. Epistasis dominates the genetic architecture of Drosophila quantitative traits. Proc Nat Acad Sci USA. 2012;109:15553–15559. doi: 10.1073/pnas.1213423109. [PubMed] [Cross Ref]
29. Vels L, Rontved CM, Bjerring M, Ingvartsen KL. Cytokine and acute phase protein gene expression in repeated liver biopsies of dairy cows with a lipopolysaccharide-induced mastitis. J Dairy Sci. 2009;92:922–934. doi: 10.3168/jds.2008-1209. [PubMed] [Cross Ref]
30. Andersen JB, Mashek DG, Larsen T, Nielsen MO, Ingvartsen KL. Effects of hyperinsulinaemia under euglycaemic condition on liver fat metabolism in dairy cows in early and mid-lactation. J Vet Med Physiol Pathol Clin Med. 2002;49:65–71. doi: 10.1046/j.1439-0442.2002.jv417.x. [PubMed] [Cross Ref]
31. Moyes KM, Sørensen P, Bionaz M. The impact of intramammary Escherichia coli challenge on liver and mammary transcriptome and cross-talk in dairy cows during early lactation using RNAseq. PLoS One. 2016;11:e0157480. doi: 10.1371/journal.pone.0157480. [PMC free article] [PubMed] [Cross Ref]
32. Liao Y, Smyth GK, Shi W. The Subread aligner: fast, accurate and scalable read mapping by seed-and-vote. Nucleic Acids Res. 2013;41:e108. doi: 10.1093/nar/gkt214. [PMC free article] [PubMed] [Cross Ref]
33. Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26:139–140. doi: 10.1093/bioinformatics/btp616. [PMC free article] [PubMed] [Cross Ref]
34. Gao H, Christensen OF, Madsen P, Nielsen US, Zhang Y, Lund MS, et al. Comparison on genomic predictions using three GBLUP methods and two single-step blending methods in the Nordic Holstein population. Genet Sel Evol. 2012;44:8. doi: 10.1186/1297-9686-44-8. [PMC free article] [PubMed] [Cross Ref]
35. Thomasen JR, Guldbrandtsen B, Su G, Brøndum RF, Lund MS. Reliabilities of genomic estimated breeding values in Danish Jersey. Animal. 2012;6:789–796. doi: 10.1017/S1751731111002035. [PubMed] [Cross Ref]
36. Brondum RF, Guldbrandtsen B, Sahana G, Lund MS, Su G. Strategies for imputation to whole genome sequence using a single or multi-breed reference population in cattle. BMC Genomics. 2014;15:728. doi: 10.1186/1471-2164-15-728. [PMC free article] [PubMed] [Cross Ref]
37. Wu X, Guldbrandtsen B, Lund MS, Sahana G. Association analysis for feet and legs disorders with whole-genome sequence variants in 3 dairy cattle breeds. J Dairy Sci. 2016;99:7221–7231. doi: 10.3168/jds.2015-10705. [PubMed] [Cross Ref]
38. Daetwyler HD, Capitan A, Pausch H, Stothard P, Van Binsbergen R, Brøndum RF, et al. Whole-genome sequencing of 234 bulls facilitates mapping of monogenic and complex traits in cattle. Nat Genet. 2014;46:858–865. doi: 10.1038/ng.3034. [PubMed] [Cross Ref]
39. Höglund JK, Sahana G, Brøndum RF, Guldbrandtsen B, Buitenhuis B, Lund MS. Fine mapping QTL for female fertility on BTA04 and BTA13 in dairy cattle using HD SNP and sequence data. BMC Genomics. 2014;15:790. doi: 10.1186/1471-2164-15-790. [PMC free article] [PubMed] [Cross Ref]
40. Fuchsberger C, Abecasis GR, Hinds DA. minimac2: faster genotype imputation. Bioinformatics. 2015;31:782–784. doi: 10.1093/bioinformatics/btu704. [PMC free article] [PubMed] [Cross Ref]
41. VanRaden PM. Efficient methods to compute genomic predictions. J Dairy Sci. 2008;91:4414–4423. doi: 10.3168/jds.2007-0980. [PubMed] [Cross Ref]
42. Johnson DL, Thompson R. Restricted maximum likelihood estimation of variance components for univariate animal models using sparse matrix techniques and average information. J Dairy Sci. 1995;78:449–456. doi: 10.3168/jds.S0022-0302(95)76654-1. [Cross Ref]
43. Madsen P, Jensen J, Labouriau R, Christensen OF, Sahana G. DMU-A package for analyzing multivariate mixed models in quantitative genetics and genomics. In: Proceedings of the 10th world congress of genetics applied to livestock production: 18–22 August 2014; Vancouver; 2014. https://asas.org/docs/default-source/wcgalp-posters/699_paper_9580_manuscript_758_0.pdf?sfvrsn=2.
44. Kang HM, Sul JH, Service SK, Zaitlen NA, Kong SY, Freimer NB, et al. Variance component model to account for sample structure in genome-wide association studies. Nat Genet. 2010;42:348–354. doi: 10.1038/ng.548. [PMC free article] [PubMed] [Cross Ref]
45. Newton MA, Quintana FA, Den Boon JA, Sengupta S, Ahlquist P. Random-set methods identify distinct aspects of the enrichment signal in gene-set analysis. Ann Appl Stat. 2007;1:85–106. doi: 10.1214/07-AOAS104. [Cross Ref]
46. Rohde PD, Demontis D, Cuyabano BCD, Børglum AD, Sørensen P, The GEMS Group Covariance association test (CVAT) identify genetic markers associated with schizophrenia in functionally associated biological processes. Genetics. 2016;203:1901–1913. doi: 10.1534/genetics.116.189498. [PubMed] [Cross Ref]
47. Xie C, Mao X, Huang J, Ding Y, Wu J, Dong S, et al. KOBAS 2.0: a web server for annotation and identification of enriched pathways and diseases. Nucleic Acids Res. 2011;39:W316–W322. doi: 10.1093/nar/gkr483. [PMC free article] [PubMed] [Cross Ref]
48. Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J R Statist Soc Series B Stat Methodol. 1995;57:289–300.
49. Pimentel Eda G, Erbe M, König S, Simianer H. Genome partitioning of genetic variation for milk production and composition traits in Holstein cattle. Front Genet. 2011;2:19. [PMC free article] [PubMed]
50. Erbe M, Hayes B, Matukumalli LK, Goswami S, Bowman PJ, Reich CM, et al. Improving accuracy of genomic predictions within and between dairy cattle breeds with imputed high-density single nucleotide polymorphism panels. J Dairy Sci. 2012;95:4114–4129. doi: 10.3168/jds.2011-5019. [PubMed] [Cross Ref]
51. Jensen J, Su G, Madsen P. Partitioning additive genetic variance into genomic and remaining polygenic components for complex traits in dairy cattle. BMC Genet. 2012;13:44. doi: 10.1186/1471-2156-13-44. [PMC free article] [PubMed] [Cross Ref]
52. MacLeod IM, Bowman PJ, Vander Jagt CJ, Haile-Mariam M, Kemper KE, Chamberlain AJ, et al. Exploiting biological priors and sequence variants enhances QTL discovery and genomic prediction of complex traits. BMC Genomics. 2016;17:144. doi: 10.1186/s12864-016-2443-6. [PMC free article] [PubMed] [Cross Ref]
53. Tiezzi F, Parker-Gaddis KL, Cole JB, Clay JS, Maltecca C. A genome-wide association study for clinical mastitis in first parity US Holstein cows using single-step approach and genomic matrix re-weighting procedure. PLoS ONE. 2015;10:e0114919. doi: 10.1371/journal.pone.0114919. [PMC free article] [PubMed] [Cross Ref]
54. Gao B, Jeong WI, Tian Z. Liver: an organ with predominant innate immunity. Hepatology. 2008;47:729–736. doi: 10.1002/hep.22034. [PubMed] [Cross Ref]
55. Bechmann LP, Hannivoort RA, Gerken G, Hotamisligil GS, Trauner M, Canbay A. The interaction of hepatic lipid and glucose metabolism in liver diseases. J Hepatol. 2012;56:952–964. doi: 10.1016/j.jhep.2011.08.025. [PubMed] [Cross Ref]
56. Hayes BJ, Bowman PJ, Chamberlain AC, Verbyla K, Goddard ME. Accuracy of genomic breeding values in multi-breed dairy cattle populations. Genet Sel Evol. 2009;41:51. doi: 10.1186/1297-9686-41-51. [PMC free article] [PubMed] [Cross Ref]
57. Gusev A, Lee SH, Trynka G, Finucane H, Vilhjálmsson BJ, Xu H, et al. Partitioning heritability of regulatory and cell-type-specific variants across 11 common diseases. Am J Hum Genet. 2014;95:535–552. doi: 10.1016/j.ajhg.2014.10.004. [PubMed] [Cross Ref]
58. Speed D, Balding DJ. MultiBLUP: improved SNP-based prediction for complex traits. Genome Res. 2014;24:1550–1557. doi: 10.1101/gr.169375.113. [PubMed] [Cross Ref]
59. Ober U, Ayroles JF, Stone EA, Richards S, Zhu D, Gibbs RA, et al. Using whole-genome sequence data to predict quantitative trait phenotypes in Drosophila melanogaster. PLoS Genet. 2012;8:e1002685. doi: 10.1371/journal.pgen.1002685. [PMC free article] [PubMed] [Cross Ref]
60. Fridley BL, Biernacka JM. Gene set analysis of SNP data: benefits, challenges, and future directions. Eur J Hum Genet. 2011;19:837–843. doi: 10.1038/ejhg.2011.57. [PMC free article] [PubMed] [Cross Ref]
61. Holmans P, Green EK, Pahwa JS, Ferreira MA, Purcell SM, Sklar P, et al. Gene ontology analysis of GWA study data sets provides insights into the biology of bipolar disorder. Am J Hum Genet. 2009;85:13–24. doi: 10.1016/j.ajhg.2009.05.011. [PubMed] [Cross Ref]
62. Medina I, Montaner D, Bonifaci N, Pujana MA, Carbonell J, Tarraga J, et al. Gene set-based analysis of polymorphisms: finding pathways or biological processes associated to traits in genome-wide association studies. Nucleic Acids Res. 2009;37:W340–W344. doi: 10.1093/nar/gkp481. [PMC free article] [PubMed] [Cross Ref]
63. O’Dushlaine C, Kenny E, Heron EA, Segurado R, Gill M, Morris DW, et al. The SNP ratio test: pathway analysis of genome-wide association datasets. Bioinformatics. 2009;25:2762–2763. doi: 10.1093/bioinformatics/btp448. [PubMed] [Cross Ref]

Articles from Genetics, Selection, Evolution : GSE are provided here courtesy of BioMed Central