Calibrating Fine-Mapping With Simulations

Why Calibrate Fine-Mapping With Simulations?

A credible set is constructed from posterior probabilities, but its nominal posterior content does not by itself tell us how often the set will contain the causal variant across repeated studies. That repeated-study behavior depends on the data-generating process: sample size, LD, causal effects, residual variation, and the priors and model used for analysis. Simulation makes these ingredients explicit and lets us evaluate a fine-mapping procedure when the causal truth is known.

The goal here is not to claim that simulated data reproduce every real genomic study. It is to create controlled scenarios, run the same analysis that would be used on data, and record performance against known truth. Earlier material on credible sets and posterior inclusion probabilities supplies the output quantities; the essential extra ingredient in a simulation is a trustworthy record of which variants were made causal.

Coverage and size answer different questions. Coverage asks whether a set captures the truth; size describes how many variants remain in that set. A procedure can achieve high coverage by returning very large sets, so neither quantity should be interpreted alone. These are simulation performance measures, not statements that an individual real-data credible set has a known probability of containing the causal variant.

Specify the Data-Generating Process

For a quantitative phenotype, a transparent starting model is \(y=X\beta+\varepsilon\). Here \(X\) is an \(n\)-by-\(p\) genotype matrix, \(\beta\) is the vector of effects in phenotype units per genotype unit, and \(\varepsilon\) is residual variation. If genotype columns are centered and scaled to sample variance one, then each coefficient is in phenotype units per genotype standard deviation. This does not mean the phenotype itself has standard deviation one.

A simple way to impose LD is to generate rows of \(X\) from a multivariate normal distribution with mean zero and covariance matrix \(R\), where \(R_{jk}\) is the intended correlation between variants \(j\) and \(k\). After generation, center and scale each column if the analysis uses standardized genotypes. This is a convenient approximation, not a literal model for discrete genotype dosages. A more realistic design can generate haplotypes or genotype dosages from specified allele frequencies and haplotype frequencies, then calculate the resulting LD and check that it is close to the target.

The matrix \(R\) must be symmetric, have diagonal entries equal to one when it represents correlations, and be positive semidefinite. Otherwise it cannot be a valid covariance matrix. The causal configuration and effect sizes must also be declared before generating outcomes. For example, choose \(K\) causal variants and set \(\beta_j=0\) for all other variants. If effects are sampled from a distribution, record the distribution and the random seed so that the simulation can be reproduced.

The analysis must be specified as carefully as the generator. Record the fine-mapping method, the prior on the number of causal variants, the effect-size prior, the variants included, and any preprocessing. In a calibration experiment, these choices are fixed across replicates unless a design factor is deliberately being compared. The analysis prior need not be identical to the mechanism used to generate effects, but any intentional difference should be stated rather than left implicit.

Run Replicates and Score the Output

A replicate is created by generating a new genotype sample, new residuals, or both, according to the design. Apply the same analysis pipeline to each dataset. Preserve variant identities and ordering throughout: the LD matrix, genotype columns, summary statistics, and causal-truth record must all refer to the same variants in the same order. For a single signal, let \(C_r\) be the reported set in replicate \(r\), and let \(c_r\) be its true causal variant. Record \(I(c_r\in C_r)\) and \(|C_r|\).

For multiple signals, define how estimated sets will be associated with true signals before looking at the results. One clean design uses distinct genomic blocks, with one simulated causal variant in each block and one corresponding reported signal per block. When signals are not separable into blocks, use a prespecified matching rule based on signal identity or the analysis model’s components. Do not choose the matching after inspecting which rule makes coverage look best. If the method returns no set for a true signal, count that signal as uncovered and record the failure; silently dropping it would make the reported coverage optimistic.

For a nominal 95% set in a single-causal-variant analysis, one common construction is to rank variants by posterior probability and include variants until their posterior probabilities sum to at least 0.95. The procedure must match the model’s definition of a signal and its posterior output; do not combine probabilities across signals as though they formed one single-causal posterior when they do not. The 95% label describes posterior mass under the fitted model, while simulated coverage is assessed separately by checking whether the known causal variant is included.

Summarize Repeated-Data Performance

Across \(B\) replicates for a single signal per replicate, the simulated coverage proportion is the average of the coverage indicators. The mean set size is the average of the set sizes. These summaries should be reported together and within clearly stated scenarios, such as sample size and LD strength. If each replicate contains several signals, define whether the target is coverage per signal or per dataset; for per-signal coverage, each signal contributes one indicator and the denominator is the total number of simulated signals.

The coverage proportion is a direct summary of the simulated outcomes. To see what it targets, let \(I_r\) be the coverage indicator in replicate \(r\), and suppose replicates are generated independently under the same design. If the probability of coverage under that design is \(q=E(I_r)\), then linearity of expectation gives

Thus the average indicator targets the coverage probability for the specified simulation design. It does not establish performance under every possible LD pattern, effect configuration, or population. The number of replicates controls how much random Monte Carlo fluctuation appears in the result; a small run can look unusually favorable or unfavorable by chance.

Design Checks, Comparisons, and Failure Modes

Useful calibration starts with a small set of prespecified scenarios, varying one factor at a time where possible. For example, hold causal count, effect sizes, sample size, and priors fixed while changing LD; then hold LD fixed while changing sample size. A factorial design can study interactions, but its scenarios and comparisons should be chosen in advance. Reusing the same generated datasets to compare two analysis procedures can make the comparison less noisy because both methods face the same simulated truth and sampling conditions.

There are several reasonable choices for generating genotypes. A multivariate-normal approximation is convenient and allows direct control of \(R\); haplotype-based simulation can better represent discrete genotypes and allele frequencies. Resampling observed genotypes preserves features of a particular reference sample, but the results then depend on that sample and may not span the LD structures of interest. None is universally best: choose the generator to match the scientific question, and state its limitations.

• Check LD validity and alignment. A non-positive-semidefinite target matrix cannot generate the stated covariance. A variant-order mismatch between \(X\) and \(R\) can yield plausible-looking output for the wrong LD structure.
• Check scaling and units. If genotype columns are standardized, state that coefficients are per genotype standard deviation. Do not call them phenotype-standardized effects unless \(y\) was also standardized.
• Check the truth record. Store causal variant identifiers and effects at generation time. Reconstructing truth from the analysis output would make coverage circular.
• Check failures and denominators. Report analysis failures and missing sets. Excluding failed or difficult replicates without a prespecified rule can inflate coverage or shrink mean set size.
• Check what the design supports. Results from one LD matrix or one effect configuration are not universal calibration claims. Include relevant alternatives and describe the scenarios actually evaluated.

Simulation complements, rather than replaces, theoretical posterior interpretation and analysis of real data. A simulation can reveal systematic undercoverage, excessive set size, or sensitivity to a chosen design, but only for the generator and conditions it represents. Comparing parametric genotype generation with haplotype-based generation, or comparing methods on the same replicates, can show whether a conclusion depends on a modeling choice. It cannot by itself prove that the generator captures all features of a target population.