Statistical power for detecting epistasis QTL effects under the F-2 design

Epistasis refers to gene interaction effect involving two or more genes. Statistical methods for mapping quantitative trait loci (QTL) with epistasis effects have become available recently. However, little is known about the statistical power and sample size requirements for mapping epistatic QTL using genetic markers. In this study, we developed analytical formulae to calculate the statistical power and sample requirement for detecting each epistasis effect under the F-2 design based on crossing inbred lines. Assuming two unlinked interactive QTL and the same absolute value for all epistasis effects, the heritability of additive × additive (a × a) effect is twice as large as that of additive × dominance (a × d) or dominance × additive (d × a) effect, and is four times as large as that of dominance × dominance (d × d) effect. Consequently, among the four types of epistasis effects involving two loci, 'a × a' effect is the easiest to detect whereas 'd × d' effect is the most difficult to detect. The statistical power for detecting 'a × a' effect is similar to that for detecting dominance effect of a single QTL. The sample size requirements for detecting 'a × d', 'd × a' and 'd × d' are highly sensitive to increased distance between the markers and the interacting QTLs. Therefore, using dense marker coverage is critical to detecting those effects.


INTRODUCTION
Epistasis refers to gene interaction effect involving two or more genes. Evidence from studies in several species, including cattle [18,21], dogs [1], rat [23], drosophila [8] and humans [10,20], indicates that epistasis can play a significant role in both quantitative and qualitative characters. Epistasis effects of quantitative trait loci (QTL) have been found in soybean [17], maize [9] and tomato [11]. Among the many models to study epistasis effects, the linear partition of genotypic values into additive, dominance, and epistasis effects by Fisher [12] is considered a classical model for epistasis [6]. Cockerham [5] and Kempthorne [15] further partitioned Fisher's epistasis effects into four components, additive × additive, additive × dominance, dominance × additive, and dominance × dominance epistasis effects with the genetic interpretation of allele × allele, allele × genotype, genotype × allele, and genotype × genotype interactions respectively. The partition of Fisher's epistasis effect by Cockerham and Kempthorne provides a necessary tool to understand the precise nature of gene interactions. The genetic modeling of epistasis by Fisher, Cockerham and Kempthorne can be applied for testing epistasis effects of candidate genes and for mapping interactive QTL using genetic markers. Several statistical methods for mapping interactive QTL have been reported, e.g., the ANOVA method [22], the randomization test [3], the mixture model likelihood analysis based on Cockerham's orthogonal contrast [14], and the Bayesian approach for an outbred population [24]. However, little is known about the statistical power and sample size requirement for detecting each epistasis effect. The purpose of this article is to derive and analyze the statistical power and sample size requirements of the F-2 design for detecting epistasis effects. Formulae for statistical power take into account major factors affecting statistical power and sample size, including separate testing and estimation of additive × additive, additive × dominance, dominance × additive, and dominance × dominance effects, various levels of epistasis effects, marker-QTL distances, type-I error and sample size; formulae for sample size requirements take into account various levels of epistasis effects, marker-QTL distances, type-I and type-II errors.

Assumptions
Throughout this paper, two unlinked quantitative trait loci (QTL) on different chromosomes, QTL 1 and QTL 2, are assumed. Let A and B denote codominant marker loci linked to QTL 1 and QTL 2, respectively, and let θ 1 , and θ 2 denote the recombination frequencies between marker A and QTL 1, and between marker B and QTL 2, respectively. Then, the marker-QTL orders are A-θ 1 -QTL1 and B-θ 2 -QTL2. Two inbred lines with gene fixation for the markers and QTL are assumed for the convenience of analytical derivations. We assume Line 1 has AAQ 1 Q 1 BBQ 2 Q 2 genotype and Line 2 has aaq 1 q 1 bbq 2 q 2 genotypes, where A, a, B and b are marker alleles, and Q 1 , q 1 , Q 2 and q 2 are QTL alleles. The cross between these two lines yields the F-1 generation with AQ 1 /aq 1 BQ 2 /bq 2 individuals. The F-2 design is a result of matings among the F-1 individuals. The least squares partitioning of genotypic values and variances [12,15] will be used to derive the common genetic modeling for QTL values, assuming Hardy-Weinberg equilibrium. The statistical testing of epistasis effects using genetic markers will use a least squares model, because analytical solutions are available from this method. From this least squares model, elements required for the calculation of statistical power and sample size requirements will be derived, including the marker contrast for testing each epistasis effect and the variance of the contrast. Theoretical results for experimental designs will be compared with simulation studies assuming various levels of epistasis effects, and marker-QTL distance.

Genetic modeling and marker contrasts
The purpose of genetic modeling of the QTL genotypic values and individual phenotypic values is to establish a theoretical foundation for defining marker contrasts for testing epistasis effects. The least squares partitioning of genotypic values by Kempthrone [15] will be used for the genetic modeling. Let g i jkl = genotypic value of individuals with genotype ij at locus 1 and kl at locus 2, (i = Q 1 and j = q 1 of locus 1, k = Q 2 and l = q 2 of locus 2). Then, using Kempthrone's partitioning of genotypic values for the case of two unlinked loci [15,19], the genotypic value can be modeled as: where µ is the population mean of QTL genotypic values, α i , α j , α k , α l are the additive effects of QTL allele Q 1 , q 1 , Q 2 , q 2 , respectively; δ i j , δ kl are the dominance effects of locus 1 and locus 2, respectively; αα ik , αα il , αα jk , αα jl are the additive × additive effects accounting for the dependency of the effect of an allelic substitution at one locus on the allele present at a second locus; αδ ikl , αδ jkl are the additive × dominance effects accounting for the interaction of single alleles at locus 1 with the genotype at locus 2; δα i jk , δα i jl are the dominance × additive effects representing the interaction of the genotype at locus 1 with single alleles at locus 2; and δδ i jkl is the dominance × dominance effect representing the interaction between the genotype at locus 1 and the genotype at locus 2. In equation (1), α i j = α i + α j , α kl = α k + α l , αα i jkl = αα ik + αα il + αα jk + αα jl , αδ i jkl = αδ ikl + αδ jkl , δα i jkl = δα i jk + δα i jl . For an F-2 population with equal allele frequencies, it can be shown that the genetic effects in equation (1) have the following symmetry property: This symmetrical property leads to simplified modeling of equation (1), as shown in Table I. More importantly, this symmetry property will greatly simplify the marker contrasts for testing epistasis effects, allowing simple analytical solutions for evaluating statistical power and sample size requirement, as to be shown later. By combining the nine equations in Table I and solving for µ, a 1 , a 2 , d 1 , d 2 , i aa , i ad , i da , and i dd , the unique solutions of the effect parameters in terms of the genotypic values are: µ = 1 16 (g iikk + 2g iikl + g iill + 2g i jkk + 4g i jkl + 2g i jll + g j jkk + 2g j jkl + g j jll ) (2) In equations (2-10), a 1 = additive effect of QTL 1, d 1 = dominance effect of QTL 1, a 2 = additive effect of QTL 2, d 2 = dominance effect of QTL 2, i aa = additive × additive epistasis effect, i ad = additive × dominance epistasis Table I. Representation of genotypic values in terms of additive, dominance, and epistasis contributions for the case of two loci with two equally frequent alleles.
Genotype Frequency Genotypic Value effect, i da = dominance × additive epistasis effect, and i dd = dominance × dominance epistasis effect between QTL 1 and QTL 2. Equations (2-10) are foundations for marker contrasts for QTL detection under the F-2 design, and can be used for testing candidate genes in an F-2 design. When a QTL genotypic value is to be predicted by linked markers, as is the case in QTL detection, the QTL genotypic value can be modeled as: (11) where m i jkl is the effect of markers, and r i jkl is the genotypic residual value due to recombination between the markers and QTL. Note that the common genetic mean (µ) term in equation (11) is dropped for convenience of derivations. The two marker models with or without the µ term are equivalent models [13,19] that achieve the same result for statistical testing. In matrix notation, the QTL genotypic value modeled by genetic markers can be expressed as: where g is the column vector of QTL genotypic values, X is the design matrix for the marker effects, r is the recombination residual of the QTL value not explained by the common mean and the markers, and m is the vector of marker effects, i.e., The normal equations for equation (12) in matrix notation are X Xm = X g, and the solution to this normal equation is m = (X X) −1 X g, where (X X) −1 is the inverse of X X, and X is the transpose of X.
A phenotypic value is modeled as the summation of a QTL genotypic value (g i jkl ) and a random residual (e) with N(0,σ 2 e ) distribution, i.e., with ε i jkl = phenotypic residual value not explained by the marker effects due to the recombination and random residuals. Using matrix notation, equation (13) can be expressed as: The normal equations are X Xm = X y, and the estimator of m is given bŷ (15), then the four marker contrasts for testing epistasis effects are: where L aa , L ad , L da and L dd are the contrasts for testing additive × additive, additive × dominance, dominance × additive, and dominance × dominance effects, respectively.

Variances of recombination and phenotypic residuals
Following the approach of Bulmer [2] (Eq. (5.1) on page 58), Table I can be expressed more succinctly as: where z 1 is the number of Q 1 alleles in a particular individual at the first putative QTL, and z 2 is the number of Q 2 alleles in a particular individual at the second putative QTL, (z 1 , z 2 = 0, 1, or 2 respectively). Equation (20) provides a convenient model for deriving variance and covariance of the QTL genotypic values (App. A). Let σ 2 g = the total QTL genotypic variance, Derivations for equations (21,22) are given in Appendix A. The total genetic variance is partitioned into eight independent components, and each variance component is a function of one effect only. This property greatly facilitates the evaluation of the contribution of an effect to the total genetic variance. Note that equivalent partitions can be obtained under alternative models, but they have different meanings in interpreting gene effects, different structures of variance components, and different properties in statistical estimation that may affect the study of QTL [14]. The population variance of recombination residuals for equation (12) is Applying equation (23) to the F-2 design, and utilizing equation (22), the recombination residual variance of QTL genotypic values is found to be: Derivation of equation (24) is given in Appendix B. The residual variance of phenotypic values for F-2 design under equation (14) is:

Mathematical formulae for statistical power and sample size
Statistical power (π) is the probability that an effect is detected when the effect is present, commonly denoted by π = 1 − β, where β is the type II error, i.e., the probability of false 'negatives'. A standardized normal distribution denoted by N(0,1) is assumed for deriving the statistical power. The normal distribution is chosen because the calculation of the exact residual degrees of freedom is unnecessary, providing analytical simplicity. Since the residual degrees of freedom are sufficiently large for the sample sizes discussed in this article (N ≥ 200), the normal distribution practically yields identical results as the t-distribution that is often used in QTL analysis. The general expression for π is: where Z is a N(0,1) normal variable, z i is the ordinate of the standardized normal curve corresponding to the type II error of β, and Φ is the cumulative distribution function of standard normal random variable. The application of equation (25) to QTL mapping designs requires two key elements: a marker contrast for detecting each epistasis effect, and the variance of the contrast. Let c be a contrast vector that defines an estimable function of m, then E(c m) = c m. Based on this result and the m vector in Appendix B, the mathematical expectation of each contrast given by equations (16)(17)(18)(19), denoted by E(L i ), i = aa, ad, da, dd, are functions of markers-QTL recombination frequencies and the QTL effects being tested, i.e., Using E(L i ) in place of L i for i = aa, ad, da, dd as defined by equations (16)(17)(18)(19), the z i value in equation (25) can be expressed as: where N i is the sample size and V i = N i var(L i ), for i = aa, ad, da, dd. For convenience, V i will be referred to as the 'kernel' of the contrast variance, meaning that V i differs from var(L i ) only by a constant of N i . Let For convenience, we will refer to h 2 a 1 and h 2 a 2 as the additive heritabilities, h 2 d 1 and h 2 d 2 as the dominance heritabilities, and h 2 aa , h 2 ad , h 2 da and h 2 dd as the additive × additive, additive × dominance, dominance × additive and dominance × dominance heritabilities, respectively. Then, the expressions of V i in terms of QTL parameters are given as follows: The derivations of equations (31-34) are given in Appendix B. Letting λ i = E(L i )/ √ V i , then z i in equation (25) can be expressed in terms of QTL parameters as: where In equations (36-39), θ 1 , θ 2 are the marker-QTL recombination frequencies.
Theoretical predictions of statistical power for various parameters using equation (25) and equations (35-39) are shown in Figures 1-4. The implications of these results will be discussed along with the results of simulation studies. Using the above results, the minimum sample size required for given levels of type I and type II errors can be expressed as: where Z α/2 and Z β are the ordinate of the standardized normal curve corresponding to the probabilities of α/2 and β. The sample size given by equation (40) is an increasing function of marker-QTL recombination frequencies, as well as type-I and type-II errors, and a decreasing function of heritability. Sample size requirements obtained from equation (40) for two type-I errors corresponding to the "suggestive" and "significant" linkages proposed by Lander and Kruglyak [16], and different levels of epistasis heritabilities and marker-QTL recombination frequencies are given in Table II, assuming a 95% statistical power.

Simulation studies on statistical power
Simulation studies were conducted using the Monte Carlo method to evaluate the theoretical results on statistical power for detecting epistasis effects under the F-2 designs. Markers and QTL genotypes were generated such that the true recombination frequencies and each QTL effect used to generate these genotypes can be obtained reversely from the data. A total of 100 sets of marker-QTL genotypes were generated. The phenotypic value of each individual is obtained as the summation of the individual QTL genotypic value and a random residual with N(0,1) distribution. For each set of genotypic data, 10 000 replicates were generated for Figures 1-4. Two interactive QTL without linkage and all epistasis effects have the same absolute value are assumed, thus the heritability of additive × additive (a × a) effect is twice as large as that of additive × dominance (a × d) or dominance × additive (d × a) effect, and is four times as large as that of dominance × dominance (d × d) effect. Heritabilities of additive × additive effect are used in the range of 0.025 to 0.25. Sample sizes of 200-2000 individuals resulting from crossing between inbred lines were generated. The significant levels (type I errors) used were those corresponding to "suggestive linkage" and "significant linkage" proposed by Lander and Kruglyak [16] with type-I errors of 0.0034 and 0.00072 respectively. Table II. Sample size required to achieve 95% power with a type I error of 5% for the F-2 design. Heritability Since the powers for a × d effect and d × a effect as expected identical in the theoretical derivations and almost identical in the simulations, we only show the results for a × d effect, and the results for d × a effect are not included. Table III shows the observed statistical power for epistasis effects for different sample sizes and heritability levels. Statistical powers observed from the simulated data agreed well with the predicted powers as shown in Figures 1-4. Among the four types of epistasis effects involving two loci, 'a × a' effect is the easiest to detect whereas 'd × d' effect is the most difficult to detect. The statistical power for detecting 'a × a' effect is similar to that for detecting dominance effect of a single QTL. The power is poor for detecting 'a × d' or 'd × a' effect and is extremely poor for detecting 'd × d' effect. This trend is consistent across a range of sample sizes (Fig. 1), heritabilities (Fig. 2), and Table III. Statistical power for detecting epistasis effects based on simulated data (θ 1 = θ 2 = 0.10).
Sample size 400 400 400 1000 1000 1000 recombination frequencies (Fig. 3). The general relationship between power and type-I error is shown in Figure 4 for a wide range of type-I errors. The sample size requirements for detecting 'a × d', 'd × a' and 'd × d' are highly sensitive to increased distance between the markers and the interacting QTLs. Therefore, using dense marker coverage is critical to detecting those effects.

DISCUSSION
Epistasis, a potentially important genetic component underlying complex traits, has not been extensively explored in QTL analysis [4]. The results obtained in this study provide some guidelines regarding the statistical power and sample size requirement for detecting epistasis effects under the F-2 design. In general, detecting epistasis effect is more difficult than detecting single QTL effect, except for additive × additive effect, which has about the same power as dominance effect. Detecting epistasis effects involving dominance effect is considerably more challenging than detecting single QTL effect. This difficulty could be reduced to some extent by decreasing marker spacing, because the statistical power for detecting epistasis effects involving dominance is highly sensitive to increased marker spacing (Tabs. I and II, Fig. 3). The statistical power and sample size requirements in this study assume the use of a single marker for detecting epistasis effects. For statistical methods using flanking    markers to detect epistasis effects, results of statistical power could be somewhat overestimates and sample size requirements somewhat underestimates.
Extending results in this study to interval mapping is straightforward theoretically but necessarily will require a lengthy development. Since the increase in power of interval mapping over a single marker analysis is only slight [7], results obtained in this study can be considered as close approximation to statistical power and sample size requirements, with statistical power being slightly underestimated and sample size slightly overestimated than those under an interval mapping.

APPENDIX A: PARTITIONING OF GENOTYPIC VALUES AND VARIANCES UNDER THE F-2 DESIGN
Under the F-2 design, the partition of QTL genotypic values is shown in Table I. The QTL alleles are assumed to have equal allele frequency and the two loci are assumed unlinked.

A.2. Analysis of variances
Since all the expectations of effect coefficients are zero, and the function of z 1 and function of z 2 are independent, we have Since z 1 is a random variable taking the values 0, 1, and 2 with probabilities 1/4, 1/2, and 1/4. Thus (z 1 − 1)(1 − 4z 1 + 2z 2 1 ) takes values −1, 0, and 1 with these probabilities, therefore, Similarly, Therefore, the covariance of any two different effect coefficients is zero. Also,