Validating Simulation Code With Invariants

Why Validate Simulation Code With Invariants?

A simulation can produce plausible-looking trial results even when its implementation is wrong. A treatment label may be assigned to the wrong outcome, a dropout may be counted twice, or an analysis may use the number randomized where it needs the number observed. Repeating such a simulation many times can make its estimates precise without making them correct. Validation therefore needs checks that are tied to properties the code must satisfy, not just to whether its output looks reasonable.

An invariant is a property that must remain true under a specified model, input, or operation. Some invariants are exact: for example, the number randomized must equal the sum of the randomized treatment-group counts. Other checks are probabilistic: if a generator produces independent Bernoulli outcomes with probability \(p\), their simulated mean should vary around \(p\), within a tolerance that reflects Monte Carlo error. Keeping those categories distinct is central to sound validation.

Definition: An exact invariant is a relation that must hold for every valid execution under stated conditions, apart from explicitly allowed numerical rounding. A stochastic check compares a random simulation result with a known probability result and allows for the variation expected from finite simulation. A limiting-case check supplies inputs under which the model or algorithm has a simple boundary result that the implementation should reproduce.

These checks are not substitutes for carefully defining the simulated trial or its decision-relevant operating characteristics. They test whether the implementation obeys consequences of that specification. The introductory material on estimating operating characteristics distinguishes a trial-level outcome, such as success or failure, from the estimated probability of that outcome across replicates; the same distinction matters here. First check the data and calculations within a replicate, then check the behavior across replicates.

Exact Invariants: Follow the Counts and Quantities

Exact checks are often the easiest to interpret because they do not require a tolerance for random variation. They are useful wherever an algorithm routes participants, records statuses, transforms data, or calculates a summary. The check should be written from the meaning of the data structure, not inferred from whatever the code happens to return.

For example, suppose each randomized participant has exactly one final disposition: completed follow-up, discontinued treatment but remained in follow-up, or lost to follow-up. If these categories are exhaustive and mutually exclusive, then their counts must sum to the number randomized. The same is true within each treatment group. An invariant of this kind can expose duplicated records, missing status assignments, and denominators that were inadvertently changed during analysis.

Worked Example: Checking Participant Accounting

A hypothetical simulation generates a trial with 120 randomized participants in a 1:1 allocation. It reports 60 in each group. In the active group, 54 complete follow-up and 6 are lost to follow-up; in the control group, 51 complete follow-up and 9 are lost to follow-up. Check the accounting and calculate the total number completing follow-up.

The allocation check is \(N_{\text{active}}+N_{\text{control}}=N_{\text{randomized}}\). Substitution gives \(60+60=120\), so allocation counts reconcile. The follow-up check is that, within each group, randomized participants equal completers plus those lost to follow-up:

60=54+6,\qquad 60=51+9.

Both group-level checks pass. The total number completing follow-up is \(54+51=105\), and the total lost is \(6+9=15\). The overall check is \(105+15=120\), which also passes.

If a report instead claimed 106 completers while retaining the stated group counts, it would violate the invariant. That is a coding or reporting failure to investigate, not Monte Carlo noise: the counts are exact for this generated trial.

Conservation checks can also be written for quantities other than participants, when the simulation represents them explicitly. For example, if a record is split into mutually exclusive categories, the counts before and after the split must reconcile. But a quantity should be called conserved only if the model and data definitions require it. A total drug amount, exposure measure, or number at risk need not be conserved across time when the modeled processes include elimination, accumulation, entry, or censoring.

Limiting Cases as Executable Thought Experiments

A limiting case is a deliberately simple input that makes a model’s expected behavior transparent. It can test whether a component is connected correctly even when realistic inputs produce complicated output. The key is to state which output is expected to simplify and which parts of the pipeline remain meaningful.

Consider a simulated response with a treatment shift and normal residual error:

Y_i=\beta_0+\beta_T I_i+\varepsilon_i,\qquad \varepsilon_i\sim N(0,\sigma^2),

where \(I_i=1\) for active treatment and \(I_i=0\) for control. If \(\beta_T=0\), the conditional mean response is the same in both groups. This does not require the two finite simulated group means to be identical when \(\sigma>0\): random residuals still differ. If \(\sigma=0\), however, every control response must equal \(\beta_0\) and every active response must equal \(\beta_0+\beta_T\). Such an input checks deterministic assignments and arithmetic without residual variation.

Worked Example: Testing a Zero-Variability Case

Use the response model above with \(\beta_0=100\) mg/L, \(\beta_T=5\) mg/L, and \(\sigma=0\). Generate three participants per group. State the expected responses and the difference in group means.

With \(\sigma=0\), each residual is exactly zero. For control, \(I_i=0\), so \(Y_i=100+5(0)+0=100\) mg/L for all three participants. For active treatment, \(I_i=1\), so \(Y_i=100+5(1)+0=105\) mg/L for all three.

Thus \(\bar{x}_{\text{control}}=(100+100+100)/3=100\) mg/L and \(\bar{x}_{\text{active}}=(105+105+105)/3=105\) mg/L. The active-minus-control difference is \(105-100=5\) mg/L, matching \(\beta_T\). Any other generated response indicates an error in treatment coding, parameter use, or residual generation. An inferential standard error based on these zero-variance data may be undefined; this test is for raw generation and arithmetic, not for forcing every downstream inferential procedure to return a valid result.

Other useful boundary inputs include event probabilities \(p=0\) and \(p=1\), where Bernoulli outcomes should always be 0 and always be 1, respectively; a zero treatment effect, where the model’s treatment contribution should vanish; and a single participant, where group-level calculations requiring multiple observations may not be defined. Tests should encode these distinctions. An undefined result can be correct behavior if the input does not support the requested calculation; silently returning a plausible number is not necessarily a success.

Known Distributions as Analytical Oracles

Some simulation components have a known probability distribution under a controlled setup. That result serves as an analytical reference, or oracle, for a stochastic check. Unlike an exact invariant, the simulated estimate need not equal its theoretical target. The comparison must allow for sampling variation, and the assumed independence and distributional conditions must match the code being tested.

For example, if each of \(N\) independent participants has a binary event with probability \(p\), the event count \(K\) has a binomial distribution with mean \(Np\) and variance \(Np(1-p)\). The sample proportion \(\hat{p}=K/N\) therefore has mean \(p\) and variance \(p(1-p)/N\). These statements follow from writing \(K\) as the sum of \(N\) independent Bernoulli variables, each with mean \(p\) and variance \(p(1-p)\): means add, and variances add under independence.

E(\hat{p})=p,\qquad \operatorname{Var}(\hat{p})=\frac{p(1-p)}{N},\qquad SE(\hat{p})=\sqrt{\frac{p(1-p)}{N}}.

Here the standard error describes the standard deviation of the estimator over repeated simulations under the stated assumptions. It is not a guarantee that one run will fall within a particular number of standard errors. It can, however, guide a tolerance and help distinguish an implausibly large discrepancy from ordinary random variation.

Worked Example: Checking a Known Bernoulli Probability

A hypothetical code path generates \(N=2{,}500\) independent binary outcomes with event probability \(p=0.20\). In one run it records 540 events. Compare the observed proportion with the analytical result and assess whether the discrepancy is a strong warning.

The observed proportion is \(\hat{p}=540/2{,}500=0.216\). Under the specified Bernoulli model, the expected proportion is \(0.20\), and its standard error is

SE(\hat{p})=\sqrt{\frac{0.20(1-0.20)}{2{,}500}} =\sqrt{0.000064}=0.008.

The difference is \(0.216-0.20=0.016\), or \(0.016/0.008=2\) standard errors. This single run is not, by itself, compelling evidence of a coding error; such finite-run departures can occur under the model. It does not prove the generator is correct either. A systematic discrepancy across repeated controlled checks, or an error much larger than the variability expected under the assumptions, would warrant investigation. The arithmetic depends on independent outcomes with the specified common probability; if the code induces dependence or uses varying probabilities, this reference calculation is not appropriate.

For a stronger controlled check, compare the simulated result with an exact calculation where feasible, or repeat the check across independent runs and account for the resulting Monte Carlo uncertainty. Avoid turning a stochastic check into an exact assertion such as “the simulated proportion must equal 0.20.” That would reject valid implementations simply because random samples fluctuate.

Seeds, Reproducibility, and Diagnostic Boundaries

A fixed random-number seed makes a run reproducible when the software, random-number generator, code path, and relevant settings are held fixed. Reproducibility is valuable because it lets a developer reproduce a failing scenario and compare outputs before and after a code change. It is a diagnostic aid, not a mathematical guarantee that the generated data are correct. A generator can return the same incorrect sequence every time.

Worked Example: Reproducing a Failing Randomized Assignment

Suppose a hypothetical randomization function is tested with a specified seed, generator, and software version. On the first run it returns the treatment labels \((1,0,1,1,0,0)\), and a repeated run from the same initialized state returns the identical sequence. The intended allocation is three active and three control participants.

The reproducibility check passes because the output is identical under the same controlled setup. The allocation invariant also passes: \(1+0+1+1+0+0=3\) active assignments, with \(6-3=3\) control assignments. These checks establish neither that the sequence was generated by the intended randomization mechanism nor that the implementation is unbiased. For example, a function that always returns three active labels followed by three control labels could pass the count check but fail to randomize positions.

A useful next test checks the function against its specified mechanism, using a suitable analytical reference or additional targeted tests. Record the seed together with the generator and relevant software details so a failure can be reproduced. Do not require different seeds to produce different outcomes in every run: random sequences can occasionally coincide, and exact seed behavior can vary across generators or software versions.

When a check fails, locate the earliest stage at which the invariant or reference result stops holding. Validate participant generation before treatment assignment, assignment before outcome simulation, and generated outcomes before analysis summaries. This staged approach helps distinguish a data-generation error from an analysis error. Keep a failing seed and the smallest useful scenario as a regression test so a later change does not reintroduce the problem.

State the claim being tested.
Write down whether it is an exact equality, a limiting-case result, or a probabilistic comparison, and specify the assumptions that make it valid.

Choose controlled inputs.
Use small, interpretable scenarios where possible, such as zero variability or a boundary event probability, without requesting undefined downstream quantities.

Check each pipeline stage.
Apply accounting checks to generated records, then verify the outcome model and analysis separately so a failure can be localized.

Interpret departures appropriately.
Treat exact violations as implementation failures; interpret stochastic departures relative to their expected Monte Carlo variability and assumptions.

Record and preserve failures.
Save the controlled inputs, seed, generator details, and failing outputs, then rerun the test after the correction.

Common Mistakes and Complementary Checks

Calling a random property an invariant. “The event proportion equals \(p\)” is not an exact rule for a finite random sample. The correction is to compare it with \(p\) using a tolerance informed by its standard error.
Using a check that does not match the assumptions. The binomial reference requires independent Bernoulli outcomes with a common probability. If probabilities vary by participant or outcomes are correlated, use a reference appropriate to that setup.
Confusing reproducibility with validity. The same erroneous code can reproducibly generate the same erroneous output. Pair seed-based reruns with mathematical and mechanism-specific checks.
Testing only totals. Correct group counts do not establish correct random allocation; verify both the invariant and the intended assignment mechanism.
Expecting every analysis to succeed in a degenerate case. A zero-variance input can be useful for testing generated values, while a standard-error calculation may be undefined. Define expected behavior for each pipeline stage.

These checks complement, rather than replace, other validation approaches. Code review and unit tests can inspect implementation logic; analytical oracles can test controlled scenarios; and broader checks of simulated data and trial behavior can reveal problems not exposed by small cases. No finite set of invariants proves a complex simulation correct. The goal is to create independent checks that make common failures visible and interpretable before operating characteristics are used to support a development decision.

Key takeaway: Use exact invariants for properties that must always hold, limiting cases for simple boundary behavior, and known distributions for stochastic references with justified tolerances. Use fixed seeds to make failures repeatable, not to certify correctness.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions.

A trial generator assigns 40 participants to active treatment and 40 to control. It records 35 and 37 participants, respectively, as completing follow-up, with the remainder lost to follow-up. Does the accounting invariant pass? What is the total number completing follow-up?
For \(N=1{,}000\) independent Bernoulli outcomes with \(p=0.10\), calculate the expected event count and the standard error of the simulated proportion.
In the response model \(Y_i=\beta_0+\beta_T I_i+\varepsilon_i\), what raw responses should be generated when \(\beta_0=80\), \(\beta_T=0\), and \(\sigma=0\)?
A fixed-seed test returns the same output on every rerun. What does this establish, and what does it not establish?

Solutions. (1) The group checks pass because \(40=35+5\) and \(40=37+3\). The total completing follow-up is \(35+37=72\). (2) The expected count is \(Np=1{,}000(0.10)=100\). The standard error of the proportion is \(\sqrt{0.10(0.90)/1{,}000}=0.00949\), approximately 0.949 percentage points. (3) With zero residual variability and no treatment shift, every response is \(80\) in both groups. (4) It establishes reproducibility under the same controlled software and generator setup; it does not establish that the output follows the intended randomization or simulation mechanism.