Volume 3 Supplement 7
Genetic Analysis Workshop 16
Genomewide association and linkage analysis of quantitative traits: comparison of likelihoodratio test and conditional score statistic
 Audrey E Hendricks^{1}Email author,
 Yanyan Zhu^{1} and
 Josée Dupuis^{1}
DOI: 10.1186/175365613S7S100
© Hendricks et al; licensee BioMed Central Ltd. 2009
Published: 15 December 2009
Abstract
Over the past decade, genetic analysis has shifted from linkage studies, which identify broad regions containing putative trait loci, to genomewide association studies, which detect the association of a marker with a specific phenotype. Because linkage and association analysis provide complementary information, developing a method to combine these analyses may increase the power to detect a true association. In this paper we compare a linkage score and association score test as well as a newly proposed combination of these two scores with traditional linkage and association methods.
Background
Improvement in genotyping technologies has led to great advances in the quest to map genes influencing complex traits. In the late 1980s came linkage studies in family samples that identified broad regions containing putative trait loci. Recently, dense singlenucleotide polymorphism (SNP) chip technology has resulted in genomewide association analysis, where the genome is queried for association with a specific phenotype. The high number of SNPs run (from 300,000 to >1,000,000) enables relatively thorough coverage of the genome, but also greatly increases the chance of falsepositive results. A low pvalue in the range of 10^{8} is often used to declare genomewide significance and finding small to moderate associations remains difficult. One advantage of association analysis is it can be carried out in samples of unrelated individuals, which may be easier to recruit. On the other hand, family samples provide extra information about segregation of the phenotype, and both linkage and association analysis may be performed when genotype and phenotype data are available on family members.
Variancecomponent analysis [1] is a commonly used approach for performing linkage analysis of quantitative traits. Its great flexibility to accommodate extended pedigrees is offset by increased type I error of the likelihoodratio test (LRT) when the trait is not normally distributed [2]. An alternate approach consists of using the efficient score statistic for linkage analysis standardized by its variancecomputed conditional on the observed phenotype [3], which renders it robust to departure from multivariate normality.
Populationbased methods can be applied to family samples provided that familial correlation is accommodated and there is no population stratification (or analyses are performed after adjustment for population stratification). Linear mixedeffect models are particularly well suited for association analysis of quantitative traits in family samples. The genotype effect and additional covariates are modeled as fixed effects while familial correlation is accommodated with a random effect with the covariance structure depending on the degree of relatedness between individuals. LRTs may be applied to determine the genotypephenotype association. Alternatively, an efficient score statistic for association analysis with the variance computed conditional on the observed phenotype may offer a robust option against departure from the normality assumption [4].
The objectives of this paper are twofold: first, we compare the LRT with conditional score statistics for both linkage and association analyses using the Genetic Analysis Workshop 16 simulated Framingham Heart Study dataset (Problem 3). Second, we explore whether a combination of the linkage and association score improves power to detect associated SNPs.
Methods
Framingham Heart Study sample
The subset of simulated Framingham Heart Study samples analyzed in our report consists of 736 families and 381 unrelated individuals. The families range in size from 2 to 291 participants with available genotypes and phenotypes, for a total of 6372 individuals.
Initial data cleaning was performed using the software PLINK [5, 6] based on call rate and sex. We removed subjects with a homozygosity rate on the X chromosome between 20% and 80% indicating ambiguous genetic sex, subjects with a conflicting reported and genetic sex, or subjects with a call rate <97%. This resulted in the exclusion of 104 participants.
Lipidrelated phenotypes were available at three exams for each participant. To evaluate methods robust to departure from normality, we concentrated on lowdensity lipoprotein (LDL) and triglyceride levels (TG) levels because of their skewed distributions. We analyzed each trait measured at the first exam, but also averaged each trait over the three exams.
Chromosome selection
Chromosome 11 was chosen because a group of 39 polygenes influencing highdensity lipoprotein (HDL) clustered in the range of 110 Mb to 134 Mb and because of the presence of two major genes, one influencing LDL and the other influencing TG. Chromosome 22 was selected for ease of computation and because it harbors a major gene explaining 1% of the LDL variance.
Linkage analysis: marker selection
For linkage analysis, we selected markers with low linkage disequilibrium, which may bias identitybydescent (IBD) estimates when parental information is unavailable [7]. We used PLINK [5, 6] to select markers on chromosomes 11 and 22 with call rate ≥ 98%, minor allele frequency of ≥ 35%, a HardyWeinberg equilibrium pvalue ≥ 0.05, and a r^{2} < 0.04. Mendelian transmission errors were detected and corrected using the program PEDCHECK [8].
Linkage analysis: variancecomponent model
The basic variancecomponent model assumes that the vector of phenotype in the k^{th} family, Y, has a multivariate normal distribution with mean E(YX) = μ_{ X }= μ + βX and covariance matrix: , where X is the covariate of interest; is the quantitative trait locus, and and are the polygenic and residuals variance components, respectively; π_{ tij }is the proportion of alleles shared IBD by relatives i and j at genome location t; and ϕ_{ ij }is the kinship coefficient. The covariance matrix unconditional on the IBD (Σ) is obtained by setting π_{ tij }= ϕ_{ ij }in the expression for Σ_{ π }.
The log likelihood for a QTL at t is , where the sum is taken over all N pedigrees. To test the null hypothesis of no linkage (H_{0}: = 0), one can use the following LRT: . The typical logarithm of odds (LOD) score is computed by dividing this LRT statistic by twice the natural logarithm of 10. The same hypothesis can be tested using the efficient score statistic, which is obtained by taking the first derivative of the log likelihood ratio with respect to , evaluated at = 0: , where A_{ π }is the matrix of centered IBD and W = (Σ^{1}(Y  μ_{ X })(Y μ_{ X })'I)Σ^{1} with elements w_{ ij }. The squared of the efficient score is standardized by an estimate of the variance conditional on the observed phenotypes to make it robust to violation of the normality assumption [3, 8].
Association analysis
The linear mixedeffects model to test for association analysis is very similar to the variancecomponent model, with the exception that the genotype effect is included in the mean rather than in the covariance matrix: E(Y g, X) = μ + βX + γg, where g is the coded genotypes. The covariance unconditional on the IBD proportions, Σ, is typically used, although one could easily substitute Σ_{ π }at the expense of added computation complexity. The null hypothesis of no association H_{0}: γ = 0 is typically tested using a LRT. Alternatively, the efficient score test may be constructed to test for association. The efficient score statistic for association reduces to , and the estimate of variance, conditional on the observed phenotypes, is , where Φ is the matrix of kinship coefficients. We refer to this statistic as the "conditional score" for association [9].
Combined linkage and association score
The conditional association score and our implementation of the linkage score are uncorrelated under the null hypothesis of no linkage and no association [10]. Therefore, we summed the chisquare form of these two statistics to produce a combined linkage and association statistic. The asymptotic null distribution of the combined statistic is an equally weighted mixture of chisquare with 1 and 2 degrees of freedom (df) for the additive genetic model and a chisquare mixture with 2 and 3 df for the general 2df genetic model.
Evaluation of methods
We computed the LRT and associationconditional score statistics using an additive genetic model and a general 2df genetic model for each SNP. For all models we adjusted for sex and age, or average age in the case of averaged phenotypes. We compared the association statistics in terms of type I error and power. To determine type I error, we selected SNPs with a call rate above 95%, a HWE pvalue above 10^{6}, and a low r^{2} (< 0.01) with all polygenes and major genes on chromosomes 11 and 22. We calculated power as the proportion of significant major gene association detected over all 200 phenotype replicates.
We computed multipoint IBD probabilities using the software LOKI [11], making full use of the entire pedigree and all available genotypes. We used the software R [12] to implement the score statistics and the KINSHIP [13] package in R to compute the LRT statistic.
Results
Type I error
LDLVisit 1  TGAverage  

0.050  0.010  0.0010  0.00010  0.050  0.010  0.0010  0.00010  
LRT (additive)  
MAF>1%  0.048  0.012  0.0019  0.00010  0.051  0.012  0.0013  0.00015 
MAF>10%  0.047  0.012  0.0018  0.00013  0.050  0.012  0.0010  0.00013 
LRT (general 2 df)  
MAF>1%  0.053  0.012  0.0012  0.00031  0.055  0.015  0.0039  0.00191 
MAF>10%  0.050  0.011  0.0011  0.00025  0.050  0.012  0.0009  0.00013 
Association conditional score (additive)  
MAF>1%  0.045  0.010  0.0013  0.00005  0.050  0.010  0.0012  0.00010 
MAF>10%  0.044  0.009  0.0014  0.00006  0.048  0.010  0.0009  0.00013 
Association conditional score (general 2 df)  
MAF>1%  0.065  0.016  0.0029  0.00118  0.073  0.023  0.0060  0.00282 
MAF>10%  0.061  0.014  0.0013  0.00019  0.067  0.017  0.0021  0.00031 
Discussion
Despite selecting phenotypes that violated the normality assumption, type I error rates for the additive associationconditional score and LRT statistics were not inflated, but slight type I error inflation was observed for the 2df LRT statistic and moderate inflation for the 2df associationconditional score. The type I error rate inflation for the 2df statistics is most likely due to a violation of asymptotic assumptions caused by SNPs with low cell counts. Note that filtering by MAF less than 10% lowers the falsepositive rate for both 2df statistics.
In our analyses, all additive statistics failed to detect the association between TG levels and rs603446, a SNP with a nonadditive effect. This result suggests, as do Lettre et al. [14], that the additive model has low power to detect associations for genetic models in which the association is limited to only one genetic group, such as recessive or overdominant genetic models. Although power is lost when the general 2df model is used to detect additive SNPs, the reduction is much smaller than the power reduction when using the additive model to detect the overdominant SNP. Because the genetic model of a SNP is always unknown a priori, the general 2df model seems to be the better model to detect an association.
Including the information from the linkage analysis using the combined conditional score did not consistently increase the power to detect an associated SNP. However, despite the lack of a large linkage signal and even though the combined conditional score statistic uses an extra degree of freedom, the power remained fairly constant when using the combined conditional score compared with the LRT and the associationonly conditional score. The minimum LOD score needed at α = 0.05 to accommodate the extra degree of freedom used by the combined conditional score is 0.467 for the additive model and 0.396 for the general 2df model. The percent of iterations that reached the minimum threshold ranged from 4.0% for the additive model for visit 1 TG levels at SNP rs603446 to 31.5% for the general 2df model for average LDL at SNP rs901824. This suggests that a larger linkage signal or a more efficient combination of the conditional association and linkage scores may ultimately increase power.
A notable difference between the LRT and the conditional score for association was computation time. Performing the analysis for a total of 35,979 SNPs on chromosomes 11 and 22 using a single processor would take more than 16 hours using the LRT and less than an hour using the associationconditional score. The computation times for the linkage analysis were more comparable. Once the IBD was computed, the LRT using SOLAR took about 45 minutes to scan chromosome 11, while the linkageconditional score took about 3.25 hours. Despite the longer linkage computation time, the combined conditional score was computed much faster than the LRT.
Although more research is needed, conditional score analysis provides an interesting possibility of a gain in power by combining linkage information using the combined conditional score. Regardless, the conditional score for association provides a fast and comparable alternative to LRTs for analysis of family data.
List of abbreviations used
 HDL:

Highdensity lipoprotein
 IBD:

Identity by descent
 LDL:

Lowdensity lipoprotein
 LOD:

Logarithm of the odds
 LRT:

Likelihood ratio test
 SNP:

Singlenucleotide polymorphism
 TG:

Triglyceride
Declarations
Acknowledgements
The Genetic Analysis Workshops are supported by NIH grant R01 GM031575 from the National Institute of General Medical Sciences. This research was conducted using the Boston University Linux Cluster for Genetic Analysis (LinGA) funded by the NIH NCRR (National Center for Research Resources) Shared Instrumentation grant (1S10RR16373601A1).
This article has been published as part of BMC Proceedings Volume 3 Supplement 7, 2009: Genetic Analysis Workshop 16. The full contents of the supplement are available online at http://www.biomedcentral.com/17536561/3?issue=S7.
Authors’ Affiliations
References
 Lange K, Westlake J, Spence MA: Extensions to pedigree analysis. III. Variance components by the scoring method. Ann Hum Genet. 1976, 39: 485491. 10.1111/j.14691809.1976.tb00156.x.View ArticlePubMedGoogle Scholar
 Allison DB, Neale MC, Zannolli R, Schork NJ, Amos CI, Blangero J: Testing the robustness of the likelihoodratio test in a variancecomponent quantitativetrait locimapping procedure. Am J Hum Genet. 1999, 65: 531544. 10.1086/302487.PubMed CentralView ArticlePubMedGoogle Scholar
 Tang HK, Siegmund D: Mapping quantitative trait loci in oligogenic models. Biostatistics. 2001, 2: 147162. 10.1093/biostatistics/2.2.147.View ArticlePubMedGoogle Scholar
 Dupuis J, Siegmund DO, Yakir B: A unified framework for linkage and association analysis of quantitative traits. Proc Natl Acad Sci USA. 2007, 104: 2021020215. 10.1073/pnas.0707138105.PubMed CentralView ArticlePubMedGoogle Scholar
 Purcell S, Neale B, ToddBrown K, Thomas L, Ferreira MAR, Bender D, Maller J, Sklar P, de Bakker PIW, Daly MJ, Sham PC: PLINK: a toolset for wholegenome association and populationbased linkage analysis. Am J Hum Genet. 2007, 81: 559575. 10.1086/519795.PubMed CentralView ArticlePubMedGoogle Scholar
 PLINK...Wholegenome Association Toolset. [http://pngu.mgh.harvard.edu/purcell/plink/]
 Schaid DJ, McDonnell SK, Wang L, Cunningham JM, Thibodeau SN: Caution on pedigree haplotype inference with software that assumes linkage equilibrium. Am J Hum Genet. 2002, 71: 992995. 10.1086/342666.PubMed CentralView ArticlePubMedGoogle Scholar
 O'Connell JR, Weeks DE: PedCheck: a program for identification of genotype incompatibilities in linkage analysis. Am J Hum Genet. 1998, 63: 259266. 10.1086/301904.PubMed CentralView ArticlePubMedGoogle Scholar
 Dupuis J, Shi J, Manning AK, Benjamin EJ, Meigs JB, Cupples LA, Siegmund D: Mapping quantitative traits in unselected families: algorithms and examples. Genet Epidemiol. 2009, 33: 617627. 10.1002/gepi.20413.PubMed CentralView ArticlePubMedGoogle Scholar
 Almasy L, Blangero J: Multipoint quantitativetrait linkage analysis in general pedigrees. Am J Hum Genet. 1998, 62: 11981211. 10.1086/301844.PubMed CentralView ArticlePubMedGoogle Scholar
 Heath SC: Markov chain Monte Carlo segregation and linkage analysis for oligogenic models. Am J Hum Genet. 1997, 61: 748760. 10.1086/515506.PubMed CentralView ArticlePubMedGoogle Scholar
 R Development Core Team: R: A language and environment for statistical computing. Vienna, Austria. 2006, ISBN 3900051070, [http://www.Rproject.org]Google Scholar
 Atkinson B, Therneau T: Kinship: mixedeffects Cox models, sparse matrices, and modeling data from large pedigrees. R package. 2007, [http://cran.rproject.org/web/packages/kinship/index.html] , version 1.1.018Google Scholar
 Lettre G, Lange C, Hirschhorn JN: Genetic model testing and statistical power in populationbased association studies of quantitative traits. Genet Epidemiol. 2007, 31: 358362. 10.1002/gepi.20217.View ArticlePubMedGoogle Scholar
Copyright
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.