Posterior Expectations As Monte Carlo Averages

From a Posterior Distribution to an Average

A posterior distribution can describe a trial parameter in detail, but some questions require a single numerical summary. What is the posterior mean response rate? What is the posterior probability that a treatment effect exceeds a clinically relevant threshold? These questions are expressed as posterior expectations. When an integral is inconvenient to evaluate directly, independent draws from the posterior provide a practical approximation.

Let $y$ denote the observed trial data and $\theta$ a parameter with posterior distribution $\pi(\theta\mid y)$. For a function $g$ whose posterior expectation exists, the quantity of interest is the posterior average of $g(\theta)$. For example, setting $g(\theta)=\theta$ gives the posterior mean of $\theta$; setting $g(\theta)=I(\theta>c)$, where $I$ is an indicator, gives the posterior probability that $\theta$ exceeds $c$.

Definition: The posterior expectation of $g(\theta)$, given data $y$, is $$ \mu_g = E[g(\theta)\mid y] = \int g(\theta)\,\pi(\theta\mid y)\,d\theta, $$ with a sum replacing the integral for a discrete parameter. If $\theta^{(1)},\ldots,\theta^{(M)}$ are independent draws from $\pi(\theta\mid y)$, the Monte Carlo estimator is $$ \widehat{\mu}_g = \frac{1}{M}\sum_{m=1}^{M}g(\theta^{(m)}). $$

The quantity being averaged is $g(\theta^{(m)})$, not necessarily the parameter itself. That distinction makes the same computation useful for means, probabilities, risks, contrasts, and other summaries. The central assumptions for the results below are that the draws are independent and identically distributed from the target posterior, and that $E[g(\theta)^2\mid y]$ is finite. The finite second moment ensures that the Monte Carlo variance is finite and supports the usual standard-error and large-sample approximations.

Why the Monte Carlo Average Works

Write $G_m=g(\theta^{(m)})$, and let $\sigma_g^2=\operatorname{Var}(g(\theta)\mid y)$. Because every draw has the posterior distribution, each $G_m$ has expectation $\mu_g$ and variance $\sigma_g^2$. Independence makes the variance of their sum the sum of their variances. Thus the average is unbiased, and its variance decreases as the number of draws increases.

E[\widehat{\mu}_g\mid y] = \frac{1}{M}\sum_{m=1}^{M}E[G_m\mid y] = \mu_g, \qquad \operatorname{Var}(\widehat{\mu}_g\mid y) = \frac{1}{M^2}\sum_{m=1}^{M}\operatorname{Var}(G_m\mid y) = \frac{\sigma_g^2}{M}.

This is the key reason the Monte Carlo average is useful: it targets the desired posterior expectation, and its simulation variance is known in terms of the variance of the sampled values. The Monte Carlo standard error (MCSE) is the standard deviation of the estimator over repeated sets of posterior draws, conditional on the observed data. If $\sigma_g$ were known, the MCSE would be $\sigma_g/\sqrt{M}$. In practice, estimate it using the sample standard deviation $s_g$ of the $M$ values $g(\theta^{(m)})$:

\widehat{\operatorname{MCSE}}(\widehat{\mu}_g) = \frac{s_g}{\sqrt{M}}, \qquad s_g^2 = \frac{1}{M-1}\sum_{m=1}^{M} \left(g(\theta^{(m)})-\widehat{\mu}_g\right)^2.

The law of large numbers says that the Monte Carlo average converges to $\mu_g$ as $M$ grows. When the variance is finite, the central limit theorem further says that, for sufficiently large $M$, its sampling distribution is approximately normal with mean $\mu_g$ and variance $\sigma_g^2/M$. This supports an approximate Monte Carlo error range such as $\widehat{\mu}_g \pm 1.96\,\widehat{\operatorname{MCSE}}$. It is an approximation to uncertainty from numerical sampling, not a posterior credible interval for $\theta$.

Worked Example: Estimating a Posterior Mean Response Rate

Suppose a hypothetical trial’s posterior distribution for the treatment-arm response rate is sampled independently five times. The five sampled response-rate values are $0.4460741,\ 0.4500000,\ 0.5000000,\ 0.5500000,\ 0.5539259$. Estimate the posterior mean response rate and its MCSE.

Here $g(\theta)=\theta$, so the estimator is the average of the sampled rates: $$ \widehat{\mu}_g = \frac{0.4460741+0.4500000+0.5000000+0.5500000+0.5539259}{5} = \frac{2.5000000}{5} = 0.5000. $$

The deviations from $0.5000$ are approximately $-0.0539259,\ -0.0500000,\ 0,\ 0.0500000,\ 0.0539259$. Their squared deviations sum to approximately $0.010816$, so $$ s_g = \sqrt{\frac{0.010816}{5-1}} = \sqrt{0.002704}=0.05200. $$ The estimated MCSE is therefore $$ \widehat{\operatorname{MCSE}} = \frac{0.05200}{\sqrt{5}} \approx 0.02326. $$

The simulation estimates the posterior mean response rate as 0.5000, or 50.00%, with estimated Monte Carlo error of 0.02326, or 2.326 percentage points. Five draws are too few for a dependable large-sample error approximation; the calculation illustrates the method, not a recommendation to stop after five draws.

Monte Carlo Averages for Posterior Probabilities

A posterior probability is itself an expectation. For a threshold $c$, use $g(\theta)=I(\theta>c)$, which equals 1 when the condition holds and 0 otherwise. Averaging these indicator values gives the proportion of posterior draws satisfying the condition. This is often a direct way to summarize whether a clinical threshold is met.

Worked Example: Probability That a Treatment Effect Exceeds a Threshold

In a hypothetical trial, let $\delta$ be the treatment benefit in percentage points, with larger values favoring treatment. The clinical threshold is 5 percentage points. Among 10 independent draws from the posterior distribution of $\delta$, 7 exceed 5. Estimate the posterior probability of exceeding the threshold and its MCSE.

Set $g(\delta)=I(\delta>5)$. There are seven values of 1 and three values of 0, so $$ \widehat{P}(\delta>5\mid y) = \frac{1}{10}\sum_{m=1}^{10}I(\delta^{(m)}>5) = \frac{7}{10} = 0.7000. $$ For indicator draws, the estimated variance of the sampled values is $\hat{p}(1-\hat{p})$, where $\hat{p}=0.7000$. The MCSE is $$ \widehat{\operatorname{MCSE}} = \sqrt{\frac{0.7000(1-0.7000)}{10}} = \sqrt{0.02100} \approx 0.1449. $$

The estimated posterior probability that the benefit exceeds 5 percentage points is 70.00%; its MCSE is about 14.49 percentage points. That large MCSE reflects the small number of draws. It does not mean there is a 70% probability that a trial result is correct, nor does it describe uncertainty in the true treatment effect in place of the posterior distribution.

How the Number of Draws Affects Error

For a fixed target $g$ and posterior distribution, the standard deviation $\sigma_g$ does not change when more draws are generated. The estimator’s MCSE is proportional to $1/\sqrt{M}$. Consequently, multiplying the number of independent draws by four halves the MCSE; multiplying it by 100 reduces the MCSE by a factor of 10. Halving an MCSE generally requires about four times as many independent draws, not twice as many.

Worked Example: Planning the Number of Independent Draws

Suppose the posterior probability of a clinically meaningful benefit is approximately $0.25$, and the aim is an MCSE near 0.025. For indicator values, the variance is $p(1-p)$. Using $p=0.25$ as a planning value, find the approximate number of independent draws required.

The MCSE is approximately $\sqrt{p(1-p)/M}$. Set this equal to the target and solve: $$ 0.025 \approx \sqrt{\frac{0.25(0.75)}{M}}, \qquad M \approx \frac{0.25(0.75)}{0.025^2} = \frac{0.1875}{0.000625} = 300. $$

About 300 independent draws are indicated under this planning approximation. The value $p=0.25$ is a planning assumption, not a guarantee of the eventual estimate or error. If the probability is unknown, a conservative planning calculation uses the maximum indicator variance, $p(1-p)\leq 0.25$, which occurs at $p=0.5$.

Choosing Computation and Checking Its Limits

When an integral has a convenient analytic solution, direct evaluation is typically preferable: it avoids Monte Carlo sampling error and may be faster. Numerical quadrature can also be effective for low-dimensional integrals. Monte Carlo averages are particularly valuable when the target is a complicated function of several posterior parameters or when direct integration is impractical. The earlier material on conjugate posterior summaries provides settings where direct calculations can be available; simulation is not automatically more accurate simply because it uses many draws.

The independence assumption matters. The formulas here apply to independent posterior draws. In later sampling procedures, successive draws may be dependent; treating such draws as independent and applying $s_g/\sqrt{M}$ can understate the actual simulation error. The number of stored draws alone does not establish independence. Methods for assessing and accounting for dependence require separate tools.

Useful checks include confirming that all draws come from the intended posterior, that the function $g$ is applied consistently, and that the estimated MCSE is small enough for the decision or report. For an estimated posterior probability, for example, inspect both the estimated probability and its MCSE. If the MCSE is large relative to the precision needed, generate more appropriate independent draws or use a suitable direct calculation. A small MCSE only says the numerical average is stable under repeated simulation from the specified posterior; it does not validate the model, remove uncertainty about the parameter, or correct a poor clinical definition of the threshold.

Confusing numerical error with posterior uncertainty. A posterior interval describes uncertainty about a parameter under the model and data; MCSE describes numerical variation in an estimate of a posterior summary.
Assuming the error falls as $1/M$. The variance falls as $1/M$, but the standard error falls as $1/\sqrt{M}$.
Using a dependent sample as if it were independent. The independent-draw formula does not account for correlation between draws.
Reporting only the number of draws. A large $M$ does not ensure adequate precision; report or assess the MCSE for the particular summary of interest.

Key takeaway: For independent posterior draws, average the sampled values of the function that answers the scientific question. The estimator’s standard error is estimated by the sample standard deviation of those values divided by $\sqrt{M}$; this quantifies Monte Carlo error, not posterior uncertainty.

Check Your Understanding

Answer each question, then compare with the solutions below.

For independent posterior draws $\theta^{(m)}$, what quantity estimates $E[g(\theta)\mid y]$?
If the sample standard deviation of sampled $g$-values is 0.12 and there are 144 independent draws, calculate the estimated MCSE.
For 25 independent draws, 18 satisfy a clinical threshold. Estimate the posterior probability and its MCSE.
If the number of independent draws is multiplied by nine, by what factor does the MCSE change, approximately?
Why is an MCSE not a credible interval for a treatment parameter?

Solutions. (1) Use $\widehat{\mu}_g=M^{-1}\sum_{m=1}^{M}g(\theta^{(m)})$. (2) $0.12/\sqrt{144}=0.12/12=0.0100$. (3) The estimate is $18/25=0.7200$; the MCSE is $\sqrt{0.72(0.28)/25}=\sqrt{0.008064}\approx0.08980$. (4) It is divided by $\sqrt{9}=3$. (5) MCSE describes the simulation variability of an estimated posterior summary across repeated draws; a credible interval describes posterior uncertainty about the parameter itself.