Checking Null Calibration by Simulation

Why Check Calibration by Simulation?

An association test is calibrated when, under its null hypothesis and the assumptions of its model, its reported significance level corresponds to the intended probability of rejecting a true null hypothesis. A test that claims a 5% type I error rate should reject about 5% of the time across repeated datasets generated under that null. In genetic association studies, this can be checked by repeatedly simulating or resampling phenotypes while holding relevant features of the study design fixed.

A simulation check is a diagnostic of a specified analysis under a specified data-generating setup. It can reveal that an implementation or modeling choice produces too many or too few small p-values, or a test-statistic distribution that differs from its reference distribution. It does not establish that the simulated setup matches every feature of the real study, nor does it prove calibration for all genetic variants or all possible alternatives.

Specify the Null Simulation

Start by stating what is null and what is held fixed. For example, a single-variant test may have a null hypothesis of no association between genotype dosage and phenotype after adjustment for specified covariates. A useful simulation can hold the observed genotypes and covariates fixed, generate a phenotype under that null model, and rerun the complete planned association analysis. This asks a conditional question: does the analysis behave as intended for these genotypes and covariates, under the chosen phenotype model?

The simulated data must reflect the features needed for the null distribution to be meaningful. These may include covariate effects, phenotype type, family or relatedness structure, ancestry or other population structure, and the test’s variance or link model. If the null model is composite—for example, it includes unknown nuisance parameters—those parameters need to be specified or estimated in a way appropriate to the simulation. The check then assesses the procedure under those choices; it is not automatically an exact guarantee for every possible nuisance-parameter value.

Two common approaches are parametric simulation and resampling. In parametric simulation, draw new phenotypes from an explicit null model, such as a normal linear model or a binary-outcome model with no genotype effect but with specified covariate effects. In resampling, rearrange observed outcomes or suitable residuals to create null replicates. Resampling requires a valid exchangeability argument: under the null, the observations being rearranged must be interchangeable under the design and conditioning assumptions. Unrestricted shuffling is generally not justified merely because the tested genetic effect is null.

For a genome-wide analysis, the test statistics at different variants are usually dependent because genotypes can be correlated through linkage disequilibrium and population structure. Simulations should preserve the genotype data when the goal is to assess the full analysis pipeline. When estimating uncertainty in a genome-wide summary, account for dependence among variants; treating every variant-by-replicate result as independent can make uncertainty appear too small.

Quantifying Monte Carlo Uncertainty

A simulation estimate is itself subject to random variation. If \(B\) independent null replicates are generated and \(R\) of them reject at the chosen threshold, the estimated rejection rate is \(\hat{r}=R/B\). If each replicate has the same true rejection probability \(r\), then \(R\) has a binomial distribution and $$ \operatorname{SE}(\hat{r})=\sqrt{\frac{r(1-r)}{B}}. $$ Because \(r\) is generally unknown, a useful approximate standard error is \(\sqrt{\hat{r}(1-\hat{r})/B}\). When comparing with a specified nominal rate \(\alpha\), the null simulation’s expected Monte Carlo standard error is \(\sqrt{\alpha(1-\alpha)/B}\). These quantities help distinguish a substantial discrepancy from an ordinary fluctuation in a finite simulation.

The binomial calculation assumes independent simulation replicates with the same rejection probability. It applies, for example, to independent generated phenotypes on a fixed design when each replicate is analyzed once. If replicates share random elements, or if the reported rate pools many correlated tests from the same simulated dataset, the simple binomial standard error may not be appropriate. In that case, define the independent simulation unit clearly and assess variability at that level.

Increasing \(B\) reduces Monte Carlo uncertainty at a rate proportional to \(1/\sqrt{B}\). It does not correct a misspecified null model, an invalid resampling scheme, or a flawed test. A very precise estimate of the wrong simulation target remains a precise answer to the wrong question.

Check Both Rejection Rates and Distributions

A rejection-rate check focuses on one or more thresholds. A distribution check uses more of the simulated output: compare the simulated test statistics with the distribution assumed by the test, or compare simulated p-values with their null expectation. For a valid test based on a continuous, correctly specified null distribution, p-values are uniformly distributed on the interval from 0 to 1 under the null. In particular, the probability of a p-value at or below \(\alpha\) is \(\alpha\). With discrete tests or estimated nuisance parameters, this uniformity may be approximate or conservative rather than exact.

For a one-degree-of-freedom score or Wald test, for example, the reference statistic may be a chi-square random variable with one degree of freedom under the null. Compare the simulated distribution with that reference across its center and tails. A test may have a broadly similar center but an excess of extreme values, or it may be shifted throughout the distribution. Looking at more than one summary can help describe the discrepancy, although the pattern alone does not uniquely identify its cause.

Resampling Must Respect the Design

A permutation-based null check is credible only when the chosen permutations preserve the structure needed under the null. If phenotype distributions differ by a covariate, unrestricted permutation of phenotypes can erase that covariate relationship. If observations are clustered or related, shuffling individuals independently can break their dependence. The resulting replicates may therefore represent a different null world from the one the test is intended to handle.

A Practical Calibration Workflow

Interpretation and Common Pitfalls

A simulation can help separate a calibration problem in the analysis from a signal pattern in observed data. If a procedure shows excess rejection under a well-justified null simulation, its observed small p-values deserve caution. Conversely, a rejection rate near nominal in one simulation is limited evidence: it supports calibration for that simulated setup and sample size, not universal validity or correct modeling of the actual study.

• Too few replicates. A small simulation can show an apparently high or low rejection rate by chance. Report \(B\) and quantify Monte Carlo uncertainty rather than treating \(\hat{r}\) as exact.
• Invalid exchangeability. Shuffling outcomes without preserving covariate or dependence structure can simulate the wrong null. State why the resampling units are interchangeable.
• Changing the target unintentionally. Simulating genotypes as well as phenotypes answers a different question from holding the observed genetic design fixed. Say which features are conditioned on and which are regenerated.
• Reading one summary as a diagnosis. An excess rejection rate identifies a discrepancy under the simulated setup, not its cause. Check the null generator, implementation, and model assumptions before assigning an explanation.
• Confusing calibration with power. A null simulation assesses false-positive behavior under the chosen null; it does not measure the probability of detecting a true genetic effect.