Stratified Simulation and Neyman Allocation

Why Stratify a Simulation?

A simulation can be more precise when its output is divided into meaningful groups before sampling. For example, a hypothetical insurer estimating average claim cost might distinguish claim types with different frequencies or severities. Rather than letting random sampling determine how many simulations fall in each group, stratified simulation fixes a sample size in each group and combines the results using the groups’ known population proportions.

This design has two related questions: how to calculate the uncertainty of the combined estimate, and how to choose the number of simulations in each stratum. The variance formula answers the first. Neyman allocation answers the second when the total number of simulations is fixed and each simulation has the same cost.

Definition: Partition the simulation population into $H$ mutually exclusive and exhaustive strata. Let $W_h$ be the known proportion of the population in stratum $h$, with $\sum_{h=1}^H W_h=1$. Let $Y$ be the output of interest, and define its conditional mean and variance in stratum $h$ as $\mu_h=E[Y\mid h]$ and $\sigma_h^2=\operatorname{Var}(Y\mid h)$. If $n_h$ independent simulations are drawn within stratum $h$, let $\bar{Y}_h$ be their sample mean. The stratified estimator of the overall mean is $\hat{\mu}_{\mathrm{str}}=\sum_{h=1}^H W_h\bar{Y}_h$.

The stratum weights are population proportions, not the proportions of simulations assigned to the strata. The estimator remains weighted by $W_h$ even when the simulation allocation $n_h$ is deliberately unequal. Each within-stratum sample should represent that stratum, and samples from different strata are assumed independent.

Variance of the Stratified Estimator

With fixed sample sizes and independent simulations within and across strata, the variance is

\operatorname{Var}(\hat{\mu}_{\mathrm{str}}) =\sum_{h=1}^H W_h^2\frac{\sigma_h^2}{n_h}.

The square on each weight matters: the estimator multiplies the random sample mean $\bar{Y}_h$ by $W_h$, and variance scales by the square of a constant. Because the stratum samples are independent, there are no covariance terms between their sample means.

Write the weighted estimator.
$\hat{\mu}_{\mathrm{str}}=\sum_h W_h\bar{Y}_h$, where the weights sum to one.

Apply variance rules.
Independence gives $\operatorname{Var}(\sum_h W_h\bar{Y}_h)=\sum_h W_h^2\operatorname{Var}(\bar{Y}_h)$.

Use the variance of a sample mean.
For $n_h$ independent draws with variance $\sigma_h^2$, $\operatorname{Var}(\bar{Y}_h)=\sigma_h^2/n_h$, yielding the stated formula.

This also shows why the estimator is unbiased: $E[\hat{\mu}_{\mathrm{str}}]=\sum_h W_h E[\bar{Y}_h]=\sum_h W_h\mu_h=E[Y]$. If the conditional variances are not known, estimate the variance of the estimator by replacing each $\sigma_h^2$ with the within-stratum sample variance $s_h^2$:

\widehat{\operatorname{Var}}(\hat{\mu}_{\mathrm{str}}) =\sum_{h=1}^H W_h^2\frac{s_h^2}{n_h}, \qquad \widehat{\operatorname{SE}}(\hat{\mu}_{\mathrm{str}}) =\sqrt{\widehat{\operatorname{Var}}(\hat{\mu}_{\mathrm{str}})}.

These formulas apply to independent simulation draws within each stratum. For sampling without replacement from a finite stratum, a finite-population correction is needed; the formula above is not the exact variance in that setting. Also, a stratum with $n_h=1$ has no within-stratum sample variance estimate, so a pilot or additional simulations are needed to estimate its variability.

Neyman Allocation

Suppose the total number of simulations $n=\sum_h n_h$ is fixed, every simulation has equal cost, and the stratum standard deviations $\sigma_h$ are known or adequately estimated. The allocation minimizing the variance is proportional to the product of stratum weight and standard deviation:

n_h =n\frac{W_h\sigma_h}{\sum_{j=1}^H W_j\sigma_j}.

The weight reflects how much the stratum contributes to the target overall mean; the standard deviation reflects how variable its simulated outcomes are. Both affect how useful another draw in that stratum will be.

To justify the result, minimize $\sum_h W_h^2\sigma_h^2/n_h$ subject to $\sum_h n_h=n$. The Lagrangian is $\mathcal{L}=\sum_h W_h^2\sigma_h^2/n_h+\lambda(\sum_h n_h-n)$. Setting its derivative with respect to $n_h$ equal to zero gives $-W_h^2\sigma_h^2/n_h^2+\lambda=0$, so $n_h=W_h\sigma_h/\sqrt{\lambda}$. Summing these allocations and using $\sum_h n_h=n$ gives $\sqrt{\lambda}=\sum_j W_j\sigma_j/n$, which yields the rule. Substitution gives the minimum variance $(\sum_h W_h\sigma_h)^2/n$.

The formula usually produces fractional sample sizes. In practice, round to integer counts whose sum is exactly $n$, adjusting the largest fractional remainders or comparing nearby feasible allocations. If every stratum must be represented, require at least one simulation per stratum and allocate the remaining budget according to the rule. An allocation based on pilot estimates of $\sigma_h$ is only as reliable as those estimates.

Worked Example: Allocating Simulations for Claim Cost

A hypothetical claim-cost model has three strata with population weights $0.50$, $0.30$, and $0.20$. Conditional mean costs are $12$, $20$, and $8$ thousand dollars, and conditional standard deviations are $2$, $4$, and $1$ thousand dollars. Allocate 120 independent simulations using Neyman allocation and find the resulting standard error, treating the standard deviations as known.

Method. First calculate the target mean: $0.50(12)+0.30(20)+0.20(8)=6+6+1.6=13.6$ thousand dollars. For allocation, the products $W_h\sigma_h$ are $0.50(2)=1.0$, $0.30(4)=1.2$, and $0.20(1)=0.2$. Their sum is $2.4$.

The continuous allocations are $120(1.0/2.4)=50$, $120(1.2/2.4)=60$, and $120(0.2/2.4)=10$. These are already integers and sum to 120.

The estimator variance, in squared thousands of dollars, is $0.50^2(2^2)/50+0.30^2(4^2)/60+0.20^2(1^2)/10 =1/50+1.44/60+0.04/10=0.020+0.024+0.004=0.048$. Thus the standard error is $\sqrt{0.048}=0.2191$ thousand dollars, or $\$219.09$, rounded to the nearest cent. The calculation checks because each allocation is proportional to its $W_h\sigma_h$, and the counts total 120. A common mistake is to allocate in proportion to $W_h$ alone; that would ignore the much greater variability in the second stratum.

Worked Example: Estimating Precision From Stratum Samples

A hypothetical portfolio is divided into three known groups with weights $0.40$, $0.35$, and $0.25$. Independent simulations produce sample sizes $25$, $35$, and $20$, sample mean costs of $10$, $14$, and $8$ thousand dollars, and sample standard deviations of $5$, $8$, and $4$ thousand dollars. Calculate the stratified estimate and its estimated standard error.

Method. Weight the sample means by population proportions, not by sample counts: $\hat{\mu}_{\mathrm{str}}=0.40(10)+0.35(14)+0.25(8)=4+4.9+2=10.9$ thousand dollars.

Estimate the variance using the squared weights: $0.40^2(5^2)/25+0.35^2(8^2)/35+0.25^2(4^2)/20 =0.16(25)/25+0.1225(64)/35+0.0625(16)/20 =0.1600+0.2240+0.0500=0.4340$ in squared thousands of dollars. Therefore the estimated standard error is $\sqrt{0.4340}=0.6588$ thousand dollars, or $\$658.79$, rounded to the nearest cent.

The result checks dimensionally: the estimate is in thousands of dollars, while the variance is in squared thousands and its square root returns to thousands. The common error avoided here is weighting each sample mean by $n_h/\sum_j n_j$; those are simulation proportions, not the portfolio’s stratum proportions.

Worked Example: Neyman Allocation for an Event Probability

A hypothetical system’s failure indicator is simulated in three strata with weights $0.80$, $0.15$, and $0.05$. The conditional failure probabilities are $0.01$, $0.20$, and $0.60$. Allocate 200 simulations approximately by Neyman allocation, estimate the variance of the resulting failure-probability estimator, and compare it with allocation proportional to stratum size.

Method. In each stratum the output is Bernoulli, so its standard deviation is $\sqrt{p_h(1-p_h)}$. The three standard deviations are $\sqrt{0.0099}=0.0995$, $\sqrt{0.16}=0.4000$, and $\sqrt{0.24}=0.4899$, to four decimals. The allocation products are approximately $0.80(0.0995)=0.0796$, $0.15(0.4000)=0.0600$, and $0.05(0.4899)=0.0245$; their sum is $0.1641$.

The continuous sample sizes are approximately $200(0.0796/0.1641)=97.1$, $200(0.0600/0.1641)=73.1$, and $200(0.0245/0.1641)=29.9$. Use integer counts $97,73,30$, which sum to 200. The variance is $0.80^2(0.0099)/97+0.15^2(0.16)/73+0.05^2(0.24)/30 =0.006336/97+0.0036/73+0.0006/30 \approx0.00006532+0.00004932+0.00002000=0.00013463$. The standard error is $\sqrt{0.00013463}=0.0116$, to four decimals.

Allocation proportional to stratum size gives $160,30,10$. Its variance is $0.006336/160+0.0036/30+0.0006/10=0.00003960+0.00012000+0.00006000=0.00021960$, with standard error $0.0148$. Thus the Neyman allocation is more precise for these assumed conditional probabilities. The total estimated failure probability, if the conditional rates are used as the stratum sample means, is $0.80(0.01)+0.15(0.20)+0.05(0.60)=0.008+0.030+0.030=0.0680$. The trap is to assign almost all simulations to the largest stratum; the smaller high-variability strata also matter to precision.

Practical Checks and Limitations

Stratification is useful when the strata are well defined, their population weights are known, and within-stratum variation can be reduced or managed through the design. Before reporting a result, check that the weights sum to one, the estimator uses those weights, the allocation counts sum to the simulation budget, and the variance uses squared weights. If an estimated conditional standard deviation is zero only because a small pilot observed identical values, it may not truly be zero; assigning no further simulations on that basis can be unsafe.

Neyman allocation is optimal under its stated assumptions, not a universal allocation rule. It assumes a fixed total budget, equal simulation costs, and suitable values for the stratum standard deviations. Its continuous solution must be translated into integer counts, and practical constraints may require minimum sample sizes. When strata have unequal simulation costs, the equal-cost allocation formula does not apply unchanged.

Key takeaway: Estimate an overall mean by weighting each stratum’s sample mean with its population weight. Its variance is $\sum_h W_h^2\sigma_h^2/n_h$; for a fixed equal-cost simulation budget, Neyman allocation assigns samples in proportion to $W_h\sigma_h$.

Check Your Understanding

For each problem, identify the relevant stratum weights and use them consistently in the estimator or allocation.

A hypothetical model has two strata with weights $0.25$ and $0.75$, conditional standard deviations $4$ and $1$, and a budget of 100 simulations. What continuous Neyman allocation does the rule give? (Answer: $57.14$ and $42.86$ simulations, approximately 57 and 43 after rounding.)
A stratified simulation has weights $0.20,0.50,0.30$, sample sizes $20,25,30$, and sample means $12,9,15$. What is the weighted estimate? (Answer: $11.4$.)
For the preceding strata, suppose the sample standard deviations are $3,2,4$, respectively. What are the estimated variance and standard error? (Answer: Variance $0.1060$; standard error $0.3256$.)
A Bernoulli event has conditional probabilities $0.04$ and $0.25$ in two equally weighted strata. For a budget of 100 simulations, what is the approximate continuous Neyman allocation? (Answer: Approximately $31.2$ and $68.8$ simulations.)

What You'll Learn