Simulation Error and Number of Trials, introduction to probability lesson

A simulation of 10,000 trials reports that an event occurred in 2484 of them, an estimate of $0.2484$ . How far from the true probability might that be, and would four times as many trials make it four times as accurate? The estimate is an average of zeros and ones, so the results about sample means answer both questions. This lesson derives the standard error of a simulation estimate, reads two actual seeded runs against a known answer, shows how to estimate or bound the standard error when the probability is unknown, and separates the random error of a finite run from mistakes in the model or the program.

Key terms

Simulation estimate

Definition Let $A$ be an event in a probability model with $P(A) = p$ . A simulation runs $n$ trials of the model and records the indicator $I_i$ , equal to $1$ if $A$ occurs in trial $i$ and $0$ otherwise. With $S_n = I_1 + \dots + I_n$ , the number of trials in which $A$ occurred, the simulation estimate of $p$ is

\hat{p}_n = \frac{S_n}{n} = \frac{1}{n} \sum_{i=1}^{n} I_i .

Each indicator is $\operatorname{Bernoulli}(p)$ , so linearity gives $E[\hat{p}_n] = p$ for every $n \ge 1$ , whether or not the indicators are independent. If they are independent, then $S_n \sim \operatorname{Bin}(n, p)$ and

\operatorname{Var}(\hat{p}_n) = \frac{p(1 - p)}{n}.

In plain words A simulation estimate is a sample mean of zeros and ones. It is unbiased, so it is right on average over possible runs, but one run gives one value, which generally differs from $p$ .

Example As a benchmark with a known answer, take an event with $p = 1/4$ . One indicator has mean $1/4$ and variance $\tfrac14 \cdot \tfrac34 = \tfrac{3}{16}$ . For this lesson the benchmark was actually simulated in Python 3.12 with the generator random.Random(20261012), counting a trial as a success when the uniform draw from random() is below $0.25$ . The recorded run had 30 successes in the first 100 trials, 249 in the first 1000 and 2484 in all 10,000, so the estimates were $0.30$ , $0.249$ and $0.2484$ .

Monte Carlo error

Definition The Monte Carlo error of a simulation estimate is the difference $\hat{p}_n - p$ that comes from using finitely many simulated trials. For independent indicators its standard deviation, called the standard error of $\hat{p}_n$ , is

\operatorname{SE}(\hat{p}_n) = \sqrt{\frac{p(1 - p)}{n}} .

Multiplying $n$ by $4$ halves the standard error, and dividing the standard error by $10$ takes $100$ times as many trials. A standard error is not a largest possible error. A rigorous bound on the probability of a large error comes from Chebyshev’s inequality, which gives, for every $\varepsilon > 0$ ,

\begin{aligned} &P(\lvert \hat{p}_n - p \rvert \ge \varepsilon) \\[0.4em] &\quad \le \frac{p(1 - p)}{n\varepsilon^2} . \end{aligned}

In plain words The standard error measures the typical size of the Monte Carlo error across possible runs with the same number of trials. Precision grows with the square root of the number of trials, so cutting the standard error by a factor of $10$ costs a hundred times as many trials. A particular run can land closer to $p$ than one standard error or farther away.

Example For $p = 1/4$ and $n = 10000$ ,

\begin{aligned} \operatorname{SE}(\hat{p}_n) &= \sqrt{\frac{3/16}{10000}} \\[0.4em] &= \frac{\sqrt{3}}{400} \approx 0.00433, \end{aligned}

about $0.43$ percentage points. With 40,000 trials it would be $\sqrt{3}/800 \approx 0.00217$ . Chebyshev’s inequality with $\varepsilon = 0.02$ bounds the probability of an error of at least $0.02$ by

\begin{aligned} &\frac{3/16}{10000 \cdot 0.0004} \\[0.4em] &\quad = \frac{3}{64} = 0.046875 . \end{aligned}

In the recorded run the errors were $0.05$ after 100 trials, $-0.001$ after 1000 and $-0.0016$ after 10,000. So the run was farther from $p$ after 10,000 trials than after 1000. That does not contradict the shrinking standard error, which describes the spread over possible runs of each length and does not promise that one run moves closer to $p$ at every step.

Plug-in standard error

Definition When $p$ is unknown, the plug-in standard error replaces $p$ by the estimate,

\widehat{\operatorname{SE}} = \sqrt{\frac{\hat{p}_n(1 - \hat{p}_n)}{n}} .

It is computed from the run, so it is itself random, and it is an estimate of the standard error, not the standard error itself. It equals $0$ whenever $S_n = 0$ or $S_n = n$ , whatever the true $p$ .

In plain words Plugging in the estimate is a practical way to report the size of the Monte Carlo error when the true probability is the very thing being estimated. It is a sensible approximation when the run contains many successes and many failures, and it can mislead for rare events, where a short run may contain no successes at all.

Example For the recorded run, $\hat{p}_{10000} = 0.2484$ gives $\widehat{\operatorname{SE}} \approx 0.00432$ , close to the true $0.00433$ . Now take a rare event with $p = 0.001$ and a run of 1000 independent trials. The chance of zero successes is $0.999^{1000} \approx 0.368$ , and such a run reports the estimate $0$ with plug-in standard error $0$ . Neither number shows that the event is impossible or that the estimate is exact, since many positive probabilities make zero successes quite likely.

Conservative standard error

Definition For every $p$ in $[0, 1]$ the identity $\tfrac14 - p(1 - p) = (p - \tfrac12)^2$ shows that $p(1 - p) \le \tfrac14$ , with equality exactly at $p = \tfrac12$ . So for independent indicators

\operatorname{SE}(\hat{p}_n) \le \frac{1}{2\sqrt{n}} ,

and the conservative standard error $1/(2\sqrt{n})$ bounds the true standard error without any knowledge of $p$ .

In plain words The worst case for the variance of an indicator is a fair coin. Planning with the worst case guarantees that the actual standard error is no larger, which makes it useful for choosing the number of trials before a run. For probabilities near $0$ or $1$ the bound can be far above the real standard error.

Example With $n = 10000$ the conservative standard error is $0.005$ , a little above the true $0.00433$ for $p = 1/4$ . A standard error of at most $0.01$ for every $p$ needs $1/(2\sqrt{n}) \le 0.01$ , which holds exactly when $n \ge 2500$ . For the rare event with $p = 0.001$ and $n = 1000$ , the true standard error is about $0.0010$ while the conservative value is about $0.0158$ , valid but very loose.

Replication

Definition A replication is a complete new run of the whole simulation experiment, with the same model, event and number of trials, whose random inputs are modeled as independent of those of the other runs. The estimates from several replications show the run-to-run variation directly. Continuing one run adds trials to the same estimate. Running the same program again with the same generator and seed reproduces the same sequence, so such a rerun is a reproducibility check, not a replication.

In plain words A fresh replication gives new information about how much the estimate varies. An exact rerun should return identical output, which confirms the computation but says nothing new about Monte Carlo error. Picking whichever seed gives the most pleasing result after seeing the output hides exactly the variability that replications are meant to show.

Example A second run of the same algorithm with seed 20261013 gave 2510 successes in 10,000 trials, an estimate of $0.2510$ . The two estimates $0.2484$ and $0.2510$ lie $0.0016$ below and $0.0010$ above $1/4$ , both within one standard error. Running the first program again with seed 20261012 gives 2484 successes once more, which checks reproducibility and adds no new estimate. Treating the streams from different seeds as independent is a working assumption about the generator.

Model error and implementation error

Definition Model error is the difference between the probability of the event in the simulated model and its probability in the real process that the model is meant to describe. Implementation error is a mismatch between the program and the stated model or event, such as counting the wrong event, using a wrong parameter or boundary, or making trials dependent where the model treats them as independent. If a program produces independent and identically distributed trials, its estimate converges in probability to the success probability of the trials it actually produces, so neither kind of error shrinks as $n$ grows.

In plain words More trials shrink the typical Monte Carlo error and leave the other two untouched. A precise estimate of the wrong quantity is still wrong. Comparing a program with an exact benchmark catches many implementation errors, while checking the model against the real process is a separate task that the simulation cannot do for itself.

Example Suppose a real event has probability $0.25$ but the model uses $0.2$ . With a million trials the standard error is $\sqrt{0.2 \cdot 0.8/10^6} = 0.0004$ , and Chebyshev’s inequality puts the estimate within $0.002$ of $0.2$ with probability at least $0.96$ . The gap of $0.05$ from the real value, about $125$ standard errors, remains. For an implementation error, take the event of two heads in two fair independent tosses, with probability $1/4$ . A program that draws one random bit and uses it for both tosses counts a success whenever that bit is heads, so it estimates $1/2$ , and more trials make that wrong answer more precise.

Check yourself

A simulation estimates an event with true probability $1/2$ from 2500 independent trials. Find the standard error of the estimate and the number of trials that would halve it.

Show the answer

The standard error is

\sqrt{\frac{0.25}{2500}} = \frac{0.5}{50} = 0.01,

one percentage point. Halving it takes four times as many trials, so 10,000 trials, with standard error $0.5/100 = 0.005$ .

A simulation of a rare event runs 400 independent trials and records no successes. Find the estimate and its plug-in standard error. Does the run show that the event has probability $0$ ? Find the probability of zero successes if the true probability were $0.005$ , and the conservative standard error for 400 trials.

Show the answer

The estimate is $0/400 = 0$ , and the plug-in standard error is $\sqrt{0 \cdot 1/400} = 0$ . The run does not show that the probability is $0$ . If the true probability were $0.005$ , zero successes would have probability $0.995^{400} \approx 0.135$ , so this outcome is not unusual for a positive probability, and the true standard error would be about $0.0035$ , not $0$ . The conservative standard error is $1/(2\sqrt{400}) = 1/40 = 0.025$ . It is a valid upper bound on the true standard error for every $p$ , but it describes the spread of the estimate, not a guaranteed limit on its error, and for a rare event it is very loose.

A colleague reruns the seed 20261012 program five times, gets 2484 successes every time, and concludes that the estimate has no Monte Carlo error. The colleague also suggests running a million trials to fix a model whose event probability is known to be wrong. Explain what is wrong with each conclusion.

Show the answer

The reruns use the same generator and seed, so they reproduce the same sequence of draws. Identical output confirms reproducibility, but the five reruns are one run repeated, not five replications. The Monte Carlo error is still described by the standard error $\sqrt{3}/400 \approx 0.00433$ , and a genuine second run with seed 20261013 gave the different estimate $0.2510$ . A million trials would shrink the standard error, but the estimate would converge to the probability in the wrong model, so the model error would remain. The model itself has to be corrected.

Simulation Error and Number of Trials

Key terms

Simulation estimate

Monte Carlo error

Plug-in standard error

Conservative standard error

Replication

Model error and implementation error

Check yourself

How this lesson connects

Builds on

Used later in