Email updates

Keep up to date with the latest news and content from BMC Bioinformatics and BioMed Central.

Open Access Methodology article

Multilocus association mapping using generalized ridge logistic regression

Zhe Liu1, Yuanyuan Shen2 and Jurg Ott3*

Author Affiliations

1 Department of Statistics, University of Chicago, 5734 S. University Avenue, Chicago, IL 60637, USA

2 Department of Biostatistics, Harvard School of Public Health, 655 Huntington Avenue, Boston, MA 02115, USA

3 Key Laboratory of Mental Health, Institute of Psychology, Chinese Academy of Sciences, 4A Datun Road, Beijing 100101, China

For all author emails, please log on.

BMC Bioinformatics 2011, 12:384  doi:10.1186/1471-2105-12-384

The electronic version of this article is the complete one and can be found online at: http://www.biomedcentral.com/1471-2105/12/384


Received:16 October 2010
Accepted:29 September 2011
Published:29 September 2011

© 2011 Liu et al; licensee 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.

Abstract

Background

In genome-wide association studies, it is widely accepted that multilocus methods are more powerful than testing single-nucleotide polymorphisms (SNPs) one at a time. Among statistical approaches considering many predictors simultaneously, scan statistics are an effective tool for detecting susceptibility genomic regions and mapping disease genes. In this study, inspired by the idea of scan statistics, we propose a novel sliding window-based method for identifying a parsimonious subset of contiguous SNPs that best predict disease status.

Results

Within each sliding window, we apply a forward model selection procedure using generalized ridge logistic regression for model fitness in each step. In power simulations, we compare the performance of our method with that of five other methods in current use. Averaging power over all the conditions considered, our method dominates the others. We also present two published datasets where our method is useful in causal SNP identification.

Conclusions

Our method can automatically combine genetic information in local genomic regions and allow for linkage disequilibrium between SNPs. It can overcome some defects of the scan statistics approach and will be very promising in genome-wide case-control association studies.

Background

In genome-wide association studies (GWAS), it is generally accepted that multilocus methods can obtain better power than single-locus approaches that test only one single-nucleotide polymorphism (SNP) at a time [1-3]. Among a large number of mathematical and statistical approaches considering many predictors simultaneously, scan statistics [4] serve as a useful multilocus analytical means for combining genetic information on multiple contiguous SNPs. In this method, the whole genome is scanned by a sliding window with estimated size, and the moving sum for each window is computed as the sum of suitable single-locus statistics. Then the scan statistic, defined as the largest moving sum, is calculated and its associated empirical p-value is evaluated by permutation tests.

Despite its remarkable advantages, the scan statistics method has two drawbacks that can restrict its practical use: (I) linkage disequilibrium (LD) within local genomic regions is not taken into account, which can result in an inflated type I error rate; (II) all contiguous SNPs within a genomic region (window) are selected simultaneously, which can bring excess noise and increase the false discovery rate (FDR).

During the last decade, various advances based on the framework of scan statistics have been developed. Zaykin et al. [5], Dudbridge et al. [6], and Yang et al. [7] proposed more effective and powerful test statistics within each sliding window and improved the sensitivity of scan statistics; Sun et al. [8] took into account the complex distribution of human genomic variation in the detection of causal chromosomal regions; Browning [9] and Li et al. [10] proposed a variable-sized sliding-window method based on Markov chain and regularized regression analysis. In their method, there is no need to specify a window size for haplotype tests, which makes it particularly useful in the investigation of association studies.

In this study, we propose a novel sliding window-based multilocus method for identifying a subset of susceptibility SNPs based on forward variable selection, using generalized ridge logistic regression (GRLR) for model fitness at each step. As a broader and generalized form of ridge logistic regression, GRLR can take advantage of prior information between any pair of SNPs and impose proper shrinkage penalty on each SNP within the genomic region of interest. Our method can automatically combine genetic information within local regions and select a subset of SNPs that best predict disease status, whose associated empirical significance level is assessed by permutation tests. We demonstrate by simulations and analysis of two published datasets that our method is highly informative and promising.

Methods

Generalized ridge logistic regression

Logistic regression (LR) is a common statistical method for looking into the dependency of a binomial response on a number of variables (predictors), whose general form is:

<a onClick="popup('http://www.biomedcentral.com/1471-2105/12/384/mathml/M1','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2105/12/384/mathml/M1">View MathML</a>

where the m-dimensional vector xi represents the ith observation and pi is the probability of observing the ith outcome yi = 1. The regression coefficients parameter <a onClick="popup('http://www.biomedcentral.com/1471-2105/12/384/mathml/M7','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2105/12/384/mathml/M7">View MathML</a> can be estimated by maximizing the following log likelihood:

<a onClick="popup('http://www.biomedcentral.com/1471-2105/12/384/mathml/M2','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2105/12/384/mathml/M2">View MathML</a>

Although logistic regression is very popular in case-control association analysis, it suffers from several shortcomings. If the number of SNPs in the regression model is larger than the number of observations, this method fails [1]. Moreover, with a large number of SNPs in the regression model, predictors can be highly correlated (high linkage disequilibrium), which can lead to further degradation of the model [11]. With the use of quadratic (L2) penalization, ridge logistic regression (RLR) [12] can overcome these disadvantages of logistic regression and increase the stability of model fitness. This method has recently been applied successfully in several biological investigations including accommodating linkage disequilibrium [13] and uncovering gene-gene interactions [11].

For the analysis discussed below, we consider generalized ridge logistic regression (GRLR), a broader and generalized form of ridge logistic regression, which amounts to maximizing the following penalized log likelihood:

<a onClick="popup('http://www.biomedcentral.com/1471-2105/12/384/mathml/M3','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2105/12/384/mathml/M3">View MathML</a>

Where p is a nonnegative definite penalty matrix and λ is a positive scale constant, that is, a tuning parameter, which can be specified by cross-validation. The regression coefficients in the model are estimated using the Newton-Raphson iterative algorithm. The effective degrees of freedom and the variance of the coefficients can be approximated by estimators introduced in [14]. Then the Wald test can be applied to assign p-values to the regression coefficients.

Due to the feature of quadratic penalization that none of the coefficient estimators would be equal to zero in the shrinkage, GRLR cannot serve as an independent tool for model selection; however, the traditional forward selection procedure can be utilized, with the use of GRLR for model fitness in each round. In forward variable selection, we start with no predictor in the model and then add the one variable that leads to the best score. We continue adding variables one at a time until the score stops improving. In this study, we choose the AIC (Akaike Information Criterion) [15] as the scoring method in the variable selection procedure, which measures goodness of fit of a statistical model.

As a proof of concept, we assume a sample of M individuals genotyped at N SNPs. For each individual, the three genotypes AA, AB, and BB on any given SNP are denoted by 0, 1, and 2, respectively, where the B allele is the minor allele, and an additive model is assumed. For the phenotype, the disease-affected status is denoted by 1, while the disease-unaffected status is denoted by 0. Let x be an N × M matrix containing the sample genotypes, where N is the number of SNPs and M is the number of individuals. Let y be an M-dimensional vector representing sample phenotypes. The procedure of our method is described below, and an intuitive flowchart is displayed in Figure 1.

thumbnailFigure 1. Summary of steps in GRLR method. Step 1: Construct search region; Step 2: SNP truncation; Step 3: Apply forward selection to the region, fitting the GRLR model in each iteration; Step 4: Calculate the p-value for final model and switch to next region.

Initially, we compute the p-values P1, P2, ..., PN for all the N SNPs, using logistic regression. For the ith SNP, where i = 1, 2, ..., N, we construct a search region with that SNP being in the center of the region. Let u and v denote the sizes (i.e. the numbers of SNPs) before and after this SNP, where u and v are always equal as long as the search region is entirely located in a chromosome. The constructed search region with the ith SNP as the central SNP is denoted by:

<a onClick="popup('http://www.biomedcentral.com/1471-2105/12/384/mathml/M4','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2105/12/384/mathml/M4">View MathML</a>

To define the penalty matrix p in the GRLR model, we use the physical position information between pairs of SNPs in the search region, in order to impose different weights on different SNPs based on their distances to the center. The larger the distance between a SNP and the center, the more penalties would be placed on it. Other weighting functions for pairwise SNPs can also be implemented, such as the LD measurements D or r2. In this study, we only consider physical position information. Assume that the positions (base pairs) of SNPs in the region Gi(u,v) are denoted by Wi(u, v) = {Di-u, Di - u + 1, ..., Di, ..., Di + v + 1, Di + v}. Then the nonnegative definite penalty matrix p is defined by:

<a onClick="popup('http://www.biomedcentral.com/1471-2105/12/384/mathml/M5','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2105/12/384/mathml/M5">View MathML</a>

To avoid introducing excess noise, we do not use the information from all SNPs in each search region, because those with single-locus p-values too large (i.e. marginal effects too low) may contribute little to the power. To do this, we denote t as the p-value truncation threshold, and those SNPs with single-locus p-values larger than t will be excluded from the search region. In practical applications, the threshold t may be set somewhere from 0.05 to 0.10. After the truncation on Gi(u,v), the adjusted search region is denoted by Ti(u,v).

With the above configuration, we apply forward selection to each search region Ti(u,v), using GRLR for model fitting in each step. Denote Bi(u,v) as the current best subset of SNPs for each step of the model selection procedure, which is empty at the beginning. We start by fitting the GRLR model using the central SNP (i.e. the ith SNP) as the only predictor, and calculating the corresponding AIC. Then we remove the central SNP from the search region, Ti(u,v), and add it to Bi(u,v). Next, the remaining SNPs in Ti(u,v) are entered into the GRLR model one by one, along with the SNPs in Bi(u,v) as predictors. We select the one SNP which can reduce the original AIC most. Then we add it to Bi(u,v), remove it from Ti(u,v), and update the current AIC. We repeat this procedure until none of the remaining SNPs in Ti(u,v) can decrease the current AIC. Finally we investigate the last model and calculate the corresponding p-value <a onClick="popup('http://www.biomedcentral.com/1471-2105/12/384/mathml/M8','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2105/12/384/mathml/M8">View MathML</a> for the model fitness, where the subset Ei stands for the selected SNPs in this search region.

Considering all the search regions over the whole genome, our test statistic is defined as the minimum of all <a onClick="popup('http://www.biomedcentral.com/1471-2105/12/384/mathml/M9','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2105/12/384/mathml/M9">View MathML</a>. To assign a global empirical p-value for the selected subset of SNPs, we use L permutations by randomly switching the case-control labels in the observed dataset. Considering the computational burden, we can construct search regions of interest using the top S SNPs whose p-values are smallest among all the p-values computed by the logistic regression. An alternative is to pre-define another truncation threshold and exclude the SNPs whose p-values exceed the threshold.

Simulations

We apply our method to simulated data and compare its performance with that of five other methods in current use. For each replication, we generate a susceptibility genomic region with 21 SNPs (consisting of the central SNP, the upper u = 10 SNPs, and the lower v = 10 SNPs) under two scenarios. In Scenario (A), two SNPs within the region are disease-causing, while in Scenario (B), three SNPs are disease-causing. In both scenarios, one causal SNP is located at the center of genomic region while other causal SNPs are randomly located. We generate 100 independent individuals (50 cases and 50 controls), using the following logistic regression model:

<a onClick="popup('http://www.biomedcentral.com/1471-2105/12/384/mathml/M6','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2105/12/384/mathml/M6">View MathML</a>

In this model, the predictor xij refers to the genotype of the ith individual at the jth causal SNP, where genotypes AA, AB, and BB are coded by 0, 1, and 2, respectively, assuming the B allele is the minor allele. In Scenario (A), k (the number of causal SNPs) equals 2, while in Scenario (B), k equals 3. The intercept βo is determined by requiring a disease prevalence of 0.05; other regression coefficients βj are set to 1.0.

Assume that each SNP within the region is derived from a multivariate normal distribution, and the variance-covariance matrix Σ is defined by: Σpq = 1, when p = q; Σpq = r, when |p - q| ≤ 5; Σpq = 0, otherwise. The correlation coefficient r can be varied in simulations. For each SNP, we determine its genotype by its corresponding value generated from a multivariate normal distribution. If the value falls into the interval (-∞, qnorm((1 - f)2)), the genotype of this SNP is set to AA; if the value falls into the interval (qnorm(1 - f 2), + ∞), the genotype of this SNP is set to BB; otherwise, the genotype of this SNP is set to AB. Here, f is the frequency of the B allele (i.e. minor allele frequency) and the function "qnorm" computes the quartile of standard normal distribution.

Results

Power calculations

We compare the performance of our method (GRLR) with that of logistic regression (LR), Fisher product method (FPM) [16], truncated product method (TPM) [5], lasso logistic regression (Lasso) [10,17,18] and elastic net (Enet) [19,20]. FPM uses a sliding window with fixed size to scan the genomic region, and for each window, it computes the sum of the logarithm of each p-value. Then the test statistic is defined as the minimum value over all windows. TPM is similar to FPM but it focuses on the p-values that are no more than a pre-fixed truncation threshold. In this simulation study, the truncation threshold for TPM method is set to 0.05 and the window size is set to 5 in both FPM and TPM methods. For the Lasso and Enet methods, all SNPs in the search region are entered into regression models and subsets of SNPs are selected. The tuning parameters in both methods are determined by cross-validation. For the GRLR method, the truncation threshold is also set to 0.05 and the tuning parameter λ is set to 1. We mimic the physical distances between different SNPs by constructing a simple mapping function that the position of the bth SNP in the genomic region is set to b.

For each combination of these parameters, with minor allele frequencies f = 0.1, 0.3, or 0.5 and correlation coefficients r = 0.0, or 0.2, we carry out 500 replications and assess the power of various methods under both scenarios. To ensure that the type I error rate of each method is equal to 0.05, we utilize the following control procedure: for the FPM, TPM, and GRLR methods, 500 permutations are applied in each simulated replication. For the Lasso and Enet methods, since they often include extra terms and tend to have a high type I error rate [21,22], we apply 10 times ten-fold cross-validations for each generated replication. For each selected SNP in the observed data, its consistency is defined as the number of times that it occurs in 100 sub-datasets that contain 90 percent of observations. The consistency thresholds are determined under the null hypothesis to ensure the normal type I error rates.

For power considerations, the testing "success" can be defined by different strategies:

Scenario (A) (two causal SNPs)

• Strategy I: at least one of the two SNPs is detected

• Strategy II: both SNPs are detected

Scenario (B) (three causal SNPs)

• Strategy I: at least one of the three SNPs is detected

• Strategy II: at least two of the three SNPs are detected

• Strategy III: all three SNPs are detected

Tables 1 and 2 show the results of power simulations for the six methods under the conditions of Scenario (A) and Scenario (B), respectively. Considering all 30 conditions in the results, the averaged power for the six methods LR, FPM, TPM, Lasso, Enet, and GRLR, are equal to 0.373, 0.400, 0.394, 0.498, 0.539, and 0.558, respectively. Our GRLR method is the overall winner and Enet is second, being slightly less powerful than our method. Lasso performs worse than GRLR and Enet, but clearly outperforms the remaining three methods. FPM and TPM have quite similar power results, and LR ranks last. It is interesting to note that our GRLR method works very well under Strategies II and III, which indicates that our method is especially powerful in uncovering all causal SNPs within a given genomic region.

Table 1. Power calculation under Scenario (A)

Table 2. Power calculation under Scenario (B)

In our GRLR method, arbitrary choices of the truncation threshold and the size of the sliding window may become a potential problem. Thus we investigate the impact of different truncation thresholds t = 0.05 or 0.10, and different sizes of regions 11, 21, or 41 (i.e. u = 5, 10, or 20). We focus on Scenario (B) (three causal SNPs), and set the minor allele frequency to 0.5. Results are shown in Figure 2, where we can see that the power fluctuation caused by different choices of the two parameters is not large. It indicates that our method is not very sensitive to the choice of truncation threshold and the size of the sliding window.

thumbnailFigure 2. Impact of different thresholds on power. The effects of different truncation thresholds t = 0.05 or 0.10, and different sizes of regions 11, 21, or 41 (i.e. u = 5, 10, or 20) on power (y-axis). Only Scenario (B) (three causal SNPs) is considered, and the minor allele frequency is set to 0.50.

We conducted additional simulations to test the performance of our method conditional on various choices of the tuning parameter λ. We consider Scenario (A) (two causal SNPs) and apply two definitions of testing success, Strategy I (requiring that at least one of the two SNPs is significant) and Strategy II (requiring significance at both SNPs). The correlation coefficient and the minor allele frequency are set to 0.0 and 0.30, respectively. The number of replication runs is 500 for all simulations. Results show that, when the tuning parameter λ equals 0.01, 0.1, 1.0, 10, and 100, the power of our proposed method under Strategy I equals 0.678, 0.650, 0.650, 0.640, and 0.672, respectively, while the power under Strategy II equals 0.600, 0.560, 0.558, 0.548, and 0.574, respectively. We can see that the power fluctuation caused by different values of λ is not large. It indicates that our method is not very sensitive to the selection of the tuning parameter, although it is unknown which selection can lead to the best result.

We further compared the performance of single-locus logistic regression and our proposed method by increasing the total number of SNPs in each simulated window to 200 and the number of causal SNPs to 10. For ease of calculation, only 100 permutations were applied in each of 200 replication runs. Simulation results show that, for logistic regression method, the proportion of times that it can successfully detect at least one causal SNP is 0.465, while the proportion of times that it can detect at least two causal SNPs is only 0.110; for our GRLR method, the proportion of times that it can successfully detect at least one causal SNP is 0.455, while the proportion of times that it can detect at least two casual SNPs is 0.445. These results indicate that our method is more powerful in discovering multiple disease variants.

Analyzing published data

To further demonstrate our GRLR method and evaluate its practical applications in real dataset analysis, we apply it to two published genome-wide datasets: (1) a case-control dataset for heroin addiction [23] with approximately 10,000 SNPs and 200 individuals (100 cases and 100 controls), and (2) a case-control dataset for age-related macular degeneration (AMD) collected in Hong Kong [24,25] with approximately 100,000 SNPs and 223 individuals (96 cases and 127 controls).

For the heroin addiction data, we initially apply logistic regression on all the SNPs using the whole-genome data analysis toolset PLINK [26]. Table 3 shows the original p-values of the top 10 SNPs, as well as the adjusted p-values after applying the Bonferroni correction and 1000 permutations. From the result, none of the SNPs are statistically significant at threshold α = 0.05 after multiple test correction. To apply our GRLR method to the dataset, we construct genomic search regions based on the SNPs whose single-locus p-values are no larger than 0.05; the truncation threshold within each search region is also set to 0.05; the size on each side of the central SNP is set to 10; the tuning parameter λ in the GRLR model is set to 1.0. We use 1000 permutations to assign empirical p-values for the selected subset of SNPs. A summary of the analysis results is shown in Table 4, where the identified subset {rs1408830, rs965972} is statistically significant with a global empirical p-value of 0.027. It is interesting to note that using logistic regression or other single-locus methods, p-values for the two SNPs are 2.23 × 10-04 and 5.06 × 10-03, neither of which is statistically significant at α = 0.05 after multiple testing corrections. In contrast, our method automatically combines the association information of the two SNPs and leads to a significant result.

Table 3. Results for heroin addiction data using logistic regression

Table 4. Results for both datasets using GRLR

For the AMD Hong Kong data, after applying logistic regression, only the SNP rs10490924 is statistically significant, whose functional significance had been established experimentally [24]. Table 5 shows the original p-values of the top 10 SNPs, as well as adjusted p-values for multiple testing corrections. For our GRLR method, we construct search regions based on the SNPs whose single-locus p-values are no larger than 0.01. Other parameter configurations are the same as those applied in the first dataset. A summary of the results is shown in Table 4. Our method shows that the identified subset {rs2736911, rs10490924, rs763720} is statistically significant with an empirical p-value of 0 (evaluated by 1000 permutations). Based on this estimate, the 95% confidence interval for p extends from 0 through 0.003, which is considerably smaller than the published p-value of 0.027 (Table 4).

Table 5. Results for AMD Hong Kong data using logistic regression

Discussion

Our method has the following advantages: (I) GRLR automatically combines association information of SNPs within each sliding window, and by truncating SNPs with low marginal effects and penalizing SNPs with a large distance from the center, our method can exclude excess noise and allow for linkage disequilibrium in local genomic regions; (II) we apply a forward model selection framework and fit the GRLR model at each step, since GRLR cannot serve as a means for variable selection. The results of power simulation and real datasets analysis indicate that this procedure works very well; (III) the global empirical p-value for the selected subset of SNPs is evaluated by permutation analysis, which properly handles the multiple testing problem and furnishes a valid type I error rate.

There are still some aspects in our GRLR method that can be improved in future work: (I) we do not consider the selection of the tuning parameter λ in this study. Cross-validation is an important tool for determining the tuning parameters in penalized regression approaches. It would be reasonable to determine λ by cross-validation, although it would introduce a higher computational burden; (II) the choices of window size and the truncation threshold in our method are arbitrary, which may reduce the stability of our method. Although simulation results show that the impact of these parameters on power is not large, a variable-size sliding-window procedure could be a better choice. Another way is to allow for the free variability of these parameters and then determine the best ones by permutations; (III) the Lasso and Enet methods can perform better than our proposed method in high-dimensional model selection, since they can deal with multiple variables simultaneously. In this case, the introduction of a truncation threshold is not quite necessary for the applications of these two methods. For power comparisons, however, we wanted to apply the same preselection procedure and the same truncation threshold to the Lasso and Enet methods. It is of our interest to investigate this further in future work; (IV) it is ideal to apply our region construction procedure to all SNPs. In this study, the search regions are constructed based on the SNPs whose marginal p-values do not exceed some pre-defined threshold; (V) it is of our interest to use the LD measurements or the combination of the LD and the physical distance information in constructing the weights for SNPs.

Conclusions

In this study, we propose a new sliding window-based multilocus approach for detecting causal SNPs, which is based on forward model selection using generalized ridge logistic regression (GRLR) for model fitness at each step. Our method can overcome some defects of the scan statistics approach and provides a powerful procedure for identifying causal genomic region and mapping susceptibility genes in routine case-control association studies. In particular, because of its capability of automatically combining association information of multiple SNPs and its advantage on variable selection, our method can be a useful technique in the analyses of human complex diseases. Our software (available by request) is written in R [27], with the use of the Design Package [28] and glmnet Package [29].

Competing interests

The authors declare that they have no competing interests.

Authors' contributions

ZL and YS participated in the study design and the programming. JO contributed to the overall organization, reviewing and editing the manuscript. All authors read and approved the final manuscript.

Acknowledgements

This work has been supported by NSFC grants from the Chinese government (project numbers 30730057 and 30700442) to JO.

References

  1. Hoh J, Ott J: Mathematical multi-locus approaches to localizing complex human trait genes.

    Nat Rev Genet 2003, 4:701-709. PubMed Abstract | Publisher Full Text OpenURL

  2. McCarthy MI, Abecasis GR, Cardon LR, Goldstein DB, Little J, Ioannidis JPA, Hirschhorn JN: Genome-wide association studies for complex traits: consensus, uncertainty and challenges.

    Nat Rev Genet 2008, 9:356-369. PubMed Abstract | Publisher Full Text OpenURL

  3. Hoggart CJ, Whittaker JC, De Iorio M, Balding DJ: Simultaneous analysis of all SNPs in genome-wide and re-sequencing association studies.

    PLoS Genet 2008, 4. OpenURL

  4. Hoh J, Ott J: Scan statistics to scan markers for susceptibility genes.

    Proc Natl Acad Sci USA 2000, 97:9615-9617. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  5. Zaykin DV, Zhivotovsky LA, Westfall PH, Weir BS: Truncated product method for combining P-values.

    Genet Epidemiol 2002, 22:170-185. PubMed Abstract | Publisher Full Text OpenURL

  6. Dudbridge F, Koeleman BP: Rank truncated product of P-values, with application to genomewide association scans.

    Genetic Epidemiology 2003, 25:360-366. PubMed Abstract | Publisher Full Text OpenURL

  7. Yang H, Hsieh H, Fann CSJ: Kernel-based association test.

    Genetics 2008, 179:1057-1068. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  8. Sun YV, Levin AM, Boerwinkle E, Robertson H, Kardia SL: A scan statistic for identifying chromosomal patterns of SNP association.

    Genetic Epidemiology 2006, 30:627-635. PubMed Abstract | Publisher Full Text OpenURL

  9. Browning SR: Multilocus association mapping using variable-Length markov chains.

    Am J Hum Genet 2006, 78:903-913. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  10. Li Y, Sung W, Liu JJ: Association mapping via regularized regression analysis of single-nucleotide-polymorphism haplotypes in variable-sized sliding windows.

    Am J Hum Genet 2007, 80:705-715. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  11. Park MY, Hastie T: Penalized logistic regression for detecting gene interactions.

    Biostat 2008, 9:30-50. OpenURL

  12. Cessie SL, Houwelingen JCV: Ridge estimators in logistic regression.

    Journal of the Royal Statistical Society Series C (Applied Statistics) 1992, 41:191-201. OpenURL

  13. Malo N, Libiger O, Schork NJ: Accommodating linkage disequilibrium in genetic-association analyses via ridge regression.

    Am J Hum Genet 2008, 82:375-385. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  14. Gray RJ: Flexible methods for analyzing survival data using splines, With applications to breast cancer prognosis.

    Journal of the American Statistical Association 1992, 87:942-951. Publisher Full Text OpenURL

  15. Akaike H: A new look at the statistical model identification.

    Automatic Control, IEEE Transactions on 1974, 19:723.

    716

    Publisher Full Text OpenURL

  16. Fisher RA: Statistical methods for research workers. 14th edition. New York: Oliver and Boyd; 1970. OpenURL

  17. Tibshirani R: Regression shrinkage and selection via the lasso.

    Journal of the Royal Statistical Society Series B (Methodological) 1996, 58:267-288. OpenURL

  18. Wu TT, Chen YF, Hastie T, Sobel E, Lange K: Genome-wide association analysis by lasso penalized logistic regression.

    Bioinformatics 2009, 25:714-721. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  19. Zou H, Hastie T: Regularization and variable selection via the elastic net.

    Journal of the Royal Statistical Society Series B 2005, 67:301-320. Publisher Full Text OpenURL

  20. Cho S, Kim H, Oh S, Kim K, Park T: Elastic-net regularization approaches for genome-wide association studies of rheumatoid arthritis.

    BMC Proceedings 2009, 3:S25. PubMed Abstract | PubMed Central Full Text OpenURL

  21. Wu J, Devlin B, Ringquist S, Trucco M, Roeder K: Screen and clean: a tool for identifying interactions in genome-wide association studies.

    Genetic Epidemiology 2010, 34:275-285. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  22. Devlin B, Roeder K, Wasserman L: Analysis of multilocus models of association.

    Genetic Epidemiology 2003, 25:36-47. PubMed Abstract | Publisher Full Text OpenURL

  23. Nielsen DA, Ji F, Yuferov V, Ho A, Chen A, Levran O, Ott J, Kreek MJ: Genotype patterns that contribute to increased risk for or protection from developing heroin addiction.

    Mol Psychiatry 2008, 13:417-428. PubMed Abstract | Publisher Full Text OpenURL

  24. DeWan A, Liu M, Hartman S, Zhang SS, Liu DTL, Zhao C, Tam POS, Chan WM, Lam DSC, Snyder M, Barnstable C, Pang CP, Hoh J: HTRA1 promoter polymorphism in wet age-related macular degeneration.

    Science 2006, 314:989-992. PubMed Abstract | Publisher Full Text OpenURL

  25. Klein RJ, Zeiss C, Chew EY, Tsai J, Sackler RS, Haynes C, Henning AK, SanGiovanni JP, Mane SM, Mayne ST, Bracken MB, Ferris FL, Ott J, Barnstable C, Hoh J: Complement factor H polymorphism in age-related macular degeneration.

    Science 2005, 308:385-389. PubMed Abstract | Publisher Full Text | PubMed Central Full Text OpenURL

  26. Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira MA, Bender D, Maller J, Sklar P, de Bakker PI, Daly MJ, Sham PC: PLINK: a tool set for whole-genome association and population-based linkage analyses.

    The American Journal of Human Genetics 2007, 81:559-575. Publisher Full Text OpenURL

  27. R Development Core Team: R: a language and environment for statistical computing. Vienna, Austria; 2010. OpenURL

  28. Harrell FE:

    Design: R package version 2.3-0. 2009. OpenURL

  29. Friedman J, Hastie T, Tibshirani R: Regularization paths for generalized linear models via coordinate descent.

    Journal of Statistical Software 2010, 33:1-22. PubMed Abstract | PubMed Central Full Text OpenURL