Volume 5 Supplement 3
Hierarchical likelihood opens a new way of estimating genetic values using genome-wide dense marker maps
© Shen et al; licensee BioMed Central Ltd. 2011
Published: 27 May 2011
Genome-wide dense markers have been used to detect genes and estimate relative genetic values. Among many methods, Bayesian techniques have been widely used and shown to be powerful in genome-wide breeding value estimation and association studies. However, computation is known to be intensive under the Bayesian framework, and specifying a prior distribution for each parameter is always required for Bayesian computation. We propose the use of hierarchical likelihood to solve such problems.
Using double hierarchical generalized linear models, we analyzed the simulated dataset provided by the QTLMAS 2010 workshop. Marker-specific variances estimated by double hierarchical generalized linear models identified the QTL with large effects for both the quantitative and binary traits. The QTL positions were detected with very high accuracy. For young individuals without phenotypic records, the true and estimated breeding values had Pearson correlation of 0.60 for the quantitative trait and 0.72 for the binary trait, where the quantitative trait had a more complicated genetic architecture involving imprinting and epistatic QTL.
Hierarchical likelihood enables estimation of marker-specific variances under the likelihoodist framework. Double hierarchical generalized linear models are powerful in localizing major QTL and computationally fast.
Genetic analyses in livestock studies are generally based on information from pedigrees and molecular markers. Traditionally, a kinship matrix can be calculated using the pedigree data, which can be used in a generalized linear mixed model (GLMM) to estimate breeding values. By including genetic marker information, genomic estimated breeding values (GEBV) can be obtained taking into account the information from these markers, and also quantitative trait loci (QTL) can be mapped by associating genotypes at a certain locus to the phenotype observations.
Dense marker genotypes along genome can now be affordably obtained due to new and efficient methods for typing single nucleotide polymorphism (SNP) markers. The dense SNP maps have made genome-wide association (GWA) studies popular for gene detection. Classic GWA methods , commonly applied to study genetic diseases in humans, are based on simple repeated single marker tests across the genome. To achieve more powerful mapping and better prediction, a unified model including all the SNPs in the genome is preferred. Such models have been estimated using Bayesian methods, implemented by Markov chain Monte Carlo (MCMC) techniques that are computationally demanding [2–5]. Lee and Nelder developed the double hierarchical generalized linear model (DHGLM) in the likelihoodist framework . DHGLM enables estimation of marker-specific variances using a fast iterative algorithm without specifying any prior distributions . The likelihoodist way of estimation is conducted through a likelihood function named hierarchical likelihood (h-likelihood) .
The aim of this paper is to map QTL and report GEBV for the simulated dataset provided by QTLMAS 2010 workshop. We employ a unified analysis via the h-likelihood and model the data using DHGLM. GEBV are calculated from the estimated marker effects, and QTL are mapped by the estimated marker-specific variances.
The dataset used in this paper was simulated for the QTLMAS 2010 workshop (Poznań, Poland). A pedigree consisting of 3226 individuals in 5 generations (F0 - F4) was simulated, where F0 contains 5 males and 15 females. Each female was mated once and gave birth to about 30 progeny. Two traits were simulated, where one is quantitative (QT), and the other is binary (BT). Young individuals in F4 (individuals 2327 to 3226) had no phenotypic records. The genome was assumed to be about 5 × 108 bp long, consisting of 5 chromosomes, each of which contained about 1 × 108 bp. Each individual was genotyped for 10031 biallelic SNPs in the genome.
DHGLM provides a unified analysis for both QTL mapping and genomic breeding value estimation. Similar to BayesA, the data are modeled on two levels, i.e. both the phenotypic mean and the variance are modeled with random effects. For a quantitative trait, the phenotype y (n × 1 vector) is postulated as a random effect model
y = X β + Zg + e (1)
where g ~ N(0, diag(λ)) are the SNP effects, λ = (λ1, λ2,..., λm)′ are the variances of the SNP effects, and the residuals e ~ N(0, σ2I). The fixed effects β included an intercept and the sex effect in our application to reduce the residual errors. The SNP variances λ are modeled as
log λ = 1 a + b (2)
with an intercept a and normally distributed random effects b. The genomic estimated breeding value (GEBV) for individual i is computed as . QTL can be scanned using the marker-specific variances λ. For a binary trait, the mean of y, is modeled by the same linear predictor X β + Zg through a logit link function.
For the marker-specific variances, the correlated random effects, b, follow a multivariate normal distribution with a mean of zero and a variance-covariance matrix , where m is the number of SNPs and k, l are the SNP indices. When ρ = 1, all the SNPs have a constant variance (GLMM); when ρ = 0, the SNPs are assumed to be independent (DHGLM); and for 0 <ρ < 1, the correlation between two SNPs is a monomial function of ρ, which is referred to as the smoothed DHGLM . We propose the use of smoothed DHGLMs since it reduces the noise in marker-specific variance estimates and highlights the signals of QTL. ρ, regarded as a spatial correlation parameter, was chosen to be 0.9 in this paper, which nicely shrank the SNPs with zero effect.
where is the variance of z . j (the j-th column of Z) across individuals. These variance values can be directly calculated from the data. The contribution (heritability) of a particular SNP is expressed by .
According to the extended likelihood principle, inference of the random SNP effects g should be drawn through the h-likelihood, fixed effects β through the marginal likelihood, and variance components λ, σ2 and through the adjusted profile likelihood . However, for efficient estimation, we propose to initialize variance components and iterate the following steps until convergence ,
Where and . The subscript M stands for ‘mean’.
• Update σ2 by fitting the deviance residuals using an intercept-only gamma GLM and prior weight w M = (1 – q M )/2, where are the residuals of (4), and are the diagonal elements of The subscript 1 and 2 stand for individuals (1 to n) and SNPs (n + 1 to n + m), respectively.
where , , z = log λ + (dM2 – λ)/λ is linearized λ in a gamma GLM with a log link, and L satisfies LL′ = A. The subscript D stands for ‘dispersion’.
• Update by fitting the deviance residuals using an intercept-only gamma GLM and prior weight w D = (1 – q D )/2, where ê D are the last m residuals of (5), and q D are the last m diagonal elements of .
Results and Discussion
Estimation of SNP effects
Estimated heritability of the detected QTL and suggestive QTL for QT and BT.
h2 of QT
h2 of BT
DHGLM were shown to be an efficient and reliable approach for both QTL mapping and genomic selection. Since DHGLM can be estimated by iterating interlinked GLMs, the execution time is greatly shortened comparing to the Bayesian computation. On a Macintosh laptop with a 2 GHz processor and 4 GB memory (1067 MHz), it took about 10-20 minutes, depending on starting values, to obtain our results using our implementation in R. No priors are required for parameters in DHGLM. Main QTL mapped via DHGLM showed very good accuracy though some QTL with small effects were shrunk or smoothed down. An R package iQTL has been implemented and is available on R-Forge: https://r-forge.r-project.org/R/?groupid=845.
XS, LR and ÖC initiated the study. XS analyzed the simulated common dataset of the QTLMAS 2010 workshop and drafted the paper. LR initiated the smoothed version of double hierarchical generalized linear models. XS, LR and ÖC worked on the revision together and approved the final manuscript.
List of abbreviations used
double hierarchical generalized linear model
genomic estimated breeding values
generalized linear model
generalized linear mixed model
hierarchical generalized linear model
Markov chain Monte Carlo
quantitative trait locus/loci
quantitative trait loci and marker assisted selection
restricted maximum likelihood
single nucleotide polymorphism
true breeding values
weighted least squares.
Xia Shen is funded by a Future Research Leaders grant from the Swedish Foundation for Strategic Research (SSF) to Örjan Carlborg. Lars Rönnegård is funded by the Swedish Research Council for Environment, Agricultural Sciences and Spatial Planning (FORMAS). François Besnier is acknowledged for sharing his IBD calculation program to validate our results by variance component methods.
This article has been published as part of BMC Proceedings Volume 5 Supplement 3, 2011: Proceedings of the 14th QTL-MAS Workshop. The full contents of the supplement are available online at http://www.biomedcentral.com/1753-6561/5?issue=S3.
- Cantor RM, Lange K, Sinsheimer JS: Prioritizing GWAS results: A review of statistical methods and recommendations for their application. Am. J. Hum. Genet. 2010, 86: 6-22. 10.1016/j.ajhg.2009.11.017.PubMed CentralView ArticlePubMedGoogle Scholar
- Meuwissen THE, Hayes BJ, Goddard ME: Prediction of total genetic value using genome-wide dense marker maps. Genetics. 2001, 157 (4): 1819-1829.PubMed CentralPubMedGoogle Scholar
- Xu S: Estimating polygenic effects using markers of the entire genome. Genetics. 2003, 163: 789-801.PubMed CentralPubMedGoogle Scholar
- Xu S: An empirical Bayes method for estimating epistatic effects of quantitative trait loci. Biometrics. 2007, 63: 513-521. 10.1111/j.1541-0420.2006.00711.x.View ArticlePubMedGoogle Scholar
- Yi N, Xu S: Bayesian LASSO for quantitative trait loci mapping. Genetics. 2008, 179: 1045-1055. 10.1534/genetics.107.085589.PubMed CentralView ArticlePubMedGoogle Scholar
- Lee Y, Nelder JA: Double hierarchical generalized linear models (with discussion). Applied Statistics. 2006, 55: 139-185.Google Scholar
- Lee Y, Nelder JA, Pawitan Y: Generalized linear models with random effects: unified analysis via h-likelihood. 2006, Chapman & Hall/CRCView ArticleGoogle Scholar
- Lee Y, Nelder JA: Hierarchical generalized linear models (with discussion). J. R. Statist. Soc. B. 1996, 58: 619-678.Google Scholar
- Yi N, Banerjee S: Hierarchical generalized linear models for multiple quantitative trait locus mapping. Genetics. 2009, 181 (3): 1101-1113. 10.1534/genetics.108.099556.PubMed CentralView ArticlePubMedGoogle Scholar
- Rönnegård L, Lee Y: Hierarchical generalized linear models have a great potential in genetics and animal breeding. Proc. WCGALP. 2010, Leipzig, GermanyGoogle Scholar
- Lee Y, Nelder JA, Noh M: H-likelihood: problems and solutions. Statistical Computing. 2007, 17: 49-55. 10.1007/s11222-006-9006-7.View ArticleGoogle Scholar
- Rönnegård L, Shen X, Alam M: hglm: a package for fitting hierarchical generalized linear models. The R Journal. 2010, 2 (2): 20-28.Google Scholar
- R Development Core Team: R: A Language and Environment for Statistical Computing. 2009, R Foundation for Statistical Computing, Vienna, Austria, [http://www.R-project.org]3-900051-07-0Google Scholar
This article is published under license to BioMed Central Ltd. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.