QTLRel: an R Package for Genome-wide Association Studies in which Relatedness is a Concern
BMC Genetics volume 12, Article number: 66 (2011)
Existing software for quantitative trait mapping is either not able to model polygenic variation or does not allow incorporation of more than one genetic variance component. Improperly modeling the genetic relatedness among subjects can result in excessive false positives. We have developed an R package, QTLRel, to enable more flexible modeling of genetic relatedness as well as covariates and non-genetic variance components.
We have successfully used the package to analyze many datasets, including F34 body weight data that contains 688 individuals genotyped at 3105 SNP markers and identified 11 QTL. It took 295 seconds to estimate variance components and 70 seconds to perform the genome scan on an Linux machine equipped with a 2.40GHz Intel(R) Core(TM)2 Quad CPU.
QTLRel provides a toolkit for genome-wide association studies that is capable of calculating genetic incidence matrices from pedigrees, estimating variance components, performing genome scans, incorporating interactive covariates and genetic and non-genetic variance components, as well as other functionalities such as multiple-QTL mapping and genome-wide epistasis.
Methods to search for quantitative trait loci (QTL) in common experimental designs are well established, and software to analyze these populations is widely available. One popular package, R/qtl , provides a comprehensive toolset for QTL mapping. Since it does not allow random effects, R/qtl is most suitable for mapping populations such as F2 and backcross where individuals are equally genetically related. Software that can model polygenic effects due to genetic relatedness includes TASSEL  and EMMA . Both allow covariates as fixed effects but are only capable of incorporating a random term to account for one genetic variance component. However, both additive and dominance modes of inheritance are common for many quantitative traits. Ignoring these variance components may result in excessive false positives. Moreover, researchers may also be interested in interactive covariates, epistasis and non-genetic random effects. We have developed an R package QTLRel that meets all these needs.
Consider the following statistical model
where y is a vector of phenotypes, β is a vector of covariate effects, γ is a vector of putative QTL effects, u is a vector of polygenic effects and ε is a vector of residual effects. X, Q and Z are incidence matrices. β can be fixed, random or a mix of both fixed and random. Assume that u~ N(0, G), ε~ N(0, I σ2), and u is independent of ε. G consists of five genetic variance components including additive and dominance components as well as three other components that model excess similarity due to inbreeding . The incidence matrices corresponding to these five genetic variance components can be obtained from condensed identity coefficients as defined in [, pp.133].
Condensed identity coefficients
While other programs are available for calculating condensed identity coefficients from pedigrees , we provide a function that is especially feasible for pedigrees with a large number of generations. Condensed identity coefficients can be derived from generalized kinship coefficients . Bottom-up and top-down are two computational strategies for calculating generalized kinship coefficients from a pedigree. The bottom-up approach starts from the target individuals and moves up the pedigree until reaching the founders. It requires minimal storage but the computational load increases approximately exponentially with the number of generations. The bottom-up approach is computationally infeasible if both the number of generations and the number of individuals are large. The top-down approach starts from founders and moves down to the target individuals. The computational load is approximately linear in the number of generations. However, the intermediate generalized kinship coefficients need to be stored, which may require extensive storage if the number of individuals in a generation is large. We have implemented a function that allows users to adjust the arguments to achieve a hybrid bottom-up/top-down approach, which can accommodate very deep pedigrees.
QTLRel can estimate variance components given the appropriate incidence matrices. QTLRel estimates these variance components using maximum likelihood. These estimates are nearly equivalent to those obtained by restricted maximum likelihood for typical sample sizes. The maximum likelihood estimates are found numerically using one of several methods. We default to Nelder-Mead since we have found it to be more numerically stable. QTLRel allows users to select variance components using a model selection procedure or perform statistical significance tests for them.
Re-estimating variance components at each marker in a genome scan may not be computationally feasible. The approach used by QTLRel is to first estimate the correlation matrix due to polygenic, residual and other random effects, which is based on the estimated variance components, and then use this matrix as known to scan the genome. Testing fixed effects conditional on estimated random effects is a general approach in mixed-effect models [2, 3, 8, 9]. Our method is most similar to the measured-genotype fixed-heritability method in Aulchenko et al. .
Empirical significance thresholds
QTLRel implements two methods for estimating genome-wide significance thresholds. The first is a permutation test in which the genotypes are permuted while the phenotypes and incidence matrices are held constant. We have previously demonstrated that when polygenic effects are ignored in the model type I error rates are inflated when a permutation is used; however, when the model is appropriate, permutation performs well . The second method is gene dropping which can appropriately control type I error rates even when polygenic variation is ignored in the model.
QTLRel has been successfully used in an AIL to identify QTL for methamphetamine sensitivity , muscle weight , prepulse inhibition  and body weight . In the analysis of body weight we calculated all five genetic incidence matrices, but only estimated additive and dominance variance components because the other variance components are negligible in general . We then performed a genome scan using 688 individual genotyped at 3,105 SNP markers. We identified 11 QTL that exceeded the .05 genome-wide significance threshold estimated from 1800 gene dropping samples. These 11 QTL were confirmed as distinct signals using a forward step-wise multiple-QTL mapping function implemented in QTLRel. We also investigated genome-wide epistatic effects but found none. The analysis was accomplished on an Linux machine equipped with a 2.40GHz Intel(R) Core(TM)2 Quad CPU. It took 295 seconds to estimate variance components and 70 seconds to perform the genome scan.
QTLRel provides a toolkit for genome-wide association studies that is capable of calculating genetic incidence matrices from pedigrees, estimating variance components, performing genome scans, and estimating significance thresholds. It can model interactive covariates and multiple genetic and non-genetic variance components. Other functions include multiple-QTL mapping and genome-wide epistasis. QTLRel can perform interval mapping based on the Haley-Knott method  for markers with 2 alleles. Because QTLRel is implemented in R users can take advantage of numerous other statistical packages; however, there is room to improve on QTLRel's speed since it makes use of many intermediate R functions. It is our intention to extend some functionalities, e.g., the Haley-Knott method to markers with more than two alleles. In summary, QTLRel provides a stand-alone, comprehensive tool to perform QTL analyses in populations in which relatedness is a concern.
Broman KW, Wu H, Sen S, Churchill GA: R/qtl: QTL mapping in experimental crosses. Bioinformatics. 2003, 19: 889-890. 10.1093/bioinformatics/btg112.
Bradbury PJ, Zhang Z, Kroon DE, Ramdoss TMCY, Buckler ES: TASSEL: software for association mapping of complex traits in diverse samples. Bioinformatics. 2003, 23: 2633-2635.
Kang HM, Zaitlen NA, Wade CM, Kirby A, Daly MJ, Eskin E: Efficient control of population structure in model organism association mapping. Genetics. 2008, 178: 1709-1723. 10.1534/genetics.107.080101.
Jackquard A: The genetics structure of populations. 1974, Springer-Verlag, NY
Lynch M, Walsh B: Genetics and analysis of quantitative traits. 1998, Sinauer Associates, Inc, 5:
Abney M: A graphical algorithm for fast computation of identity coefficients and generalized kinship coefficients. Bioinformatics. 2009, 25: 1561-1563. 10.1093/bioinformatics/btp185.
Karigl G: A recursive algorithm for the calculation of identity coefficients. Ann Hum Genet. 1981, 45: 299-305. 10.1111/j.1469-1809.1981.tb00341.x.
Pinheiro JC, Bates DM: Mixed-effects models in S and S-PLUS. 2000, Springer-Verlag, New York
Aulchenko YS, de Koning DJ, Haley C: Genomewide rapid association using mixed model and regression: a fast and simple method for genome-wide pedigree-based quantitative trait loci association analysis. Genetics. 2009, 177: 577-585.
Cheng R, Lim JE, Samocha KE, Sokoloff G, Abney M, Skol AD, Palmer AA: Genome-wide association studies and the problem of relatedness among advanced intercross lines and other highly recombinant populations. Genetics. 2010, 185: 1033-1044. 10.1534/genetics.110.116863.
Lionikas A, Cheng R, Lim JE, Palmer AA, Blizard DA: Fine-mapping of muscle weight QTL in LG/J and SM/J intercrosses. Physiol Genomics. 2010, 42A: 33-38. 10.1152/physiolgenomics.00100.2010.
Samocha KE, Lim JE, Cheng R, Sokoloff G, Palmer AA: Fine mapping of QTL for prepulse inhibition in LG/J and SM/J mice using F2 and advanced intercross lines. Genes Brain Behav. 2010, 9: 759-767. 10.1111/j.1601-183X.2010.00613.x.
Parker CC, Cheng R, Lim JE, Skol AD, Abney M, Palmer AA: Fine-mapping quantitative trait loci for body weight in LG/J × SM/J F2 and F34 advanced intercross lines. Mamm Genome. 2011
Abney M, McPeek MS, Ober C: Estimation of variance components of quantitative traits in inbred populations. Am J Hum Genet. 2000, 141: 629-650.
Haley CS, Knott SA: A simple regression model for interval mapping in line crosses. Heredity. 1992, 69: :315-324. 10.1038/hdy.1992.131.
We appreciate helpful comments from the anonymous reviewers. This project was supported by NIH grants R01DA021336, R01MH079103 and R21DA024845.
RC wrote the program, this paper and the tutorial and has made seminal intellectual contributions throughout this project. MA worked with the other authors and lent his considerable experience addressing similar issues in human populations such as the Hutterites. AAP initiated the projects that lead to the development of the software and has interacted extensively with all other authors. ADS worked with the other authors and lent his considerable experience with the analysis of animal breeding designs and helped to develop earlier versions of both this paper and the tutorial. All authors read and approved the final manuscript.
About this article
Cite this article
Cheng, R., Abney, M., Palmer, A.A. et al. QTLRel: an R Package for Genome-wide Association Studies in which Relatedness is a Concern. BMC Genet 12, 66 (2011). https://doi.org/10.1186/1471-2156-12-66