Abstract
Background
In genomewide association studies, it is widely accepted that multilocus methods are more powerful than testing singlenucleotide 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 windowbased 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 genomewide casecontrol association studies.
Background
In genomewide association studies (GWAS), it is generally accepted that multilocus methods can obtain better power than singlelocus approaches that test only one singlenucleotide polymorphism (SNP) at a time [13]. 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 singlelocus statistics. Then the scan statistic, defined as the largest moving sum, is calculated and its associated empirical pvalue 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 variablesized slidingwindow 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 windowbased 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:
where the mdimensional vector x_{i }represents the ith observation and p_{i }is the probability of observing the ith outcome y_{i }= 1. The regression coefficients parameter can be estimated by maximizing the following log likelihood:
Although logistic regression is very popular in casecontrol 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 (L_{2}) 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 genegene 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:
Where p is a nonnegative definite penalty matrix and λ is a positive scale constant, that is, a tuning parameter, which can be specified by crossvalidation. The regression coefficients in the model are estimated using the NewtonRaphson 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 pvalues 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 diseaseaffected status is denoted by 1, while the diseaseunaffected 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 Mdimensional vector representing sample phenotypes. The procedure of our method is described below, and an intuitive flowchart is displayed in Figure 1.
Figure 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 pvalue for final model and switch to next region.
Initially, we compute the pvalues P_{1}, P_{2}, ..., P_{N }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:
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 r^{2}. In this study, we only consider physical position information. Assume that the positions (base pairs) of SNPs in the region G_{i}(u,v) are denoted by W_{i}(u, v) = {D_{iu}, D_{i  u + 1}, ..., D_{i}, ..., D_{i + v + 1}, D_{i + v}}. Then the nonnegative definite penalty matrix p is defined by:
To avoid introducing excess noise, we do not use the information from all SNPs in each search region, because those with singlelocus pvalues too large (i.e. marginal effects too low) may contribute little to the power. To do this, we denote t as the pvalue truncation threshold, and those SNPs with singlelocus pvalues 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 G_{i}(u,v), the adjusted search region is denoted by T_{i}(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 B_{i}(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, T_{i}(u,v), and add it to B_{i}(u,v). Next, the remaining SNPs in T_{i}(u,v) are entered into the GRLR model one by one, along with the SNPs in B_{i}(u,v) as predictors. We select the one SNP which can reduce the original AIC most. Then we add it to B_{i}(u,v), remove it from T_{i}(u,v), and update the current AIC. We repeat this procedure until none of the remaining SNPs in T_{i}(u,v) can decrease the current AIC. Finally we investigate the last model and calculate the corresponding pvalue for the model fitness, where the subset E_{i }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 . To assign a global empirical pvalue for the selected subset of SNPs, we use L permutations by randomly switching the casecontrol labels in the observed dataset. Considering the computational burden, we can construct search regions of interest using the top S SNPs whose pvalues are smallest among all the pvalues computed by the logistic regression. An alternative is to predefine another truncation threshold and exclude the SNPs whose pvalues 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 diseasecausing, while in Scenario (B), three SNPs are diseasecausing. 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:
In this model, the predictor x_{ij }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 variancecovariance 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 pvalue. Then the test statistic is defined as the minimum value over all windows. TPM is similar to FPM but it focuses on the pvalues that are no more than a prefixed 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 crossvalidation. 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 b^{th }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 tenfold crossvalidations 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 subdatasets 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.
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.
Figure 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 (yaxis). 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 singlelocus 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 genomewide datasets: (1) a casecontrol dataset for heroin addiction [23] with approximately 10,000 SNPs and 200 individuals (100 cases and 100 controls), and (2) a casecontrol dataset for agerelated 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 wholegenome data analysis toolset PLINK [26]. Table 3 shows the original pvalues of the top 10 SNPs, as well as the adjusted pvalues 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 singlelocus pvalues 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 pvalues 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 pvalue of 0.027. It is interesting to note that using logistic regression or other singlelocus methods, pvalues 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 pvalues of the top 10 SNPs, as well as adjusted pvalues for multiple testing corrections. For our GRLR method, we construct search regions based on the SNPs whose singlelocus pvalues 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 pvalue 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 pvalue 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 pvalue 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. Crossvalidation is an important tool for determining the tuning parameters in penalized regression approaches. It would be reasonable to determine λ by crossvalidation, 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 variablesize slidingwindow 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 highdimensional 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 pvalues do not exceed some predefined 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 windowbased 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 casecontrol 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

Hoh J, Ott J: Mathematical multilocus approaches to localizing complex human trait genes.
Nat Rev Genet 2003, 4:701709. PubMed Abstract  Publisher Full Text

McCarthy MI, Abecasis GR, Cardon LR, Goldstein DB, Little J, Ioannidis JPA, Hirschhorn JN: Genomewide association studies for complex traits: consensus, uncertainty and challenges.
Nat Rev Genet 2008, 9:356369. PubMed Abstract  Publisher Full Text

Hoggart CJ, Whittaker JC, De Iorio M, Balding DJ: Simultaneous analysis of all SNPs in genomewide and resequencing association studies.

Hoh J, Ott J: Scan statistics to scan markers for susceptibility genes.
Proc Natl Acad Sci USA 2000, 97:96159617. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Zaykin DV, Zhivotovsky LA, Westfall PH, Weir BS: Truncated product method for combining Pvalues.
Genet Epidemiol 2002, 22:170185. PubMed Abstract  Publisher Full Text

Dudbridge F, Koeleman BP: Rank truncated product of Pvalues, with application to genomewide association scans.
Genetic Epidemiology 2003, 25:360366. PubMed Abstract  Publisher Full Text

Yang H, Hsieh H, Fann CSJ: Kernelbased association test.
Genetics 2008, 179:10571068. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Sun YV, Levin AM, Boerwinkle E, Robertson H, Kardia SL: A scan statistic for identifying chromosomal patterns of SNP association.
Genetic Epidemiology 2006, 30:627635. PubMed Abstract  Publisher Full Text

Browning SR: Multilocus association mapping using variableLength markov chains.
Am J Hum Genet 2006, 78:903913. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Li Y, Sung W, Liu JJ: Association mapping via regularized regression analysis of singlenucleotidepolymorphism haplotypes in variablesized sliding windows.
Am J Hum Genet 2007, 80:705715. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

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

Cessie SL, Houwelingen JCV: Ridge estimators in logistic regression.
Journal of the Royal Statistical Society Series C (Applied Statistics) 1992, 41:191201.

Malo N, Libiger O, Schork NJ: Accommodating linkage disequilibrium in geneticassociation analyses via ridge regression.
Am J Hum Genet 2008, 82:375385. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Gray RJ: Flexible methods for analyzing survival data using splines, With applications to breast cancer prognosis.
Journal of the American Statistical Association 1992, 87:942951. Publisher Full Text

Akaike H: A new look at the statistical model identification.
Automatic Control, IEEE Transactions on 1974, 19:723.
716
Publisher Full Text 
Fisher RA: Statistical methods for research workers. 14th edition. New York: Oliver and Boyd; 1970.

Tibshirani R: Regression shrinkage and selection via the lasso.
Journal of the Royal Statistical Society Series B (Methodological) 1996, 58:267288.

Wu TT, Chen YF, Hastie T, Sobel E, Lange K: Genomewide association analysis by lasso penalized logistic regression.
Bioinformatics 2009, 25:714721. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Zou H, Hastie T: Regularization and variable selection via the elastic net.
Journal of the Royal Statistical Society Series B 2005, 67:301320. Publisher Full Text

Cho S, Kim H, Oh S, Kim K, Park T: Elasticnet regularization approaches for genomewide association studies of rheumatoid arthritis.
BMC Proceedings 2009, 3:S25. PubMed Abstract  PubMed Central Full Text

Wu J, Devlin B, Ringquist S, Trucco M, Roeder K: Screen and clean: a tool for identifying interactions in genomewide association studies.
Genetic Epidemiology 2010, 34:275285. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Devlin B, Roeder K, Wasserman L: Analysis of multilocus models of association.
Genetic Epidemiology 2003, 25:3647. PubMed Abstract  Publisher Full Text

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:417428. PubMed Abstract  Publisher Full Text

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 agerelated macular degeneration.
Science 2006, 314:989992. PubMed Abstract  Publisher Full Text

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 agerelated macular degeneration.
Science 2005, 308:385389. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Purcell S, Neale B, ToddBrown K, Thomas L, Ferreira MA, Bender D, Maller J, Sklar P, de Bakker PI, Daly MJ, Sham PC: PLINK: a tool set for wholegenome association and populationbased linkage analyses.
The American Journal of Human Genetics 2007, 81:559575. Publisher Full Text

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

Friedman J, Hastie T, Tibshirani R: Regularization paths for generalized linear models via coordinate descent.
Journal of Statistical Software 2010, 33:122. PubMed Abstract  PubMed Central Full Text