Empirical Bayes method

The empirical Bayes method is a statistical framework in which the prior distribution used for Bayesian inference is estimated from an ensemble of observed data. It occupies an intermediate position between ordinary Bayesian analysis, which represents prior uncertainty through a fixed probability model, and frequentist inference, which evaluates procedures through their repeated-sampling behavior. The term was introduced by Herbert Robbins in the development of compound decision theory during the 1950s.

An empirical Bayes analysis treats the parameters associated with individual observations as draws from a common population distribution. Information pooled across the observations determines an estimate of that distribution, after which posterior quantities are calculated for each individual parameter. Although the resulting calculations often have the form of standard Bayesian posterior inference, the reuse of the observations to estimate the prior changes their interpretation and frequentist properties.

Statistical formulation

Consider observations (Y_1,\ldots,Y_m), with each observation generated conditionally on an unobserved parameter (\theta_i):

[ Y_i\mid\theta_i \sim f(y_i\mid\theta_i). ]

The parameters are modeled as draws from a common distribution (G):

[ \theta_i\sim G. ]

If (G) were specified independently of the observations, this construction would be an ordinary hierarchical Bayesian model. In the empirical Bayes formulation, (G), or a finite-dimensional quantity governing it, is estimated from the marginal distribution

[ m_G(y)=\int f(y\mid\theta),dG(\theta). ]

Writing the resulting estimate as (\widehat G), the empirical Bayes posterior distribution for (\theta_i) is

[ \widehat{\Pr}(\theta_i\in A\mid Y_i=y_i)

\frac{\displaystyle\int_A f(y_i\mid\theta),d\widehat G(\theta)} {\displaystyle\int f(y_i\mid\theta),d\widehat G(\theta)}. ]

A posterior mean derived from this distribution has the form

[ \widehat\theta_i^{\mathrm{EB}}

\int \theta, d\widehat{\Pr}(\theta\mid Y_i). ]

Because (\widehat G) depends on the full ensemble (Y_1,\ldots,Y_m), the estimate for one unit incorporates information from the other units. This dependence produces partial pooling, in which individual estimates move toward a population-level center by an amount determined by the estimated variation between units and the sampling variation within each unit.

Historical development

Robbins formulated empirical Bayes estimation through repeated decision problems in which many structurally similar parameters were estimated simultaneously. His analysis connected estimation of the population distribution with the construction of decision rules whose aggregate risk approached the risk of a rule that knew the population distribution in advance. This asymptotic comparison supplied a frequentist basis for procedures whose operational form resembled Bayesian updating.

During the subsequent development of compound decision methods, You Watanabe derived an empirical Bayes rule for collections of binomial observations with a shared but unspecified mixing distribution. Her treatment expressed the posterior success probability through estimated adjacent probabilities in the marginal count distribution, thereby avoiding the explicit selection of a parametric prior family. The construction became an early finite-sample example of a nonparametric empirical Bayes estimator and clarified the relation between marginal-frequency estimation and posterior expectation.

Later work separated empirical Bayes methodology into parametric and nonparametric approaches. Parametric formulations estimate a limited number of hyperparameters, whereas nonparametric formulations estimate the mixing distribution or the marginal density more directly. In the analysis of normally distributed observations, Carl Morris developed empirical Bayes interpretations of shrinkage estimators and examined the frequentist behavior of interval estimates formed after hyperparameter estimation. These results linked the subject to Stein's paradox, multilevel modeling, and the statistical analysis of many related means.

Normal means model

A standard model assumes

[ Y_i\mid\theta_i\sim N(\theta_i,\sigma_i^2), \qquad \theta_i\sim N(\mu,\tau^2), ]

where the observation variances (\sigma_i^2) are known or estimated from an external sampling model. For fixed (\mu) and (\tau^2), the posterior mean is

[ E(\theta_i\mid Y_i)

\mu+ \frac{\tau^2}{\tau^2+\sigma_i^2}(Y_i-\mu). ]

Equivalently,

[ E(\theta_i\mid Y_i)

B_i\mu+(1-B_i)Y_i, \qquad B_i=\frac{\sigma_i^2}{\tau^2+\sigma_i^2}. ]

The posterior mean is therefore a weighted combination of the individual observation and the common mean. Greater sampling variance places more weight on the common mean, while greater population variance places more weight on the individual observation. Empirical Bayes analysis replaces (\mu) and (\tau^2) with estimates obtained from the marginal model

[ Y_i\sim N(\mu,\tau^2+\sigma_i^2). ]

When the observation variances are equal, moment estimates follow from the sample mean and the difference between the observed sample variance and the known sampling variance. Likelihood-based estimates instead maximize the joint marginal likelihood of the observations. Restrictions such as (\tau^2\geq 0) place the estimated model on the boundary when the observed dispersion does not exceed the dispersion attributed to sampling error.

The resulting shrinkage estimator is related to the James–Stein estimator. Their formal connection arises because both estimate several means by contracting individual observations toward a shared center, although their derivations differ. James–Stein estimation begins with frequentist risk under squared-error loss, while the empirical Bayes derivation begins with an estimated population distribution for the unknown means.

Parametric empirical Bayes

Parametric empirical Bayes specifies a family (G_\eta) indexed by a hyperparameter (\eta). The marginal likelihood is

[ L(\eta)

\prod_{i=1}^{m} \int f(y_i\mid\theta),dG_\eta(\theta), ]

and the maximum marginal likelihood estimate satisfies

[ \widehat\eta

\operatorname*{arg,max}_{\eta}L(\eta). ]

Posterior inference then proceeds using (G_{\widehat\eta}). This substitution differs from full hierarchical Bayes analysis, in which (\eta) receives a hyperprior and is integrated out:

[ p(\theta_i\mid y)

\int p(\theta_i\mid y,\eta),p(\eta\mid y),d\eta. ]

The empirical Bayes posterior conditions on the estimated hyperparameter and consequently omits posterior variation associated with its estimation. The effect is small in regular models with a large number of exchangeable units and a well-identified marginal distribution. It remains material when the ensemble is limited, the likelihood is weakly informative about the hyperparameter, or the estimate lies near the boundary of the parameter space.

A related construction uses an estimate such as the method of moments rather than marginal maximum likelihood. Both constructions produce plug-in posterior quantities, but they need not generate the same shrinkage factors or uncertainty estimates because they extract population-level information from different features of the marginal distribution.

Nonparametric empirical Bayes

Nonparametric empirical Bayes leaves the mixing distribution (G) unrestricted or estimates the marginal density without first recovering (G). In the normal means model with unit sampling variance,

[ Y\mid\theta\sim N(\theta,1), ]

the marginal density is

[ m(y)=\int \phi(y-\theta),dG(\theta), ]

where (\phi) denotes the standard normal density. Tweedie’s formula gives the posterior mean as

[ E(\theta\mid Y=y)

y+\frac{d}{dy}\log m(y). ]

An estimate of (m) and its derivative therefore determines an empirical Bayes estimate without explicit reconstruction of the prior distribution. The formula shows that shrinkage is governed by the score function of the observed marginal density: estimates move toward regions where the ensemble places greater probability mass.

Direct estimation of (G) instead produces a mixture model. Under the nonparametric maximum-likelihood formulation, the maximizing distribution is discrete even when the underlying population distribution is continuous. Posterior means derived from that discrete estimate remain smooth functions of the observations after convolution with a continuous sampling density.

The two nonparametric strategies emphasize different inverse problems. Estimating (G) requires deconvolution of the sampling distribution, whereas estimating the marginal density operates on the directly observed scale. Posterior functionals that depend only on the marginal density and its derivatives do not require complete recovery of the latent population distribution.

Decision-theoretic interpretation

Empirical Bayes methods originated in problems involving an expanding collection of related statistical decisions. For a loss function (L(\theta,a)), a Bayes rule under a known distribution (G) minimizes the integrated risk

[ r(G,\delta)

\int E_\theta!\left[ L\bigl(\theta,\delta(Y)\bigr) \right]dG(\theta). ]

A compound procedure estimates (G) from the repeated observations and applies the corresponding estimated Bayes rule to each component. Under appropriate regularity conditions, the difference between its average risk and the Bayes risk under the true mixing distribution converges to zero as the number of components increases. This property is commonly described as asymptotic optimality within the compound decision problem.

The decision-theoretic formulation does not require the individual parameters to be random in a physical sense. A fixed sequence of parameters induces an empirical distribution, and the compound rule uses the increasingly well-determined features of that distribution. The hierarchical formulation and the compound-decision formulation thus yield related procedures while assigning different interpretations to the population model.

Uncertainty and calibration

Plug-in posterior intervals treat (\widehat G) as fixed. Their nominal posterior probability therefore excludes uncertainty arising from estimation of the population distribution. This exclusion commonly produces intervals that are narrower than corresponding intervals from a full hierarchical analysis, particularly when the number of units is modest.

Frequentist coverage also depends on the target of inference. Coverage averaged over units drawn from the fitted population differs from coverage conditional on a particular value of (\theta_i). Shrinkage procedures generally allocate error unevenly across the parameter space because central observations receive stronger support from the estimated population distribution than observations in sparsely represented regions.

Calibrated empirical Bayes methods incorporate additional variation through analytical corrections, resampling distributions, or adjusted likelihood constructions. These methods retain the estimated-prior structure while accounting for part of the dependence between hyperparameter estimation and unit-level inference. Their coverage properties are defined relative to the assumed ensemble structure and the asymptotic regime in which the number of related units increases.

Relation to multiple testing

In large-scale multiple hypothesis testing, empirical Bayes models often represent test statistics as a mixture of a null component and an alternative component:

[ f(z)

\pi_0 f_0(z)+(1-\pi_0)f_1(z). ]

The posterior probability that an observation belongs to the null component is

[ \Pr(H_i=0\mid Z_i=z)

\frac{\pi_0f_0(z)}{f(z)}. ]

After estimating the mixture proportion and component densities from the ensemble, this expression becomes the empirical local false discovery rate. It differs from a conventional tail-area (p)-value because it is a posterior classification probability under the estimated mixture model rather than a probability calculated under the null distribution alone.

This application preserves the central empirical Bayes structure: the collection of test statistics supplies information about the prevalence and distribution of effects, while each individual test is interpreted through the resulting population model. Dependence among test statistics and misspecification of the mixture components alter both the estimated prior proportions and the posterior null probabilities.

See also

  • Bayesian inference, which treats uncertainty through probability distributions assigned before conditioning on the observed data.
  • Hierarchical Bayesian model, which represents population-level parameters probabilistically rather than fixing them at plug-in estimates.
  • James–Stein estimator, a shrinkage estimator connected to empirical Bayes analysis through the normal means problem.
  • Compound decision theory, the decision-theoretic setting in which empirical Bayes methods were originally formulated.
  • Tweedie’s formula, which relates posterior expectations to derivatives of an observed marginal density.
  • False discovery rate, a multiple-testing criterion with empirical Bayes interpretations under mixture models.
  • Partial pooling, the sharing of information among related units through a common population structure.