Combining identity by descent and association in genetic case-control studies

Background In human case-control association studies, one of the chi-square tests typically carried out is based on a 2 × 3 table of genotypes (homogeneity of three genotype frequencies in case and control individuals). We formulate the two degrees of freedom associated with a given genotype distribution in terms of two biologically relevant parameters, (1) the probability F that an individual's two alleles are identical by descent (IBD) and (2) the frequency p of one of the alleles. Results Imposing the restriction, F ≥ 0, makes some of the genotype frequencies invalid thereby reducing noise. We propose a new statistical association test, the FP test, by focusing on allele frequency differences between case and control individuals while allowing for suitable IBD probabilities. Power calculations show that (1) the practice of generally carrying out two association tests (allele and genotype test) has an increased type I error and (2) our test is more powerful than conventional genotype and allele tests under recessive trait inheritance, and at least as powerful as these conventional tests under dominant inheritance. Conclusion For dominant and recessive modes of inheritance, any apparent power gain by an allele test when carried out in conjunction with a genotype test tends to be purchased entirely by an increased rate of false positive results due to omission of a multiple testing correction. As an alternative to these two standard association tests, our FP test represents a convenient and more powerful alternative.


Background
In their well-known paper on homozygosity mapping published 20 years ago, Lander and Botstein [1] showed that fewer than a dozen unrelated inbred children should suffice to map a recessive trait given a dense map of genetic markers. They recommended that one should search for extended regions of homozygosity shared by a set of inbred individuals. The power of this approach was demonstrated by the mapping of a gene for a rare recessive trait in a genome-wide investigation of only three distantly related patients [2]: On chromosome 18, five of the six disease chromosomes shared a region of 19 cM in length.
Even individuals seemingly collected at random from the population tend to exhibit extended regions of allele sharing [3], which may be interpreted as the consequence of the mating of (distantly) related individuals. Such homozygosity represents autozygosity, that is, the sharing of two alleles that are copies of a single ancestral allele [3].
In recent years, researchers have shown renewed interest in extended segments of homozygosity and have generally done so by focusing on segments of specific lengths [4][5][6][7]. Our purpose here is to focus on individual SNPs rather than on genomic segments of arbitrary lengths, and to allow for biologically reasonable values of IBD while testing for allele frequency differences. The reason for this is not so much biological plausibility but rather statistical power: By disregarding parameter values that are unlikely to be of importance, we reduce the effect of statistical noise and thereby gain accuracy and power. It will be seen below that our approach is somewhat analogous to the "possible triangle" method [8] in affected sib-pair linkage analysis and is expected to lead to similar increases in power.

Statistical model for SNP genotype frequencies
Consider a SNP marker with two alleles, A and B, with p = P(A). Let F denote the probability that the two alleles in an individual are identical by descent (IBD) or autozygous (implies homozygous) or, equivalently, F = proportion of individuals who are IBD at that marker, where F is the inbreeding coefficient. Given autozygosity (F = 1), the frequencies of genotypes AA, AB, and BB are given by p, 0, and 1 -p, respectively, and given allozygosity (F = 0) we assume the corresponding frequencies to be p 2 , 2p(1 -p), and (1 -p) 2 (Hardy-Weinberg equilibrium, HWE). Then the genotype probabilities may be formulated as given in Table 1.
As there are equal numbers (i.e., 2) of free (independent) genotype classes and parameters, finding maximum likelihood estimates (MLEs) amounts to simply equating genotype frequencies with their probabilities of occurrence (Table 1). Solving the expressions for q 1 and q 2 in Table 1 for F and p leads to the following solutions: The MLEs of q 1 , q 2 , and q 3 are simply the proportions of individuals with given genotypes. Because of the invariance property of MLEs, the functions p and F of q 1 and q 2 (equation 1) are also MLEs. The inbreeding coefficient F may be viewed as an indicator of how far the genotype frequencies deviate from HWE. If q 2 is smaller than expected under HWE then F > 0.
The inbreeding coefficient F in human populations is known to be small and positive. In North America, it is generally much smaller than 0.01 [9] but in isolated populations may reach value of around 0.10 [10]. However, in samples of individuals affected with a heritable trait, the inbreeding coefficient may be even higher. For example, for the SNP most strongly associated with age-related macular degeneration (AMD) [11], rs380390, application of equations (1) leads to an estimated inbreeding coefficient of 0.13 in case individuals. Such increases occur because of enrichment of a disease genotype due to ascertainment of case individuals, and/or they may be a consequence of the fact that individuals affected with a heritable trait may be distantly related. Conversely, the estimated (unrestricted) F value for rs380390 was -0.07 in control individuals.
So far, the expressions for the genotype frequencies in Table 1 simply amount to a specific formulation of the two df's associated with the three frequencies. Several other such parametrizations have previously been proposed [12]. The current transformation (1) of genotype frequencies does not by itself yield any new insights. However, in the next section we impose restrictions on the range of parameter values, which will result in a new test.

Statistical test
Parametrization of genotype frequencies (F = probability of autozygosity, p = frequency of A allele).
( Table 1), p a and F a are the parameter values in case individuals. For control individuals, the log likelihood log is obtained in an analogous manner. The test where the subscript c refers to the combined data (case and control individual), that is, common parameter values. For unrestricted parameter values, the test statistic T has an asymptotic chi-square distribution under H 0 .
However, because of the conditions imposed (F ≥ 0), T does not follow a chi-square distribution. Therefore, as outlined in the Methods section, we compute associated significance levels numerically by computer-based permutation testing. Another reason for applying permutation tests is as follows. This single test may replace the conventional two association tests (allele and genotype test) in genome-wide association studies, where 100,000s of markers are used. Test results for these markers are not independent, which is optimally taken into account in permutation tests but would be difficult to capture analytically.

Power calculations
To evaluate the performance of our new test with existing tests, we carry out power calculations under a recessive and a dominant model of disease inheritance, where we assume a functional SNP fully associated with the disease variant. Model parameters (penetrances and disease allele frequencies) are calibrated to predict a trait prevalence of 5% for each model. The proportion of affected individuals in the population whose disease is due to the given gene is fixed at 10%. The "strength" of a model is measured by the penetrance ratio, γ, where γ = 1 corresponds to the null hypothesis.
We compare three tests, our FP test, the chi-square genotype test based on a 2 × 3  Figure 1 shows power curves for recessive disease models. Clearly, our FP test has better power than any of the other association tests. For dominant traits (Figure 2), if researchers carried out only the allele test, that test has higher power than the FP and genotype tests. As both the allele and genotype tests are usually carried out at the same time, the proper comparison is between the MaxGA statistic and our FP test, where the latter exhibits a slight advantage over the former. Thus, it is fair to say that the FP test is more powerful than conventional test statistics (allele and genotype tests) under our recessive and dominant inheritance models. We have implemented it in a computer program, FPtest, which is freely available on our website [13].
Under both dominant and recessive models, the MaxGA statistic has essentially the same power as the genotype test. Thus, any superiority of the allele test, if carried out at the same time as the genotype test, is wiped out by the multiple testing correction. As the usual practice is to do both, an allele and a genotype test, the potentially higher power of the allele test is fallacious as it is entirely purchased by an increased rate of false positive results.
We also considered Risch's genotype relative risk model [14], in which the penetrances for the three disease geno-

Application to published data
We applied our FP test to published data on age-related macular degeneration (AMD) [11,15] and Parkinson disease (PD) [16,17]. In these genome-wide studies, the number of SNPs ranged from 100,000 (AMD) to 393,000 (PD). The empirical (experiment-wise) significance level associated with each association test was estimated in 5,000 randomization samples (labels "case" and "control" randomly permuted) as the proportion of such samples exhibiting a test statistic at least as large as the observed test statistic for any one of the SNPs. Table 2 shows the results of the genotype and FP tests in these published data. For completeness, results for the allele test are also shown. Significance levels in a given row of Table  2 were computed on the basis of the same set of randomization samples but different samples were used for different rows. For these calculations, we used the permutation procedure [18] implemented in our sumstat program [19], which employs the most recently published 64 bit random number generator [20].
For each of the three studies in Table 3, all SNPs are listed that achieved an experiment-wise significance level of 0.05 or less in either one of the three association tests. Of the 13 resulting SNPs, 11 show a smaller p-value for the FP test than the genotype test and one SNP shows the same p-value. These results clearly demonstrate the usefulness of our new association test. As expected, observed unrestricted F values are larger in case than control individuals and in the latter are often negative.

Discussion
As mentioned in the introduction, researchers often look for genomic regions of increased homozygosity or autozygosity by sliding a window of fixed length across the genome. Our test offers an elegant alternative to such windows of fixed and arbitrary lengths. We propose to work with scan statistics as previously developed [21]. This method also employs a window of a fixed length (fixed number of SNPs) and determines the maximum of the sum of test statistics for all such windows in the genome, which is the scan statistic of the given length. What sets this approach apart from ad hoc approaches is that it applies different window sizes from 1 up to a specified maximum length and estimates optimal window length by maximum likelihood. An updated version of our scan-stat program is available that incorporates the FP test statistic [22]. This implementation allows users to determine the most significant stretch of continguous markers with high values of the FP statistic, which we interpret as a genomic region of high IBD.
It is interesting to note unrestricted values of F predicted by our disease models. For example, for γ = 5, the recessive model predicts 0.316 in cases and -0.025 in controls. This explains why our FP test has higher power than conventional tests for recessive traits: Inbreeding coefficients tend to be strongly positive in cases and only slightly negative in controls (this implies strong deviations from HWE in cases). Our parameter restrictions disallow negative F values, which reduces "noise" in the determination of significance levels. On the other hand, the dominant model predicts unrestricted F values of -0.041 in cases and 0.002 in controls. These values are only slightly different from 0 and will become non-negative in the FP test, that is, the FP essentially defaults to the genotype test with only a slight advantage over it.
The null distribution of the FP test statistic (under F ≥ 0) is not known and would be difficult to obtain, particularly for large numbers of markers whose test results are nonindependent.

Conclusion
At least for the recessive and dominant models considered here, our FP test is more powerful than allele and genotype tests. Thus, it represents an attractive alternative to Power for recessive disease models these conventional tests. A potential disadvantage of the FP test might be that it requires permutation testing for an appropriate determination of p-values. However, permutation testing is one of the best approaches to correct for multiple testing in genome-wide association studies and is often carried out anyway, so the FP test does not represent an additional burden.

Methods
The restriction of F and p to non-negative values reduces the parameter space of genotype frequencies, q 1 and q 3 . Simple algebraic manipulation (not shown here) of equations (1) demonstrates that 0 ≤ p ≤ 1 is always satisfied and so is F ≤ 1. However, F ≥ 0 is only satisfied for 2p(1p) ≥ q 2 . Figure 3 shows the parameter space for the three SNP genotypes. Because of q 1 + q 3 ≤ 1, the unrestricted parameter space corresponds to the lower triangle in Figure 1. The requirement that F ≥ 0 translates into the inequality, The resulting restricted parameter space corresponds to the area marked "F > 0" in Figure 3, that is, the area between the convex solid line and the diagonal. The restricted parameter space is only 2/3 of the surface of the unrestricted parameter space.
For power calculations, we assume disease models with two alleles and three genotypes, dd, Dd, and DD. The respective penetrances are f 1 , f 2 , and f 3 , where we set f 2 = f 1 for recessive models and f 2 = f 3 for dominant models. The "strength" of a model is measured by the penetrance ratio, γ = f 3 /f 1 , which is very approximately equal to the odds ratio. Thus, we have three genetic parameters, p, f 1 , and γ, which predict trait prevalence as K = γf 1 p 2 + f 1 (1 -p 2 ) for recessive traits and as K = q 2 f 1 + (1-q 2 )γf 1 , q = 1 -p, for dominant traits.
Fixing K = 0.05 and Q = 0.10 leaves one free parameter, which we vary to generate power curves. For a fixed set of genetic parameter values, Bayesian calculations yield conditional genotype frequencies in case and control individuals, from which random samples (replicates) are drawn.
For each of dominant and recessive models, with a value of the penetrance ratio, γ = 1, critical limits for test statistics are chosen so as to make the type I error for each of them equal to 0.05. That is, the critical limits are chosen such that the proportion of randomization samples exceeding this limit is equal to 0.05 (in other words, we q q q q q Parameter space for SNP genotype frequencies Figure 3 Parameter space for SNP genotype frequencies.
Parameters q 1 and q 3 are frequencies for genotypes AA and BB, respectively; F = probability of autozygosity, inbreeding coefficient.  Results of the genotype, FP, and allele tests applied to three published datasets: AMD in Caucasians [11], AMD HK in data from Hong Kong [15], and PD [16]. Estimated F values based on observed numbers of genotypes in case and control individuals are unrestricted, that is, may be negative.
are using the 95 th percentile of the computer-generated null distribution of the test statistic as the critical limit). Then power is determined for penetrance ratios ranging from 1 through 8. All power calculations were carried out based on 5,000 replicates.