---
title: Empirical Partially Bayes Methods Review
url: https://www.emergentmind.com/topics/empirical-partially-bayes-methods
type: topic
---

# Empirical Partially Bayes Methods Review

Empirical partially Bayes methods are procedures in which a Bayesian hierarchy is used only in part: a prior, hyperprior, or nuisance-parameter distribution is estimated from data and then treated as fixed for subsequent posterior calculation, decision making, or testing. In the canonical form, one starts from a hierarchical model
\[
y_i\mid \theta_i\sim p(y_i\mid \theta_i),\qquad \theta_i\mid \eta\sim \pi(\theta_i\mid \eta),\qquad \eta\sim \rho(\eta),
\]
but replaces full posterior inference on \((\theta,\eta)\) with a plug-in step based on the marginal likelihood \(m(y\mid \eta)\). This places empirical Bayes, and the broader class of “partially Bayes” procedures, between classical frequentist methods and fully hierarchical Bayes: the final inference is posterior-shaped, but uncertainty in the learned prior component is not fully propagated [2605.21535]. Several recent works use “partially Bayes” to describe both classical type-II maximum-likelihood hyperparameter estimation and nuisance-only hierarchies in which parameters of interest are treated as fixed while nuisance distributions are learned across many related problems [1612.01403; 2512.08847].

## 1. Canonical formulation and the meaning of partial Bayes

In the standard empirical-Bayes construction, the joint density is
\[
p(y,\theta,\eta)=\rho(\eta)\prod_{i=1}^n\bigl[\pi(\theta_i\mid \eta)\,p(y_i\mid \theta_i)\bigr],
\]
and the fully Bayesian posterior is
\[
p(\theta,\eta\mid y)=\frac{p(y\mid \theta)\,\pi(\theta\mid \eta)\,\rho(\eta)}
{\iint p(y\mid \theta')\,\pi(\theta'\mid \eta')\,\rho(\eta')\,d\theta'\,d\eta'}.
\]
Empirical Bayes replaces integration over \(\eta\) by a type-II maximum-likelihood estimate,
\[
\hat\eta=\arg\max_\eta m(y\mid \eta),\qquad
m(y\mid \eta)=\prod_{i=1}^n\int p(y_i\mid \theta_i)\,\pi(\theta_i\mid \eta)\,d\theta_i,
\]
and then forms
\[
p_{\rm EB}(\theta_i\mid y_i)=p(\theta_i\mid y_i,\eta=\hat\eta).
\]
This plug-in posterior ignores posterior uncertainty in \(\eta\) and uses the data twice [2605.21535].

The same structural idea appears in broader partial-Bayes settings. One formulation assumes a partially specified prior family \(\pi(\theta)=g(\theta;\eta)\) with unknown hyperparameter \(\eta\), together with sampling model \(X\mid \theta\sim L(x\mid \theta)\), and seeks exact, frequentist-calibrated inference on \(\theta\) even though \(\eta\) is unspecified [1802.04050]. In another formulation, the prior itself is treated as an unknown object estimated from a population \(Z_1,\dots,Z_M\), after which Bayesian updating is carried out for a new individual. This is explicitly described as “empirical Bayes” or “partially Bayes,” and as a bridge between purely frequentist and purely Bayesian paradigms [1612.01403].

A related compound-decision perspective replaces the unknown prior by the empirical distribution \(G_n=n^{-1}\sum_{i=1}^n\delta_{\theta_i}\). Under exchangeability, the compound risk of a separable rule equals the Bayes risk under \(G_n\), so empirical partially Bayes methods can be viewed as attempts to estimate an oracle mixing distribution and then apply the usual Bayes rule [2404.03422].

## 2. Inferential targets, double use of data, and coherence

A central criticism of empirical Bayes is that it does not target the same inferential object as a fully hierarchical model. Dennis Lindley’s quip that “there is only one thing worse than a frequentist, and that is an empirical Bayesian” is presented not merely as caricature but as a technical objection: empirical Bayes uses the same data twice, conflates levels of a hierarchy, and produces posterior-shaped summaries whose uncertainty quantification differs from what a fully hierarchical model delivers [2605.21535].

This distinction has become more salient as empirical-Bayes ideas have been extended. Recent variants include empirical Bayes via probabilistic symmetries, empirical Bayes with implicit likelihoods through simulation-based inference, and empirical Bayes for combining experimental and observational data through calibration studies. The criticism advanced against these developments is that they often target inferential objects distinct from the posterior conditional on the realized data. A fully hierarchical alternative retains
\[
p(\theta,\eta\mid y)\propto \Bigl[\prod_{i=1}^n p(y_i\mid \theta_i)\,\pi(\theta_i\mid \eta)\Bigr]\rho(\eta),
\]
so credible intervals for \(\eta\) propagate directly into inference for \(\theta\), preserving coherence [2605.21535].

The computational argument that historically favored empirical Bayes has also weakened. Modern hierarchical Bayes routinely uses MCMC, including Gibbs sampling, Metropolis within Gibbs, and HMC; variational inference, including mean-field, structured VI, and black-box VI; neural amortization with normalizing flows and invertible networks; and sequential Monte Carlo or particle methods for streaming data. In this setting, the extra cost of sampling a hyperparameter or adding an extra variational dimension is described as modest relative to the benefit of avoiding data-twice use, obtaining correct coverage, and preserving coherence [2605.21535].

A closely related point appears in nuisance-only partial-Bayes testing. Partially Bayes \(p\)-values condition on ancillary nuisance-parameter statistics and pool information across units only through those ancillaries, thereby avoiding the posterior-predictive reuse of the same test statistic in both fitting and checking. This is explicitly contrasted with posterior predictive \(p\)-values, which double-use the data and are not uniform even in a fully Bayesian frame [2512.08847].

## 3. Tweedie identities, f-modeling, and coherent shrinkage

Much of empirical partially Bayes methodology is organized around identities that express posterior summaries in terms of marginal densities. In the normal-means model \(X\mid \theta\sim N(\theta,\sigma^2)\) with marginal density
\[
m(x)=\int \phi_\sigma(x-\theta)\,\pi(\theta)\,d\theta,
\]
Tweedie’s formula states
\[
\mathbb E[\theta\mid X=x]=x+\sigma^2\frac{d}{dx}\log m(x).
\]
Efron’s \(f\)-modeling estimates \(\hat m(x)\) by smoothing the observed \(\{x_i\}\) with splines, kernels, or histograms, differentiates \(\log \hat m(x)\), and plugs the estimated score into Tweedie’s identity. The important caveat is that a generic smoothed \(\hat m\) need not have the form \(\int \phi_\sigma(x-\theta)\,\hat\pi(\theta)\,d\theta\) for any probability density \(\hat\pi\); consequently, the resulting shrinkage estimate need not equal the posterior mean under any joint model. The horseshoe Tweedie formula is presented as a coherent alternative because it does arise from an explicit hierarchical prior [2605.21535].

For the hierarchical horseshoe prior,
\[
\theta_i\mid \lambda_i,\tau\sim N(0,\lambda_i^2\tau^2),\qquad
\lambda_i\sim C^+(0,1),\qquad
\tau\sim C^+(0,1),
\]
the shrinkage weight
\[
\kappa_i=\frac{1}{1+\lambda_i^2\tau^2/\sigma^2}
\]
yields the conditional posterior-mean identity
\[
\mathbb E[\theta_i\mid x_i,\tau]=\bigl(1-\mathbb E[\kappa_i\mid x_i,\tau]\bigr)x_i.
\]
The exact score function satisfies
\[
\frac{d}{dx}\log m_{\rm HS}(x)=\frac{\mathbb E[\theta\mid x]-x}{\sigma^2}.
\]
The half-Cauchy mixture produces aggressive shrinkage near zero and near-zero shrinkage for large \(|x|\), a combination emphasized as one reason to prefer coherent global-local priors over ad hoc score smoothing [2605.21535].

The same logic has been generalized to heteroscedastic normal means with unequal and unknown variances. If
\[
x_i\mid \theta_i,\sigma_i^2\sim N(\theta_i,\sigma_i^2),\qquad
s_i^2\mid \sigma_i^2\sim (\sigma_i^2/k)\chi_k^2,
\]
and \(f(x,s^2)\) denotes the joint marginal density under an unspecified prior \(g(\theta,\sigma^2)\), then under weighted loss \(L(\hat\theta,\theta)=(\hat\theta-\theta)^2/\sigma^2\) the Bayes rule is
\[
\hat\theta_B=\frac{E(\theta/\sigma^2\mid x,s^2)}{E(1/\sigma^2\mid x,s^2)},
\]
and can be written as
\[
\hat\theta_B
=x+\frac{k\,s^2\,f_x(x,s^2)}
{(k-2)\,f(x,s^2)-2\,s^2\,f_{s^2}(x,s^2)}.
\]
That framework also introduces a moment-generating-function representation of the full posterior, so that once \(f\) is known the posterior law is determined without explicitly specifying or estimating the prior \(g\). The method is proposed as a unified \(f\)-modeling framework for point estimation, uncertainty quantification, and hypothesis testing in heterogeneous data environments [2604.21865].

## 4. Major methodological variants

The contemporary literature contains several distinct empirical partially Bayes constructions, differing mainly in what is estimated and what inferential target is retained.

| Variant | Estimated object | Resulting target |
|---|---|---|
| Population empirical Bayes | Empirical population distribution or latent dataset prior | Population posterior or POP predictive density |
| Bayesian empirical Bayes via probabilistic symmetries | Directing measure \(G\) from an ergodic decomposition | Plug-in posterior conditional on \(\hat G\) |
| Implicit-likelihood empirical Bayes | Simulator-induced posterior approximation \(\hat p(\theta\mid x)\) | Amortized “EB posterior” or “population posterior” |
| Nuisance-only partially Bayes | Mixing law of nuisance parameters across units | Conditional tail-area \(p\)-values and moderated tests |

Population empirical Bayes (POP-EB) introduces a latent dataset \(D\) with prior \(F_{\rm emp}\), the empirical distribution of the observed sample. The POP predictive density is
\[
p_{\rm POP}(y^*\mid X)=\int p(y^*\mid D)\,p(D\mid X)\,dD,
\]
with bootstrap approximations yielding a MAP version based on the single best bootstrap and a full-Bayes approximation based on a weighted mixture over bootstraps. For complex models, POP-EB uses bumping variational inference (BUMP-VI), which interleaves bootstrap-specific ELBO gradients with predictive-score selection on the original data. Reported applications include linear regression on health data, a Gaussian mixture model of natural images, and latent Dirichlet allocation on scientific documents, where POP-EB improves held-out predictive scores relative to classical Bayesian inference [1411.0292].

A second family of methods estimates a directing measure directly. Under probabilistic symmetries, exchangeable arrays, spatial fields, and covariate-indexed hierarchies admit an ergodic decomposition over a directing measure \(G\), which is then estimated by maximizing
\[
\hat G=\arg\max_G\prod_i\int p(y_i\mid \theta_i)\,dG(\theta_i).
\]
Variational approximations with neural networks are used for \(G\). Closely related simulation-based methods replace explicit likelihood evaluation by a simulator \(p(x\mid \theta)\), then use neural posterior estimation, neural likelihood estimation, or neural ratio estimation to produce \(\hat p(\theta\mid x)\), sometimes treating it directly as an empirical-Bayes posterior even when no joint prior on \(\theta\) is written down [2605.21535].

In nuisance-only partially Bayes inference, the primary parameters remain fixed while only nuisance parameters are hierarchically modeled. One formulation observes summaries \((T_i,U_i)\) with \(T_i\) a test statistic and \(U_i\) ancillary for the primary parameter, posits
\[
(T_i,U_i)\mid \theta_i,\phi_i\sim p(t,u\mid \theta_i,\phi_i),\qquad
\phi_i\mid G\overset{\rm iid}{\sim}G,\qquad
G\sim \Pi,
\]
and defines the partially Bayes \(p\)-value by conditioning on all ancillaries but leaving out \(U_i\) when learning the nuisance law for unit \(i\). Under the fully Bayes frame this yields exact uniformity, and under large-\(K\) or large-\(n\) asymptotics it yields frequentist uniformity as well [2512.08847].

This nuisance-only perspective has become prominent in high-throughput biology. In the normal-means problem with unknown variances, conditional \(p\)-values based on a prior \(G\) over variances satisfy an Eddington/Tweedie-type formula, and nonparametric maximum-likelihood estimation of \(G\) yields plug-in \(p\)-values that can be combined with Benjamini–Hochberg for asymptotic FDR control [2303.02887]. Limma-trend is analyzed as computing approximate partially Bayes \(p\)-values that condition on residual sample variance and a unit-level summary, while nonparametric generalizations estimate the residual variance prior by NPMLE and retain asymptotic FDR control even when the trend is misspecified or inconsistently estimated [2605.21304]. Two-sample extensions with unequal variances use either a prior on the variance ratio \( \lambda_i=\sigma_{iA}^2/\sigma_{iB}^2 \) or a joint prior on \((\sigma_{iA}^2,\sigma_{iB}^2)\), both estimated by NPMLE, yielding asymptotically uniform \(p\)-values as the number of features grows with fixed replicate counts [2510.00432].

## 5. Frequentist calibration, asymptotics, and computation

A major line of work studies whether empirical Bayes and Bayes eventually agree. For a prior family \(\Pi(\cdot\mid \lambda)\) and plug-in posterior \(\Pi(\cdot\mid \hat\lambda_n,X_{1:n})\), Bayesian weak merging with every fixed-\(\lambda\) posterior is equivalent to weak consistency of the empirical-Bayes posterior under exchangeability and injectivity of \(\theta\mapsto P_\theta\). Under regularity conditions, marginal-likelihood empirical Bayes asymptotically selects the hyperparameter value whose prior density most favors the truth, and in regular parametric settings the plug-in posterior can merge strongly in total variation with the Bayesian posterior at the oracle hyperparameter. In nonparametric Dirichlet-process mixture settings, by contrast, only weak merging is generally possible because priors with different concentration parameters are mutually singular [1204.1470].

A different route to calibration is to avoid plug-in posteriors altogether. For partial-Bayes problems with prior family \(g(\theta;\eta)\) and unspecified \(\eta\), the inferential-model approach builds an association \(X=a(\theta,U)\), adds a prior mapping \(\theta=\theta(\eta,V)\), and then propagates a predictive random set through a reduced association for \(\eta\). The resulting plausibility interval
\[
C_\alpha(x)=\{\eta:\,Pl_{T(x)\mid h}(\eta)\ge \alpha\}
\]
satisfies
\[
P_{X,\eta,\eta^*}\{C_\alpha(X)\ni \eta\}\ge 1-\alpha
\]
under mild regularity. When the missing part of the prior is known, the construction recovers the Bayesian \(100(1-\alpha)\%\) credible interval; in large samples it converges to the optimal one [1802.04050].

For compound decisions and denoising, nonparametric maximum-likelihood estimation gives sharper oracle comparisons. In exchangeable models \(y_i\mid \theta_i\sim f(y\mid \theta_i)\), the Bayes rule under squared-error loss uses the posterior mean under the unknown mixing distribution \(G\). NPMLE estimates \(G\) by convex optimization without an explicit tuning parameter, and the resulting empirical-Bayes posterior mean enjoys oracle inequalities: for sub-Gaussian \(G\), regret relative to the oracle Bayes rule is
\[
O\bigl(n^{-1}(\log n)^5\bigr),
\]
while heavier-tailed classes achieve the stated nonparametric rates up to logarithmic factors [2404.03422]. In multivariate heteroscedastic Gaussian mixtures,
\[
X_i\mid \theta_i\sim N_d(\theta_i,\Sigma_i),
\]
finite-dimensional approximations to the NPMLE still yield average Hellinger bounds for the fitted marginals and regret bounds showing that the empirical-Bayes posterior means nearly match the oracle posterior means [2109.03466].

Empirical partially Bayes computation has also developed into a distinct technical area. One strategy runs a single MCMC chain under a reference prior \(\nu_{h_1}\), then reweights to estimate the entire curves
\[
m_y(h)=\int p_\theta(y)\,\nu_h(d\theta),\qquad
I_g(h)=E[g(\theta)\mid Y=y;h]
\]
for all \(h\in\mathcal H\). Under geometric ergodicity and empirical-process conditions, this yields uniform strong consistency, functional central limit theorems, argmax consistency for empirical-Bayes hyperparameter estimates, and simultaneous confidence bands [1807.02191]. A related multi-chain method chooses skeleton priors, estimates normalizing-constant ratios by reverse logistic regression, and then uses importance sampling with control variates to estimate Bayes factors and posterior expectations across a continuum of hyperparameters [1202.5160].

## 6. Domains of application and contemporary practice

In systems medicine, empirical Bayes is used to construct informative priors from many patients before carrying out patient-specific inference. One comparison considered four priors: a non-informative prior, a nonparametric maximum-likelihood prior, a maximum penalized likelihood prior, and a doubly-smoothed MLE prior. In a harmonic-oscillator example with \(M=300\), the reported posterior success-rate error standard deviations were \(0.1634\) for the uniform prior, \(0.1215\) for NPMLE, \(0.1097\) for DS-MLE, and \(0.1003\) for MPLE. In a \(d=115\) ODE model of the human menstrual cycle, estimated priors produced much sharper marginal posteriors than the non-informative baseline, while DS-MLE and MPLE avoided the atomic collapse of NPMLE [1612.01403].

In high-dimensional prediction, empirical Bayes is used to learn from many variables and from external “co-data.” In Gaussian ridge models, hyperparameters can be estimated by marginal likelihood or moment matching; in group-regularized logistic regression, hybrid empirical Bayes–full Bayes models retain EB-estimated group multipliers while placing a hyperprior on a global shrinkage parameter. In spike-and-slab variable selection, co-data enter through
\[
\logit \pi_j=C_j^\top \alpha,
\]
with \(\alpha\) estimated by Gibbs-EB or moment-based EB. These constructions are described as especially useful when the prior has multiple parameters that encode a priori information on variables [1709.04192].

In matrix completion, empirical Bayes has been adapted to binary observation models. For a partially observed binary matrix \(Y\), the method places a row-wise Gaussian prior \(m_i\sim N_q(\mu,\Sigma)\), estimates \((\mu,\Sigma)\) by MCEM from the marginal likelihood, and samples from the posterior using Albert–Chib augmentation and Gibbs updates. In the continuous Gaussian limit, the construction connects to the Efron–Morris singular-value shrinker; in the 1-bit setting the Gibbs sampler performs adaptive singular-value shrinkage through \(\Sigma\). Reported real-data results show that the heterogeneous-column model EB2 attained classification accuracy \(73.05\%\), cross-entropy \(0.5326\), and \(ECE\approx 0.0198\) on Jester, and accuracy \(67.8\%\), cross-entropy \(0.5983\), and \(ECE\approx 0.036\) on MovieLens 100K [2605.09509].

Large-scale inference in genomics and related fields has become a principal domain for empirical partially Bayes methods. Limma-trend can be interpreted as a parametric partially Bayes procedure that places a trend-dependent prior on residual variances and then computes moderated-\(t\) \(p\)-values conditional on the residual sample variance and a unit-level summary. Nonparametric generalizations estimate either the residual multiplier distribution or the full conditional variance distribution and obtain asymptotic FDR control under dense signals [2605.21304]. In unequal-variance two-sample testing, VREPB and DVEPB estimate respectively a prior on \(\sigma_{iA}^2/\sigma_{iB}^2\) or a joint prior on \((\sigma_{iA}^2,\sigma_{iB}^2)\), again by NPMLE, to produce partially Bayes \(p\)-values with asymptotic type-I control [2510.00432].

Pharmacovigilance provides another compound-decision setting. In spontaneous reporting systems, counts satisfy
\[
N_{ij}\mid \lambda_{ij}\sim \mathrm{Poisson}(E_{ij}\lambda_{ij}),
\]
with signal strengths \(\lambda_{ij}\) drawn from an unknown prior \(g(\lambda)\). Nonparametric empirical Bayes methods estimate \(g\) by discrete-support convex optimization, spline-based exponential-family deconvolution, \(K\)-gamma mixtures, or sparse overfitted general-gamma mixtures. Posterior summaries include the empirical Bayes geometric mean
\[
EBGM_{ij}=\exp\{E[\log \lambda_{ij}\mid N_{ij}]\}
\]
and the lower-bound \(EB05\), the \(5\)th posterior percentile. The package pvEBayes implements these methods and graphical summaries; in one statin example, approximately \(104\) AE-drug signals with \(EB05>1.001\) were reported among \(6\) statins [2512.01057].

Across these applications, a consistent contemporary recommendation is to redeploy modern variational, neural, and simulation-based machinery in service of fully hierarchical Bayes rather than to bypass the hierarchy. In that view, one should write down a hyperprior for \(\eta\) or for the directing measure \(G\), use scalable inference to approximate the full posterior, replace incoherent marginal-score smoothing by g-modeling or continuous global-local priors, and prefer online VB or sequential Monte Carlo for streaming data. The broad implication is not that empirical partially Bayes methods are obsolete, but that they are best understood as approximations with distinct inferential targets, useful theoretical guarantees in some regimes, and increasingly explicit alternatives when coherence is the primary objective [2605.21535].

Source: https://www.emergentmind.com/topics/empirical-partially-bayes-methods