Monte Carlo Error for Event Probabilities

Introduction

A simulation often answers a question about whether a particular event occurs: does an annual portfolio loss exceed a threshold, does a claim payment reach a policy limit, or does a simulated reserve fall below zero? Each simulation replication either satisfies the event condition or does not. That binary outcome makes an event probability a special case of estimating a mean.

The key consequence is practical: the observed proportion of event occurrences estimates the event probability, and its Monte Carlo standard error depends on both that proportion and the number of independent replications. The standard error measures simulation variability; it does not measure uncertainty about whether the model itself is appropriate.

Definition: Let $A$ be a specified event under the simulation model, with probability $p=P(A)$. For replication $j$, define the indicator $I_j=1$ if $A$ occurs and $I_j=0$ otherwise. If there are $N$ independent, identically distributed replications, the event count is $K=\sum_{j=1}^N I_j$, and the simulated event-probability estimate is $$ \widehat{p}=\frac{K}{N}=\frac{1}{N}\sum_{j=1}^N I_j. $$

This setup treats one complete replication as one trial. For example, if a replication simulates an entire year's claims for a portfolio, the trial is the year's portfolio outcome, not each individual claim within that year. The trials must follow the same simulation model and be independent for the variance formula below to apply directly.

Deriving the Variance and Standard Error

For a fixed event $A$, each indicator $I_j$ has a Bernoulli distribution: it equals 1 with probability $p$ and 0 with probability $1-p$. Its expectation is $p$. Since $I_j^2=I_j$, its variance is $$ \operatorname{Var}(I_j)=E[I_j^2]-(E[I_j])^2=p-p^2=p(1-p). $$

The estimator averages $N$ independent indicators. The variance of a sum of independent random variables is the sum of their variances; dividing that sum by $N^2$ gives the variance of the average:

E[\widehat{p}]=p,\qquad \operatorname{Var}(\widehat{p}) =\operatorname{Var}\left(\frac{1}{N}\sum_{j=1}^N I_j\right) =\frac{1}{N^2}\sum_{j=1}^N p(1-p) =\frac{p(1-p)}{N}.

Thus the standard deviation of the estimate over repeated runs of the same simulation is $\sqrt{p(1-p)/N}$. It is the theoretical Monte Carlo standard error, but $p$ is the quantity being estimated and is ordinarily unknown. Replacing it with the observed proportion produces the usual plug-in estimate:

Result: Under independent, identically distributed Bernoulli replications, the estimated event probability is $\widehat{p}=K/N$, its estimated Monte Carlo variance is $\widehat{p}(1-\widehat{p})/N$, and its plug-in standard error is $$ \widehat{\operatorname{SE}}(\widehat{p}) =\sqrt{\frac{\widehat{p}(1-\widehat{p})}{N}} =\sqrt{\frac{K(N-K)}{N^3}}. $$

The formula describes the typical random variation of $\widehat{p}$ across repeated simulations of size $N$. It is not a guaranteed maximum error, and it does not say that the unknown probability lies within one standard error of the estimate. Constructing an interval for the probability is a separate inferential task.

Worked Example: Exceeding an Annual Loss Threshold

A hypothetical insurer simulates $N=4{,}000$ independent annual portfolio outcomes under one fixed model. In each replication, the event is that aggregate annual loss exceeds $ \$250{,}000$. The event occurs in $K=164$ replications. Estimate the event probability and its plug-in Monte Carlo standard error.

Method. Record one indicator for each simulated year, then use the event count as the sum of those indicators. The method assumes the 4,000 simulated annual outcomes are independent and generated under the same model.

Step 1: Estimate the probability. The number of non-event replications is $4{,}000-164=3{,}836$, and $$ \widehat{p}=\frac{164}{4{,}000}=0.041. $$

Step 2: Estimate the variance. $$ \widehat{\operatorname{Var}}(\widehat{p}) =\frac{0.041(1-0.041)}{4{,}000} =\frac{0.041(0.959)}{4{,}000} =\frac{0.039319}{4{,}000} =0.00000982975. $$

Step 3: Take the square root. $$ \widehat{\operatorname{SE}}(\widehat{p}) =\sqrt{0.00000982975} \approx 0.0031. $$

The estimated probability is $0.0410$, or 4.10%, and the estimated standard error is 0.0031, or about 0.31 percentage points. Probabilities and standard errors are rounded to four decimal places. Check: $164/4{,}000$ is 4.1%, consistent with the reported estimate. The common mistake avoided here is using the number of events, 164, as the probability instead of dividing by all 4,000 replications.

Interpreting the Calculation Correctly

The variance formula is specific to the way the probability is estimated. Each replication contributes one zero-or-one result to the estimate. The fact that an underlying model may contain many policies, claims, or time periods does not change the formula if the final output for each independent replication is just whether the target event occurred.

The binomial model also clarifies what belongs in the denominator. If $N$ is the number of independent simulated portfolios and $K$ is the number whose annual loss breaches a threshold, use $K/N$. Do not instead divide by the number of policies in those portfolios. Conversely, if each simulated policy is genuinely an independent trial for the event being studied, then the trial count is the number of simulated policies.

An alternative calculation uses the sample variance of the indicators. With the usual sample-variance denominator $N-1$, that variance is $s_I^2=N\widehat{p}(1-\widehat{p})/(N-1)$, so $s_I/\sqrt{N}=\sqrt{\widehat{p}(1-\widehat{p})/(N-1)}$. The plug-in formula above substitutes $\widehat{p}$ into the Bernoulli variance and divides by $N$; it is the standard plug-in Monte Carlo standard error for an event proportion. For large $N$, the two estimates are close, but they are not algebraically identical.

Worked Example: A Relatively Infrequent Claim Event

A hypothetical actuary simulates $1{,}600$ independent policy-years using the same model. The event is that the claim payment in a policy-year exceeds $ \$80{,}000$. This event occurs in 24 replications. Find the estimated probability, plug-in variance, and plug-in standard error.

Method. Treat each policy-year as one Bernoulli trial because the event is defined for a policy-year, and the problem stipulates independent, identically distributed trials.

Step 1: Calculate the event proportion. There are $1{,}600-24=1{,}576$ non-event trials, so $$ \widehat{p}=\frac{24}{1{,}600}=0.015. $$

Step 2: Calculate the estimated variance. $$ \widehat{\operatorname{Var}}(\widehat{p}) =\frac{0.015(1-0.015)}{1{,}600} =\frac{0.015(0.985)}{1{,}600} =\frac{0.014775}{1{,}600} =0.000009234375. $$

Step 3: Calculate the estimated standard error. $$ \widehat{\operatorname{SE}}(\widehat{p}) =\sqrt{0.000009234375} \approx 0.0030. $$

The estimated probability is 0.0150 and its plug-in standard error is 0.0030, both reported to four decimal places. Check: $24/1{,}600=1.5\%$, and the standard error is smaller than the estimate but is not negligible relative to it. The common mistake avoided is treating a small observed percentage as if its Monte Carlo error must also be small in absolute terms.

Worked Example: Defining a Combined Portfolio Event

For each of $5{,}000$ independent simulated years of a hypothetical portfolio, an actuary records whether either (i) aggregate retained loss exceeds $ \$900{,}000$, or (ii) the year has at least 12 claims. The combined event occurs in 137 years. Estimate the probability and its plug-in standard error.

Method. Define one indicator for the combined event in each year. A year meeting both conditions still contributes just one event indicator, not two.

Step 1: Determine the event proportion. $$ \widehat{p}=\frac{137}{5{,}000}=0.0274. $$

Step 2: Calculate the estimated variance. $$ \widehat{\operatorname{Var}}(\widehat{p}) =\frac{0.0274(1-0.0274)}{5{,}000} =\frac{0.0274(0.9726)}{5{,}000} =\frac{0.02664924}{5{,}000} =0.000005329848. $$

Step 3: Take the square root. $$ \widehat{\operatorname{SE}}(\widehat{p}) =\sqrt{0.000005329848} \approx 0.0023. $$

The estimate is 0.0274 and the standard error is 0.0023, to four decimal places. Check: $137/5{,}000$ is 2.74%, matching the estimate. The common mistake avoided is counting years meeting both criteria twice: the target is the probability of the union of the conditions, and each year has one combined-event indicator.

What the Standard Error Does—and Does Not—Measure

Monte Carlo error arises because a finite set of random replications is used in place of the full probability calculation. If the model and event definition are held fixed, a different set of random draws will generally produce a different event count and estimate. The standard error summarizes the estimated scale of that variation.

This is distinct from model uncertainty. For instance, the simulation may use an assumed claim-frequency distribution, severity distribution, or dependence structure. The Bernoulli standard error does not assess whether those assumptions are suitable, nor does it account for uncertainty in parameters supplied to the simulation unless that uncertainty is explicitly represented in the simulation design. It quantifies sampling variation from the simulated trials under the specified setup.

The independence condition matters. If simulated years are dependent, the variance of their average includes covariance terms as well as individual indicator variances. In that case, $p(1-p)/N$ need not be the correct variance. Likewise, if trials have unequal weights or are generated by a different design, the simple event-proportion formula may not apply. State what counts as one replication and verify the design before using the calculation.

Worked Example: No Events Observed

A hypothetical catastrophe model is run for $N=600$ independent simulated years. The event is that annual loss exceeds $ \$5{,}000{,}000$, and no simulated year meets the condition. Calculate the plug-in estimate and standard error, then explain what the zero result means.

Method. Apply the event-proportion formula to the observed count, then distinguish the formula's output from a claim that the event is impossible.

Step 1: Calculate the event proportion. Here $K=0$, so $$ \widehat{p}=\frac{0}{600}=0. $$

Step 2: Calculate the plug-in standard error. $$ \widehat{\operatorname{SE}}(\widehat{p}) =\sqrt{\frac{0(1-0)}{600}}=0. $$

The plug-in standard error is zero because the plug-in estimate of $p$ is zero. This is a boundary effect of the formula, not proof that the true event probability is zero. To see why zero occurrences can happen even when the probability is positive, suppose solely for illustration that the true probability were $0.005$. The probability of observing no events in 600 independent trials would be $$ (1-0.005)^{600}=0.995^{600}\approx 0.0494. $$

Thus, even under that illustrative positive probability, a zero count would occur about 4.94% of the time. The actual probability is unknown, but the example shows why a zero plug-in standard error must not be interpreted as certainty. Check: the reported estimate and plug-in standard error follow directly from $K=0$; the separate calculation uses the stated illustrative probability, not an estimate from the observed zero count. The common mistake avoided is concluding that an event cannot occur simply because it did not occur in a finite simulation.

Exam Traps and a Reliable Workflow

A careful solution keeps the event definition, trial unit, count, and denominator aligned. Before calculating, identify the event as a precise condition and decide what constitutes one independent replication. Then check that $K$ counts replications satisfying that condition and that $0\leq K\leq N$.

Define one trial.
Specify the simulated unit, such as one policy-year or one complete portfolio-year, and confirm that the trials are identically distributed and independent.

Define the event and count.
Apply the stated event condition once per trial. Record $K$, the number of trials satisfying it, and $N$, the total number of trials.

Calculate the estimate and error measure.
Use $\widehat{p}=K/N$ and $\sqrt{\widehat{p}(1-\widehat{p})/N}$. Keep the variance and standard error distinct: the latter is the square root of the former.

Interpret within the design.
Report the standard error as an estimate of simulation variability for this event proportion. Do not treat it as a model-validation measure or as a guaranteed error bound.

A frequent incorrect calculation is $\sqrt{\widehat{p}/N}$, which omits the factor $1-\widehat{p}$. For the 164 events in 4,000 trials above, this incorrect calculation gives $\sqrt{0.041/4{,}000}\approx0.0032$, whereas the correct plug-in standard error is $\sqrt{0.041(0.959)/4{,}000}\approx0.0031$. The omitted factor is part of the Bernoulli variance and must be retained. The difference can be substantial when $\widehat{p}$ is not close to zero.

Key takeaway: A simulated event probability is a sample mean of zero-one indicators. Under independent, identically distributed replications, its variance is $p(1-p)/N$, and its usual plug-in standard error is $\sqrt{\widehat{p}(1-\widehat{p})/N}$. Correct interpretation depends on defining the event and replication unit consistently.

Check Your Understanding

For each problem, use the independent, identically distributed replication assumption stated in the question. Report probabilities and standard errors to four decimal places.

A risk model generates 2,000 independent portfolio-years, and 90 have an annual loss above the stated threshold. Find the estimated event probability and plug-in standard error. (Answer: $\widehat{p}=0.0450$; standard error $=0.0046$.)
In 1,200 independent simulated policy-years, 18 have a claim payment above the policy's stated attachment amount. Calculate the event-probability estimate and plug-in standard error. (Answer: $\widehat{p}=0.0150$; standard error $=0.0035$.)
A simulation of 3,500 independent years records 77 years in which the reserve is below zero at year-end. Find the estimate and plug-in standard error. (Answer: $\widehat{p}=0.0220$; standard error $=0.0025$.)
In 800 independent simulated seasons, no season has more than 20 claims. Calculate the plug-in estimate and standard error, and state the important qualification on interpreting the result. (Answer: $\widehat{p}=0.0000$; plug-in standard error $=0.0000$; neither value proves the event is impossible.)
A hypothetical portfolio simulation has 5,000 independent years. In 125 years, either the aggregate loss exceeds its threshold or the claim count is at least 15. A year satisfying both conditions is counted once. Calculate the estimated probability and plug-in standard error. (Answer: $\widehat{p}=0.0250$; standard error $=0.0022$.)

What You'll Learn