Mixture Chi-Square Null Distribution for SKAT

Why the SKAT Null Is a Mixture

The SKAT statistic is quadratic in a vector of variant-level scores. Even when each score is approximately normal under the null hypothesis, the quadratic form is not generally a single chi-square variable. Its null distribution depends on how variable the scores are, how they covary, and how the variants are weighted. The eigenvalues of the resulting quadratic form summarize those features: each positive eigenvalue becomes a weight on an independent one-degree-of-freedom chi-square variable.

The derivation starts with the quadratic score statistic from the earlier material on deriving the SKAT quadratic score statistic. This tutorial focuses on its null distribution and on using that distribution to calibrate a test. The result is exact when the score vector is normal with known covariance; in many association models, normality is an asymptotic approximation.

Assumptions and the Eigenvalue Derivation

For the result, the null score must have mean zero and covariance \(\Sigma_U\), and \(W\) must be symmetric positive semidefinite. If \(U\sim N(0,\Sigma_U)\), these conditions give the distribution exactly, even if \(\Sigma_U\) is singular. If the score is not exactly normal, the same result is an asymptotic null approximation when a multivariate central limit theorem applies and the null covariance is consistently estimated. In association analyses, nuisance covariates must be accounted for in the score and its covariance; otherwise the stated null distribution may not apply.

To see the result, write \(U=\Sigma_U^{1/2}Z\), where \(Z\sim N(0,I)\). Substitution into the statistic gives \(Q=Z^\mathsf{T}\Sigma_U^{1/2}W\Sigma_U^{1/2}Z\). The matrix in this expression is symmetric and positive semidefinite. Orthogonally diagonalizing it rotates \(Z\) without changing its standard multivariate normal distribution. Its rotated coordinates remain independent standard normal variables, and their squares are independent \(\chi^2_1\) variables.

The eigenvalues in the definition use \(A=W^{1/2}\Sigma_UW^{1/2}\). This has the same nonzero eigenvalues as \(\Sigma_U^{1/2}W\Sigma_U^{1/2}\), so either matrix gives the mixture weights. In the common linear-model formulation, if \(G\) is the \(n\times m\) genotype matrix and \(R\) is a residual vector with covariance \(\sigma^2 M\), then \(U=G^\mathsf{T}R\) has covariance \(\sigma^2G^\mathsf{T}MG\). The score-space matrix is therefore \(W^{1/2}(\sigma^2G^\mathsf{T}MG)W^{1/2}\). Equivalently, the nonzero weights can be obtained from the observation-space matrix \(\sigma^2M^{1/2}GWG^\mathsf{T}M^{1/2}\). This explains the reference to kernel eigenvalues: they are the eigenvalues of the relevant weighted score or observation-space quadratic form, not necessarily the eigenvalues of an unweighted genotype matrix.

Zero eigenvalues contribute nothing to \(Q\); they do not add degrees of freedom. A positive eigenvalue \(\lambda_j\) scales one \(\chi^2_1\) component. Thus, the mixture is not generally determined just by the number of variants or by a nominal degrees-of-freedom count.

Worked Mixture Calculations

Moments, Calibration, and Numerical Evaluation

The mixture representation gives useful checks on a calculation. Because a \(\chi^2_1\) variable has mean 1 and variance 2, independence implies that the null mean and variance of \(Q\) are

For the first worked example, these are \(4+1=5\) and \(2(4^2+1^2)=34\). These moments can detect some mistakes—for example, an incorrect eigenvalue list may fail to reproduce the trace and sum of squared eigenvalues—but matching moments alone does not establish that a tail probability is correct.

For a general mixture, the upper-tail probability \(P(Q\geq q_{\mathrm{obs}}\mid H_0)\) is computed numerically. Methods such as the Davies and Imhof procedures evaluate the distribution of a weighted chi-square sum; they do not replace the mixture with an ordinary chi-square. Approximate moment-matching methods can be faster, but may be less accurate in extreme tails. For a genome-wide or other stringent threshold, tail accuracy matters: a small absolute error can be large relative to the true probability.

The weights must be nonnegative for this representation. A materially negative computed eigenvalue usually signals numerical error or an incorrectly constructed matrix, since the target matrix is positive semidefinite. Tiny negative values caused by floating-point rounding may be treated as zero only after checking the scale and numerical precision. If all eigenvalues are zero, then \(Q=0\) under the stated model and there is no nondegenerate test statistic to calibrate.

Diagnostics, Alternatives, and Limitations

A reliable calibration requires more than diagonalizing a matrix. Check that the score is centered under the null, that nuisance covariates were handled consistently in both the score and its covariance, that the same weights were used to form \(Q\) and \(A\), and that the eigenvalues are nonnegative up to numerical tolerance. Also compare the computed sum of eigenvalues with \(\operatorname{tr}(A)\), and compare the sum of their squares with \(\operatorname{tr}(A^2)\). Disagreement indicates an implementation or numerical problem.

The normal-score approximation can be weak when data are highly discrete, observations are influential, or the effective information is limited. In those settings, the mixture based on an asymptotic covariance may not give accurate tail probabilities. This is a limitation of the approximation, not a change in the eigenvalue derivation. Diagnostics should include null calibration in suitable simulations or other checks appropriate to the data and model.

A single scaled chi-square is justified when all positive eigenvalues are equal; otherwise it is only an approximation. A burden-score test instead tests a collapsed, directional score and has a different null statistic and alternative emphasis. It is not a substitute for evaluating the SKAT mixture when the SKAT statistic is the chosen test. Resampling or permutation can offer an alternative calibration when its exchangeability assumptions hold, but may be computationally demanding and must respect the study design and covariate structure.

A common mistake is to count variants and compare \(Q\) with a chi-square distribution having that many degrees of freedom. This reasoning assumes equal mixture weights, which need not hold even when all variants are included. The correction is to calculate the eigenvalues of \(W^{1/2}\Sigma_UW^{1/2}\) and evaluate the resulting weighted sum. Another mistake is to interpret a mixture p-value as evidence that a particular variant is associated: the statistic tests the set-level null, and its calibration does not identify which variant contributes.