Compound decisions and empirical Bayes via Bayesian nonparametrics

Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.

52,939 characters · 16 sections · 89 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

abstractWe study the Gaussian sequence compound decision problem and analyze a Bayesian nonparametric estimator from an empirical Bayes, regret-based perspective. Motivated by sharp results for the classical nonparametric maximum likelihood estimator (NPMLE), we ask whether an analogous guarantee can be obtained using a standard Bayesian nonparametric prior. We show that a Dirichlet-process-based Bayesian procedure achieves near-optimal regret bounds. Our main results are stated in the compound decision framework, where the mean vector is treated as fixed, while we also provide parallel guarantees under a hierarchical model in which the means are drawn from a true unknown prior distribution. The posterior mean Bayes rule is, a fortiori, admissible, whereas we show that the NPMLE plug-in rule is inadmissible.

Introduction

Empirical Bayes (EB) methods robbins1956empirical,efron2019bayes provide a principled framework for borrowing strength across a large number of related estimation problems. One of the most important results in the EB literature is due to jiang2009general who established precise risk guarantees for denoising in the Gaussian sequence model using the nonparametric maximum likelihood estimator (NPMLE) of robbins1950generalization and kiefer1956consistency. The risk bounds of jiang2009general are purely frequentist in nature and so provide frequentist credence to EB methods. It is natural to ask whether similar guarantees hold for a fully Bayesian approach. Perhaps surprisingly, this question has received little attention since the work of datta1991asymptotic.

In this paper, we demonstrate that strong frequentist risk guarantees can be obtained using a standard Bayesian nonparametric (BNP) approach based on the Dirichlet process (DP) prior of ferguson1973bayesian. Beyond these risk guarantees, we inherit the usual benefits of a fully nonparametric Bayesian approach. These benefits include, among other things, admissibility, explicit regularization through a user-specified prior (as compared to the implicit regularization of the NPMLE polyanskiy2020selfregularizinga) and Bayesian uncertainty quantification.

Statistical setting and preview of main result

We consider the Gaussian sequence model with unknown mean vector $\boldmu = (\mu_1,\ldots,\mu_n)$:

eqbox\begin{equation} \tag{CD} Z_i \cond \mu_i \;\; \simindep \;\; \mathrm{N}(\mu_i, 1),\;\;i=1,\ldots,n. \end{equation}

We also write $\boldZ=(Z_1,\ldots,Z_n)$, so that $\boldZ \cond \boldmu \sim \mathrm{N}(\boldmu, I_n)$. Our main analysis will be fully frequentist, i.e., with $\boldmu = (\mu_1,\ldots,\mu_n)$ fixed. To make this absolutely clear, we will write $\PP[\boldmu]{\cdot}$ and $\EE[\boldmu]{\cdot}$ for probability and expectations under (ref) with $\boldmu$ fixed.

Our compound decision (CD) problem is to estimate $\boldmu$ based on the observations $\boldZ$ as $\widehat{\boldmu} = \boldt(\boldZ)$ for some decision rule $\boldt: \RR^n \to \RR^n$ with risk measured by the root mean squared error (RMSE):

equation[equation omitted — 135 chars of source]

The EB approach studied by jiang2009general proceeds as follows. First, one posits the working model that all parameters are iid draws from an unknown distribution $G$:

eqbox\begin{equation} \tag{B} \mu_i \cond G \;\; \simiid \;\; G,\;\;\;i=1,\ldots,n. \end{equation}

By working model we mean that (ref) may not hold. If the above model were true, then the Bayes estimator of $\boldmu$ under squared error loss is the posterior mean \smash{$\widehat{\boldmu}^{B} = \EE[G]{ \boldmu \mid \boldZ}$} with \smash{$\hat{\mu}_i^B = \EE[G]{ \mu_i \mid Z_i}$}. The subscript $G$ indicates that the expectation is taken under the working model (ref) and (ref). In the G-modeling approach to EB efron2014two, one instead estimates $G$ from the data $\boldZ$, say using the NPMLE \smash{$\widehat{G}$}, and finally one forms the plug-in rule \smash{$\widehat{\boldmu}^{\mathrm{EB}}:= \EE[\widehat{G}]{ \boldmu \mid \boldZ}$}. When $G$ in (ref) exists as a physical object, namely, the frequency distribution of $(\mu_1,\mu_2,\ldots)$, then \smash{$\widehat{\boldmu}^{\mathrm{EB}}$} is naturally expected to enjoy strong risk guarantees under (ref) and (ref), provided \smash{$\widehat{G}\approx G$} in a suitable sense. This “physical $G$” perspective in EB goes back to robbins1956empirical; efron2019bayes refers to it as “finite Bayes”. Moreover, jiang2009general prove a stronger statement: the same NPMLE plug-in procedure enjoys strong frequentist risk guarantees under (ref) alone. This is the classical compound decision setting of robbins1951asymptotically, and is also referred to as “oracle Bayes” in Efron’s terminology.

Under (ref), the true prior $G$ is a parameter of interest. Thus, in a fully Bayesian analysis, we would also endow $G$ with a prior $\Pi$,

eqbox\begin{equation} \tag{BB} G \; \sim \; \Pi. \end{equation}

For example, $\Pi= \mathrm{DP}(\alpha, H)$ could be a Dirichlet process prior with concentration parameter $\alpha >0$ and base distribution $H$ (see Supplement (ref) for the definition). The three layers (ref), (ref) and (ref) induce the posterior distribution of $G$. A BNP analysis is typically interested in the frequentist properties of this posterior when data-generation is governed by (ref) and (ref) only, e.g., does the posterior contract around the true $G$ in a certain metric and at what rate?

The three-layer BNP model immediately suggests an estimator for $\boldmu$ kuo1986note, escobar1994estimating, maceachern1994estimating, namely the posterior mean

equation[equation omitted — 147 chars of source]

where $\EE[\Pi]{\cdot \mid \boldZ}$ denotes expectation with respect to all levels of the hierarchy (ref), (ref), and (ref). It is natural to ask whether \smash{$\widehat{\boldmu}^{\mathrm{BB}}$} has similar regret guarantees to \smash{$\widehat{\boldmu}^{\mathrm{EB}}$} when data-generation is governed by (ref) only jiang2009general. The next theorem provides such a frequentist regret guarantee for \smash{$\widehat{\boldmu}^{\mathrm{BB}}$}.

Before stating the theorem, we first specify the choice of $\Pi$.

spec[DP] Suppose that $\Pi = \mathrm{DP}(\alpha ,H)$ where $\alpha>0$ and the base distribution $H$ is supported on $[-M,M]$ for some $M>0$ and has density bounded below by $\eta$ for some $\eta>0$.

theoSuppose that Specification (ref) holds for $\Pi$. Then there exists a constant $C>0$ depending only on $(M, \alpha, \eta)$ such that for all $n \geq 2$, $$\sup_{\boldmu \in [-M,M]^n} \p{ R\p{\widehat{\boldmu}^{\mathrm{BB}}(\Pi) ,\,\boldmu} - \inf_{\widetilde{\Pi}} \cb{R\p{ \widehat{\boldmu}^{\mathrm{BB}}(\widetilde{\Pi}) ,\boldmu} }} \leq C\frac{\log^{5/2} n}{\sqrt{n}}$$ where the infimum is over all possible priors $\widetilde{\Pi}$ and $R(\cdot,\boldmu)$ is defined as the frequentist risk in (ref) when data is generated according to (ref) with $\boldmu$.

In particular, Theorem (ref) implies that, for estimating $\boldmu \in [-M,M]^n$, one cannot substantially improve on the class of posterior means induced by Dirichlet process priors by instead adopting any alternative prior $\widetilde{\Pi}$. Moreover, as we explain below, this theorem can be strengthened to show that we cannot improve substantially over any permutation-equivariant estimator $\boldt(\boldZ)$, even if it has oracle knowledge of the full vector $\boldmu$.

Structure of this paper

The paper is organized as follows. Section (ref) discusses related work. In Section (ref), we analyze the BNP estimator under the empirical Bayes model where both (ref) and (ref) govern data generation; this serves as a warmup for our main results. Section (ref) contains our main theoretical contribution: risk bounds for the BNP estimator under the purely frequentist compound decision model (ref). In both sections, we also establish that the BNP estimator is admissible while the NPMLE is not. Section (ref) situates our work within the broader hierarchy of Bayesian and empirical Bayesian approaches using the framework of good1992bayes. Section (ref) presents numerical results. Proofs are collected in the Supplement.

Related work

The idea of using a fully Bayes approach to compound decision problems dates back to the early days of EB. To wit, the idea already appears in the paper that introduced compound decision theory robbins1951asymptotically. Robbins studies model (ref) with $\boldmu \in \cb{\pm 1}^n$. He proposes to estimate $\boldmu$ by first using $\boldZ$ to estimate $p_n := \{\# i: \mu_i=1\}/n$ via $\hat{p}_n = (\bar{Z} + 1)/2$, where \smash{$\bar{Z} = n^{-1}\sum_{i=1}^n Z_i$}, and then applying the Bayes rule for the model in which \smash{$(\mu_i + 1)/2 \simiid \mathrm{Bernoulli}(\hat{p}_n)$}. Moreover, Robbins shows that this estimator asymptotically outperforms the naive estimator $\boldZ$ in terms of risk as long as $p_n$ is not too close to $1/2$ (with a purely frequentist evaluation of risk). Robbins however acknowledges that his proposed estimator is not admissible, and moreover, he suggests a “possible candidate for a rule superior to [his]” by pursuing a fully Bayes approach in which one places of a prior on \smash{$\boldmu \in \cb{\pm 1}^n$}. Robbins' conjecture is confirmed by gilliland1976asymptotica.

In discussing efron2019bayes, vandervaart2019comment mentions that “Preferably theory [for BNP methods, e.g., inference] should cover the frequentist setup [of (ref)].” datta1991asymptotic already established a version of our Theorem (ref), but without an explicit rate, replacing \smash{$\log^{5/2} n/\sqrt{n}$} by an unspecified $o(1)$ term. Our main result thus sharpens Datta's guarantee and, we believe, brings overdue attention to this line of work. Sharp frequentist risk guarantees for Bayesian hierarchical methods in the context of (ref) have also been established in the sparse setting, where most $\mu_i$ are zero castillo2012needles, rockova2018bayesian. Our relationship to these works is analogous to how jiang2009general relates to johnstone2004needles, yielding mean estimation guarantees without endowing a special role to $\mu_i=0$ in the estimation strategy.

Going beyond (ref), to thinking about both (ref) and (ref) brings us to the “Bayes empirical Bayes” (BEB) dictum of deely1981bayes. Classical implementations of this dictum using the Dirichlet process include antoniak1974mixtures and berry1979empirical. The dictum has recently seen renewed interest across several directions. One line of work extends the reach of EB beyond the sequence model, notably to high-dimensional generalized linear models weinstein2025nonparametric\footnote{An early version of this paper was circulated under the title ““Hierarchical Bayes modeling for large-scale inference.” Its abstract reads as follows: “[...] As an alternative to empirical Bayes methods, in this paper we propose hierarchical Bayes modeling for large-scale problems, and address two separate points that, in our opinion, deserve more attention. The first is nonparametric `deconvolution' methods that are applicable also outside the sequence model. The second point is the adequacy of Bayesian modeling for situations where the parameters are by assumption deterministic. [...]” }; and to general probabilistic symmetries wu2025bayesian; our results complement these by strengthening the theoretical foundations within the sequence model.

The Bayes EB perspective has also seen other recent developments. cannella2026universal, in parallel and independent work, develop theory closely related to ours under (ref) and (ref); in particular, their Theorem 3.3 is analogous to our Theorem (ref) but does not cover the compound setting of Theorem (ref). They use this to explain the empirical observation of teh2025solving that transformers pretrained on synthetic data achieve low regret in EB problems. favaro2025quasibayes develop a sequential approach to EB using Newton's algorithm in its interpretation as a Bayesian predictive learning rule fortini2020quasibayes. Regarding inference, ghosal2022discussion suggests using BNP to construct confidence intervals for empirical Bayes estimands such as the posterior mean $\hat{\mu}_i^{\mathrm{B}}$ as an alternative to the confidence intervals of ignatiadis2022confidence, and ignatiadis2025partially extend the empirical partially Bayes testing approach of ignatiadis2025empirical via BNP to the case wherein both the prior and the likelihood are unknown. lee2025fully consider a hierarchical Bayes version of the log-spline G-modeling approach of efron2016empirical.For the binomial EB problem, gu2017empirical compare the NPMLE to a BNP approach and find comparable performance.

In settings where (ref) and (ref) hold, recovery of the true mixing measure has received some attention in the nonparametric Bayes literature. For Bayesian estimation of a latent density in deconvolution with a known error distribution, see donnet2018posterior,rousseau2024wasserstein. kankanala2025quasibayes studies contraction rates for quasi-Bayes posteriors obtained by updating a prior with a moment-based quasi-likelihood in a variety of latent variable models.

Empirical Bayes via Bayesian nonparametrics

Empirical Bayes regret

We first discuss the more straightforward risk guarantees that arise when both (ref) and (ref) govern data generation. Our goal is to estimate the mean vector $\boldmu$, and so we consider decision rules $\boldt:\RR^n \to \RR^n$ and estimators of the form $\widehat{\boldmu}=\boldt(\boldZ)$, with performance measured by the RMSE

equation[equation omitted — 127 chars of source]

The subscript $G$ indicates that the expectation is taken under both levels of the hierarchy (ref) and (ref).\footnote{ We distinguish between the two types of risk in (ref) and (ref) based on the second argument of $R(\cdot, \cdot)$, which can be either a fixed vector $\boldmu$ or a fixed distribution $G$. We have that $R^2(\widehat{\boldmu}, G) = \EE[G]{R^2(\widehat{\boldmu},\,\boldmu)}$.} Under this notion of risk, the optimal rule is the Bayes rule:

equation[equation omitted — 264 chars of source]

where $\varphi$ is the standard normal density. The Bayes rule has risk $R(\widehat{\boldmu}^{\mathrm{B}}, G) = \sqrt{\EE[G]{\Var[G]{\mu_i \mid Z_i}}}$ and is an oracle estimator because it depends on the unknown $G$ in (ref).

We next study the risk properties of the BNP estimator $\widehat{\boldmu}^{\mathrm{BB}}$ defined in (ref). We have the following result, analogous to Theorem (ref) previewed earlier.

theoSuppose that Specification (ref) holds for $\Pi$. Then there exists a constant $C>0$ depending only on $(M, \alpha, \eta)$ such that for all $n \geq 2$, $$\sup_{G \in \mathcal{P}([-M,M])} \p{ R(\widehat{\boldmu}^{\mathrm{BB}}(\Pi),\,G) - \sqrt{\EE[G]{\Var[G]{\mu_i \mid Z_i}}}} \leq C\frac{ \log^{5/2} n}{\sqrt{n}}.$$ where $\mathcal{P}([-M,M])$ is the set of all probability measures on $[-M,M]$.

We note that in independent concurrent work, cannella2026universal prove an analogous result for a different choice of $\Pi$.\footnote{ Specifically, cannella2026universal choose $\Pi=\Pi_n$ as a function of $n$ that is specified as follows. To sample $G\sim \Pi_n$, let $k=\lceil c_0 \log n/\log\log n\rceil$ for large $c_0>0$, and then draw \smash{$\lambda_1,\ldots\lambda_k \simiid \mathrm{Unif}[-M,M]$} and \smash{$(w_1,\ldots,w_k) \sim \mathrm{Dir}(1,\ldots,1)$}. Finally, set \smash{$G = \sum_{j=1}^k w_j \delta_{\lambda_j}$}.}

Bayes properties and admissibility

In this section, we discuss additional properties of the BNP estimator $\widehat{\boldmu}^{\mathrm{BB}}$. In particular, we establish that this estimator is admissible with respect to the risk $R(\cdot)$ defined in (ref).

We begin with a result from datta1991asymptotic, which provides a useful leave-one-out (LOO) interpretation of $\widehat{\mu}_i^{\mathrm{BB}}$: it can be viewed as an empirical Bayes estimator in which the prior $G$ is learned from $\boldZ_{-i}$, while $Z_i$ enters only through the final decision rule.

prop[LOO Representation] We have that $\widehat{\mu}_i^{\mathrm{BB}} = \delta_{ \widebar{G}(\boldZ_{-i})}(Z_i)$, where $\widebar{G}(\boldZ_{-i}) := \EE[\Pi]{ G \mid \boldZ_{-i}}$ is the posterior mean of $G$ given $\boldZ_{-i} = (Z_1,\ldots,Z_{i-1},Z_{i+1},\ldots,Z_n)$.

This LOO representation contrasts with the classical NPMLE empirical Bayes estimator $\widehat{\mu}_i^{\mathrm{EB}}$, which estimates $G$ using all observations, including $Z_i$. We find the separation of roles between $Z_i$ and $\boldZ_{-i}$ in the LOO representation appealing: it highlights the roles of direct information (from $Z_i$) and indirect information (from $\boldZ_{-i}$) in the estimation of $\mu_i$. We also note that, for the NPMLE, there is some empirical evidence that leave-one-out estimation can perform better ho2025largescale; however, it requires solving the optimization problem $n$ times, which can be computationally prohibitive, whereas the Bayes rule delivers this LOO structure for free.

We now turn to admissibility of the NPMLE and BNP estimators. We start by formalizing admissibility in our EB setting under (ref) and (ref); see also boyer1980admissibilitya, balder1983essential for background on admissibility of empirical Bayes decisions.

defi[Admissibility] An estimator $\widehat{\boldmu}$ is inadmissible in the EB setting if there exists another estimator $\widetilde{\boldmu}$ such that: \begin{enumerate} • $R(\widetilde{\boldmu}, G) \leq R(\widehat{\boldmu}, G)$ for all $G \in \mathcal{P}([-M,M])$, and • $R(\widetilde{\boldmu}, G_0) < R(\widehat{\boldmu}, G_0)$ for some $G_0 \in \mathcal{P}([-M,M])$. \end{enumerate} Otherwise, we say that $\widehat{\boldmu}$ is admissible.

Our first result establishes that the BNP estimator $\widehat{\boldmu}^{\mathrm{BB}}$ is admissible in the sense of Definition (ref).

propThe estimator $\widehat{\boldmu}^{\mathrm{BB}}$ is admissible in the EB setting with data generating process given by both (ref) and (ref).

We next show that the NPMLE-based EB estimator $\widehat{\boldmu}^{\mathrm{EB}}$ is inadmissible. To state this result precisely, we first define the NPMLE of the mixing distribution $G$. Let

equation[equation omitted — 201 chars of source]

where $\mathcal{P}(\Theta)$ denotes the set of all probability distributions supported on $\Theta$, with $\Theta=\RR$ in the unconstrained case and $\Theta=[-M,M]$ in the constrained case. Here, $f_G(\cdot)$ is the marginal density of $Z_i$ induced by the hierarchical model (ref)--(ref). The EB estimator for the $i$-th coordinate is the plug-in Bayes rule $\widehat{\mu}_i^{\mathrm{EB}}(\boldZ) := \delta_{\widehat{G}}(Z_i)$.

propLet $\widehat{\boldmu}^{\mathrm{EB}}$ denote the NPMLE EB estimator defined above, where $\widehat{G}$ in (ref) is either constrained to $\mathcal{P}([-M,M])$ or taken over $\mathcal{P}(\RR)$. Then $\widehat{\boldmu}^{\mathrm{EB}}$ is inadmissible.

We note that datta2025polynomial also consider admissibility questions in a closely related EB setting, but from a different perspective. They focus on F-modeling strategies efron2014two, in which one directly estimates the marginal density $f_G$ by $\hat f$ and then plugs $\hat f$ into the Eddington/Tweedie formula dyson1926method,efron2011tweedie,

equation[equation omitted — 78 chars of source]

yielding an estimator of the form $\hat{\mu}_i = Z_i + \hat{f}'(Z_i)/\hat{f}(Z_i)$. For several popular F-modeling procedures, such as the polynomial log-marginal proposal of efron2011tweedie, they ask whether there exists an implicit (data-driven) prior \smash{$\tilde{G}$} such that \smash{$\hat{\mu}_i = \delta_{\tilde{G}}(Z_i)$}. However, even when such an implicit prior \smash{$\tilde{G}$} exists, this does not by itself imply that the corresponding estimator is Bayes, or admissible, for the full decision problem.

Posterior contraction and proof of Theorem (ref)

We now describe the main results used in the proof of Theorem (ref). Our starting point is a posterior contraction result for the Dirichlet process mixture model under the data generating process (ref)--(ref). Under (ref)--(ref), there exists a true mixing distribution $G_{\star}$ (and hence a true marginal density $f_{G_{\star}}$). This object will be important for the posterior contraction results that we state below. For two densities $f$ and $h$ define their Hellinger distance as:

equation[equation omitted — 107 chars of source]

The following result is a minor strengthening of ghosal2001entropies: in addition to the contraction rate itself, it tracks the corresponding high-probability events.

theo[Posterior contraction] Suppose that Specification (ref) holds, the true prior $G_{\star}$ is supported on $[-M,M]$. Then there exist constants $C,c>0$ depending only on $(M,\alpha,\eta)$ such that for all sufficiently large $n$, $$ \PP[G_{\star}]{ \Pi\p{G\,:\, \Dhel( f_{G_{\star}}, f_G) \geq C\frac{\log n}{\sqrt{n}} \; \Big | \; \boldZ} \leq \exp\p{-c \log^2 n}} \geq 1-\frac{1}{n}. $$ Moreover, the posterior mean marginal density $\bar{f} = \int f_{G}\, \dd \Pi(G \mid \boldZ)$ satisfies $$ \PP[G_{\star}]{ \Dhel( f_{G_{\star}}, \bar{f}) \geq C' \frac{ \log n}{\sqrt{n}}} \leq \frac{1}{n}.$$

To leverage Theorem (ref) in the proof of Theorem (ref), we next connect Hellinger contraction of the marginal densities to the discrepancy between the BNP estimator $\widehat{\boldmu}^{\mathrm{BB}}$ and the Bayes oracle $\widehat{\boldmu}^{\mathrm{B}}$. To that end, note that for any priors $G$ and $Q$ we have

equation[equation omitted — 208 chars of source]

where $\Dfisher{f_G}{f_Q}$ is the Fisher divergence; see ghosh2025steins for discussion of its relationship to empirical Bayes estimation.

Next, by combining Lemma 6.1 of zhang2005general with Theorem E.1 of saha2020nonparametric (which builds on Theorem 3 of jiang2009general), we can relate Fisher divergence to Hellinger distance as follows.

lemmThere exists a universal constant $C>0$ such that for any $\rho \in (0, (2\pi e)^{-1/2})$, $M>0$ and any two distributions $G$ and $Q$ supported on $[-M,M]$, we have that $$ \Dfisher{f_{G}}{f_Q} \leq C\p{ \Dhel^2(f_G, f_Q) \cdot \max\cb{ \abs{\log \rho}^3,\; \abs{ \log \Dhel(f_G,f_Q)}} \,+\, (M+1)\rho \abs{\log \rho}}. $$

We now have all the ingredients to prove Theorem (ref). Since it is short, we give the proof here.

proof[Proof of Theorem (ref)] Fix $G_{\star}$ with $\mathrm{supp}(G_{\star})\subseteq[-M,M]$ and write $\hat\mu_i^{\mathrm B}:=\delta_{G_{\star}}(Z_i)$, so that \smash{$R(\widehat{\boldmu}^{\mathrm B},G_{\star})=\sqrt{\EE[G_{\star}]{\Var[G_{\star}]{\mu\mid Z}}}$}. By the triangle inequality, \begin{equation} R(\widehat{\boldmu}^{\mathrm{BB}},G_{\star})-\sqrt{\EE[G_{\star}]{\Var[G_{\star}]{\mu\mid Z}}} \le \left\{\frac1n\sum_{i=1}^n \EE[G_{\star}]{\big(\hat\mu_i^{\mathrm{BB}}-\hat\mu_i^{\mathrm B}\big)^2}\right\}^{1/2}. \end{equation} By Proposition (ref), \smash{$\hat\mu_i^{\mathrm{BB}}=\delta_{\widebar G(\boldZ_{-i})}(Z_i)$}. Let \smash{$\bar f_{-i}:= f_{\widebar G(\boldZ_{-i})} = \int f_{G} \,\dd\Pi(G\mid \boldZ_{-i})$}. Conditioning on $\boldZ_{-i}$ and using (ref) and (ref), we have that $$ \EEInline[G_{\star}]{(\hat\mu_i^{\mathrm{BB}}-\hat\mu_i^{\mathrm B})^2} = \EE[G_{\star}]{\EEInline[G_{\star}]{(\hat\mu_i^{\mathrm{BB}}-\hat\mu_i^{\mathrm B})^2\mid \boldZ_{-i}}} = \EE[G_{\star}]{\Dfisher{f_{G_{\star}}}{\bar f_{-i}}}. $$ Consider the event $A_i:= \cb{ \Dhel(f_{G_{\star}},\bar f_{-i})< C' \log n/\sqrt{n}}$, where $C'$ is the constant from Theorem (ref) and note that $\PP[G_{\star}]{A_i^c} < 1/n$. On the event $A_i$, Lemma (ref) with $\rho=n^{-1}$ yields \smash{$\Dfisher{f_{G_{\star}}}{\bar f_{-i}}\lesssim (\log^5 n)/n$}, while always \smash{$\Dfisher{f_{G_{\star}}}{\bar f_{-i}}\le 4M^2$} since \smash{$\hat\mu_i^{\mathrm{BB}}$} and \smash{$\hat\mu_i^{\mathrm B}$} are both in $[-M,M]$. Hence \smash{$\EE[G_{\star}]{(\hat\mu_i^{\mathrm{BB}}-\hat\mu_i^{\mathrm B})^2}\lesssim (\log^5 n)/n$} uniformly in $i \in \cb{1,\ldots,n}$. Plugging into (ref) and then taking the supremum over $G_{\star} \in \mathcal{P}([-M,M])$ completes the proof.

Compound decisions via Bayesian nonparametrics

Compound decisions regret

In this section, we suppose that data-generation is given by (ref) only and measure risk by the RMSE $R(\widehat{\boldmu},\boldmu)$ defined in (ref). In Section (ref) we previewed our main regret result, namely Theorem (ref). Below we state a slightly stronger result by benchmarking $\widehat{\boldmu}^{\mathrm{BB}}$ against an even broader class of estimators, that of all permutation equivariant decision rules weinstein2021permutation,

equation[equation omitted — 238 chars of source]

where $\mathcal{S}_n$ is the set of permutations of $\{1,\ldots,n\}$. In the absence of further information on $\boldmu$, permutation equivariant rules are natural in the setting of (ref).

theoSuppose that Specification (ref) holds for $\Pi$. Then there exists a constant $C>0$ depending only on $(M, \alpha, \eta)$ such that for all $n \geq 2$, $$\sup_{\boldmu \in [-M,M]^n} \p{ R\p{\widehat{\boldmu}^{\mathrm{BB}}(\Pi) ,\,\boldmu} - \inf_{\boldt \in \mathcal{T}^{\mathrm{PE}} } \cb{R\p{ \boldt(\boldZ) ,\boldmu} }} \leq C\frac{\log^{5/2} n}{\sqrt{n}}$$ where the infimum is over all permutation equivariant decision rules defined in (ref) and $R(\cdot,\boldmu)$ is defined as the frequentist risk in (ref) when data is generated according to (ref) with $\boldmu$.

The reason that Theorem (ref) is a strengthening of Theorem (ref) is that for any prior $\widetilde{\Pi}$, it holds that $\widehat{\boldmu}^{\mathrm{B}}(\widetilde{\Pi})$ can be written as $\boldt(\boldZ)$ for a certain $\boldt \in \mathcal{T}^{\mathrm{PE}}$. The regret bound of Theorem (ref) matches the regret bound for the NPMLE shown by jiang2009general.

Bayes properties and admissibility

We now turn to the question of admissibility in the compound setting, complementing the results of Section (ref). Here admissibility is understood in the usual sense for multivariate normal mean estimation under squared error loss. We consider two parameter spaces: the compact set $[-M,M]^n$ (which aligns with our regret bounds) and the full space $\RR^n$.

defi[Admissibility] Let $\Theta \subseteq \RR^n$ be the parameter space (either $\Theta = [-M,M]^n$ or $\Theta = \RR^n$). An estimator $\widehat{\boldmu}$ is inadmissible over $\Theta$ if there exists another estimator $\widetilde{\boldmu}$ such that: \begin{enumerate} • $R(\widetilde{\boldmu}, \boldmu) \leq R(\widehat{\boldmu}, \boldmu)$ for all $\boldmu \in \Theta$, and • $R(\widetilde{\boldmu}, \boldmu_0) < R(\widehat{\boldmu}, \boldmu_0)$ for some $\boldmu_0 \in \Theta$. \end{enumerate} Otherwise, we say that $\widehat{\boldmu}$ is admissible over $\Theta$.

Our first result is that the BNP estimator $\widehat{\boldmu}^{\mathrm{BB}}$ is admissible in the above sense.

propThe estimator $\widehat{\boldmu}^{\mathrm{BB}}$ is admissible over $\RR^n$ and over $[-M,M]^n$ in the compound decision setting with data generating process given by (ref).

By contrast, the classical NPMLE-based empirical Bayes rule is not admissible.

propLet $\widehat{\boldmu}^{\mathrm{EB}}$ be the NPMLE-based EB estimator defined above with the NPMLE $\widehat{G}$ in (ref) either constrained to distributions supported on $[-M,M]$ or unconstrained. Then: \begin{enumerate}[label=(\roman*)] • For the compact parameter space $[-M,M]^n$: $\widehat{\boldmu}^{\mathrm{EB}}$ is inadmissible for all $n \geq 1$, regardless of whether $\widehat{G}$ is constrained or unconstrained. • For the full parameter space $\RR^n$: $\widehat{\boldmu}^{\mathrm{EB}}$ is inadmissible for $n \geq 2$ when $\widehat{G}$ is unconstrained, and for all $n \geq 1$ when $\widehat{G}$ is constrained to $[-M,M]$. \end{enumerate}

remaFor $n=1$, the unconstrained NPMLE is admissible over $\RR$. The reason is that the NPMLE \smash{$\widehat{G}$} is a point mass at $Z_1$, and so \smash{$\widehat{\mu}_1^{\mathrm{EB}} \equiv Z_1$} which is admissible for $n=1$.

proof[Proof sketch.] There are four inadmissibility claims above. Let us sketch the proof of the inadmissibility of the unconstrained NPMLE-based EB estimator over $\RR^n$ for $n \geq 2$. The key idea is to observe the following. Fix $n \geq 2$. For $\boldz\in \RR^n$, let $\bar{z} := n^{-1}\sum_{i=1}^n z_i$ and define the set \begin{equation} \mathcal U := \left\{\boldz\in\RR^n:\ \max_{1\le i\le n}|z_i-\bar z| \leq 1 \right\}. \end{equation} By the first order optimality condition of the NPMLE in (ref) over $\mathcal{P}(\RR)$, we have that $\widehat{G}(\boldz) = \delta_{\bar z}$ for all $\boldz \in \mathcal{U}$. Thus, $$ \widehat{\boldmu}^{\mathrm{EB}}(\boldz) = \bar{z} \cdot \mathbf{1}_n = (\bar{z},\ldots,\bar{z}) \quad \text{for all } \boldz \in \mathcal{U}. $$ Now suppose that \smash{$\widehat{\boldmu}^{\mathrm{EB}}$} is admissible. Then by Theorem 3.1.1. of brown1971admissible, \smash{$\widehat{\boldmu}^{\mathrm{EB}}$} is a generalized Bayes estimator and so it is analytic as a function of $\boldz$. Hence, \smash{$\widehat{\boldmu}^{\mathrm{EB}}(\boldz) = \bar{z} \cdot \mathbf{1}_n$} for almost all $\boldz \in \RR^n$ by the identity theorem. However, this can easily be shown to be a contradiction. Thus, $\widehat{\boldmu}^{\mathrm{EB}}$ is inadmissible.

remaThe proof technique above has been used before to show inadmissibility of estimators in the multivariate normal mean estimation problem. Recall the James-Stein james1961estimation estimator and its positive-part, $$ \widehat{\boldmu}^{\mathrm{JS}} := \left(1 - \frac{(n-2)}{\Norm{\boldZ}^2_2}\right) \boldZ,\;\;\;\;\;\;\; \widehat{\boldmu}^{\mathrm{JS+}} := \left(1 - \frac{(n-2)}{\Norm{\boldZ}^2_2}\right)_+ \boldZ, $$ where $(x)_+ = \max\{x,0\}$. For $n\geq3$, the James-Stein estimator dominates the MLE $\boldZ$, and thus the latter is inadmissible. However, \smash{$\widehat{\boldmu}^{\mathrm{JS}}$} is itself inadmissible since it is dominated by \smash{$\widehat{\boldmu}^{\mathrm{JS+}}$} efron1973stein, and \smash{$\widehat{\boldmu}^{\mathrm{JS+}}$} is inadmissible as well. To show inadmissibility of the latter, brown1986fundamentals observes that \smash{$\widehat{\boldmu}^{\mathrm{JS+}}(\boldz) =0$} for \smash{$\Norm{\boldz}^2_2 \leq n-2$}. By the same continuation argument as above, if \smash{$\widehat{\boldmu}^{\mathrm{JS+}}$} were admissible, then \smash{$\widehat{\boldmu}^{\mathrm{JS+}}(\boldz) = 0$} for almost all $\boldz \in \RR^n$, which is a contradiction. Thus, \smash{$\widehat{\boldmu}^{\mathrm{JS+}}$} is inadmissible.

Posterior contraction and proof of Theorem (ref)

The overall proof strategy has four main elements, each of which we outline in turn:

enumerate[noitemsep,leftmargin=*] • reduction to separable decision rules; • the fundamental theorem of compound decision; • posterior contraction under (ref); • control of posterior means following jiang2009general.

To begin, we define the class of simple separable decision rules:

equation[equation omitted — 192 chars of source]

As a preliminary step, we note the following quantitative bound, due to greenshtein2009asymptotic, that sharpens earlier results of hannan1955asymptotic; see also han2025besting and liang2025sharp.

theoThere exists a constant $C_M >0$ that depends only on $M$ such that $$\sup_{\boldmu \in [-M,M]^n} \p{ \inf_{\boldt \in \mathcal{T}^{\mathrm{S}} } \cb{R(\boldmu, \boldt(\boldZ)) } - \inf_{\boldt \in \mathcal{T}^{\mathrm{PE}} } \cb{R(\boldmu, \boldt(\boldZ)) }} \leq C_M \frac{1}{\sqrt{n}}.$$

The upshot of Theorem (ref) is that it suffices to compete with the best simple separable estimator, rather than with the best permutation equivariant estimator. Moreover, the optimal simple separable rule admits an explicit characterization via the fundamental theorem of compound decision theory, which will be useful in later steps of the proof.

A key idea in compound decision theory is that the compound model with fixed $\boldmu$ in (ref) behaves, in certain respects, similarly to the Bayes model in (ref)--(ref) with \smash{$\mu_i \simiid G_n$}, where $G_n$ is the empirical distribution of $\boldmu$,

equation[equation omitted — 109 chars of source]

and where $\delta_u$ denotes the Dirac mass at $u$. This similarity also underlies the result of Theorem (ref) above and can be formalized by comparing the marginal distribution of a randomly permuted version of $\boldZ$ under the two models; see, e.g., han2025approximate. As we show below, many useful consequences of this connection follow already from the linearity of expectation. We illustrate this by giving a self-contained proof of the fundamental theorem of compound decisions robbins1951asymptotically, zhang2003compound.

theo[Fundamental theorem of compound decisions] Let $t:\RR\to\RR$ be measurable and consider the separable estimator with $\hat\mu_i=t(Z_i)$ for all $i$. Then under (ref), \begin{equation} R( \widehat{\boldmu},\,\boldmu) =R(\widehat{\boldmu},\, G_n(\boldmu)), \end{equation} where the left-hand side risk is defined as in (ref) and the right-hand side risk is defined in (ref). Consequently, $t^\star_{\boldmu}(z)=\delta_{G_n(\boldmu)}(z)$ yields the separable estimator that minimizes $R( \boldt(\boldZ),\,\boldmu)$ over all $\boldt \in \mathcal{T}^{\mathrm{S}}$.

proofLet $\widehat{\boldmu}$ be a separable estimator as above with coordinate-wise function $t(\cdot)$. Writing $G_n = G_n(\boldmu)$ as in (ref), then we can write $R^2(\widehat{\boldmu},\,\boldmu)$ as follows, $$ \frac{1}{n}\sum_{i=1}^n \EE[\boldmu]{(t(Z_i)-\mu_i)^2}= \frac{1}{n}\sum_{i=1}^n \int (t(z)-\mu_i)^2\,\varphi(z-\mu_i)\,\dd z = \int\!\int (t(z)-\nu)^2\,\varphi(z-\nu)\,\dd z\,\dd G_n(\nu). $$ Analogously, we can write $R^2( \widehat{\boldmu},\,G_n)$ as $$ \frac{1}{n} \sum_{i=1}^n \EE[G_n]{(t(Z_i)-\mu_i)^2} = \EE[G_n]{(t(Z_i)-\mu_i)^2} = \int\!\int (t(z)-\nu)^2\,\varphi(z-\mu)\,\dd z\,\dd G_n(\nu). $$ Thus indeed $R(\widehat{\boldmu},\,\boldmu) = R( \widehat{\boldmu},\,G_n)$. As explained before (ref), $R( \widehat{\boldmu},\,G_n)$ is minimized over all estimators $\widehat{\boldmu}$ by $\widehat{\boldmu}^{\mathrm{B}}(G_n)$, which is in fact a separable estimator with coordinate-wise function $\delta_{G_n}(\cdot)$ (also defined in (ref)). Thus the best separable estimator is also the same, and also minimizes the equivalent objective $R(\boldt(\boldZ), \boldmu)$ over all separable estimators.

Given Theorem (ref) and the results of Section (ref), it is natural to expect that the BNP posterior $\Pi(\cdot \mid \boldZ)$ may concentrate around the empirical mixing distribution $G_n$ in (ref). The next result makes this intuition precise. To that end, define the marginal density induced by $G_n$ (equivalently, by (ref) with $\boldmu$ fixed) as

equation[equation omitted — 76 chars of source]

A result similar to the next one (but without a rate) was shown by datta1991consistency.

theo[Posterior contraction] Assume $\max_{1\le i\le n}|\mu_i|\le M$ and suppose Specification (ref) holds. Then there exist constants $C,c>0$ depending only on $(M,\alpha,\eta)$ such that for all sufficiently large $n$, $$ \PP[\boldmu]{ \Pi\p{G\,:\, \Dhel( f_{\boldmu}, f_G) \geq C\frac{\log n}{\sqrt{n}} \; \Big | \; \boldZ} \leq \exp\p{-c \log^2 n}} \geq 1-\frac{1}{n}. $$

proof[Proof sketch.] Fix $C>0$ (we will choose it later), let $\varepsilon_n := \log n/\sqrt{n}$, and consider the set $ \mathcal{U} := \cb{ G \,:\, \Dhel(f_{\boldmu}, f_G) \geq C\varepsilon_n }. $ By the simplest Schwartz posterior contraction argument, \begin{equation} \begin{aligned} \Pi( \mathcal{U} \mid \boldZ) \;= \; \frac{ \int_{\mathcal{U} } \prod_{i=1}^n f_G(Z_i) \dd\Pi(G)}{\int \prod_{i=1}^n f_G(Z_i)\dd\Pi(G)} \;=\; \frac{ \int_{\mathcal{U} } \prod_{i=1}^n \frac{f_G(Z_i)}{f_{\boldmu}(Z_i)}\dd\Pi(G)}{\int \prod_{i=1}^n \frac{f_G(Z_i)}{f_{\boldmu}(Z_i)}\dd\Pi(G)}. \end{aligned} \end{equation} We have to show two facts: that with high probability under $\PP[\boldmu]{\cdot}$, both the numerator is “small” and the denominator is “large.” Note that $\prod_{i=1}^n f_{\boldmu}(Z_i) \neq \prod_{i=1}^n \varphi(Z_i-\mu_i),$ where the latter is the likelihood under the data generating process in (ref), and so the ratio introduced in (ref) is not a likelihood ratio under $\PP[\boldmu]{\cdot}$. We can show that the posterior contracts despite this misspecification.\footnote{Similar arguments also appear in the analysis of nonparametric Bayesian procedures that utilize quasi-likelihoods kato2013quasi,kankanala2025generalized.} To see why, consider the following argument. Fix $f_G$ for $G \in \mathcal{P}([-M,M])$. To understand both the numerator and denominator in (ref), we must understand the behaviour of the log-likelihood ratio $\sum_{i=1}^n \log\cb{f_G(Z_i)/f_{\boldmu}(Z_i)}$. Taking expectations under $\PP[\boldmu]{\cdot}$, we get $$ \begin{aligned} \EE[\boldmu]{ \sum_{i=1}^n \log\p{\frac{f_G(Z_i)}{f_{\boldmu}(Z_i)}} } &= \sum_{i=1}^n \int \log\p{\frac{f_G(z)}{f_{\boldmu}(z)}} \varphi(z-\mu_i)\, \dd z \\ &= n \int \log\p{\frac{f_G(z)}{f_{\boldmu}(z)}} \underbrace{\frac{1}{n}\sum_{i=1}^n \varphi(z-\mu_i)}_{=f_{\boldmu}(z)}\, \dd z \,=\, - n \DKL{f_{\boldmu}}{f_G}, \end{aligned} $$ where $\DKL{\cdot}{\cdot}$ is the Kullback-Leibler divergence defined as follows for two densities $f,h$, \begin{equation} \DKL{f}{h} := \int f(z) \log\frac{f(z)}{h(z)} \dd z. \end{equation} By comparison, if data generation were given by both (ref) and (ref) with true prior $G_{\star}$ (as in the setting of Theorem (ref)), then we would have $\EE[G_{\star}]{ \sum_{i=1}^n \log\cb{f_G(Z_i)/f_{G_{\star}}(Z_i)}} = - n \DKL{f_{G_\star}}{f_{G}}$, i.e., effectively $f_{\boldmu}$ plays the role of $f_{G_\star}$. Given suitable concentration uniformly over $f_G$, we can now surmise the following. In the numerator of (ref), fix $G \in \mathcal{U}$ (within the integral). Note that $-n\DKL{f_{\boldmu}}{f_G} \leq -2n\Dhel^2(f_{\boldmu}, f_G) \leq -2 C^2 n \varepsilon_n^2$, and so the numerator should be small. Lemma (ref) in the supplement makes this argument precise. The argument closely tracks the arguments showing Hellinger rates for the NPMLE by zhang2009generalized under (ref). Similarly, for the denominator, as long as the prior $\Pi$ puts enough mass around $f_{\boldmu}$ in KL neighborhoods, the denominator will be sufficiently large with high probability. Lemma (ref) in the supplement makes this argument precise and closely tracks the arguments of ghosal2001entropies under the standard frequentist BNP setup where both (ref) and (ref) determine the data generating process. The key difference here is that the “true” density $f_{\boldmu}$ depends on $\boldmu$ in a non-iid way, but this does not impact the prior mass calculations.

The fourth (and last) step of the proof amounts to controlling $\EEInline[\boldmu]{ \NormInline{ \widehat{\boldmu}^{\mathrm{B}}(G_n(\boldmu)) - \widehat{\boldmu}^{\mathrm{B}}}^2}$. This term can be controlled using Theorem (ref), Lemma (ref) as well as results of jiang2009general on regularized Bayes rules. Supplement (ref) explains this step and puts together the overall proof of Theorem (ref).

I. J. Good's staircase

The hierarchy set forth through (ref), (ref), and (ref) can be continued further, say as,

eqbox\begin{equation} \tag{BBB} \Pi \; \sim \; \Gamma. \end{equation}

where $\Gamma$ is a hyperhyperprior on $\Pi$, e.g., on the parameters $\alpha$ and/or $H$ of the Dirichlet process and so forth. At each stage of the hierarchy, say, at (ref), we have three options: (i) fix that level, treating the distribution at that level as known, (ii) estimate the distribution at that level by empirical Bayes (EB), or (iii) proceed to the next level of the hierarchy. good1992bayes uses the terminology “B”, “EB”, “BB”, “EBB”, “BBB”, etc., to denote these various options.

Thus, standard EB as introduced by robbins1956empirical is also called EB in Good's staircase notation. Meanwhile, EB that estimates \smash{$\widehat{\Pi}$} (e.g., \smash{$\widehat{\alpha}$} and/or \smash{$\widehat{H}$ when $\Pi=\mathrm{DP}(\alpha,H)$}) in (ref) as in e.g., liu1996nonparametric, mcauliffe2006nonparametric, donnet2018posteriora is called EBB. Fully Bayesian nonparametric approaches that put hyperpriors on $\alpha$ and/or $H$ as in e.g., escobar1995bayesian are called BBB. Although imperfect and not always applicable,\footnote{The hierarchy is not always clear-cut. For instance, one could treat (ref) and (ref) as a single level given by the composition of the two. Or if one is not interested in the $\mu_i$ in (ref) but only in the induced marginal densities $f_G(\cdot)$, then it makes sense to merge (ref) with (ref) and call (ref) as BB instead of BBB. } we find this terminology useful in communicating the various levels of hierarchy and estimation.

In our setting, we could analyze any of these schemes under (ref) (as in Section (ref)), or under both (ref) and (ref) (as in Section (ref)). For instance, an EBB scheme could be analyzed by using techniques in e.g., petrone2014empirical, rousseau2017asymptotic, donnet2018posteriora.

We also refer to vandervaart2023frequentism and ignatiadis2025partially for further discussion on the various levels of the hierarchy and how frequentist guarantees can be obtained by treating randomness up to a certain level of the hierarchy in a frequentist manner.

Numerical results

Remarks on implementation

Throughout, to compute the BNP estimator of $\boldmu$, we make the following choices for the hyperparameters of the Dirichlet Process in (ref): $$H=\mathrm{Unif}[-10,\,10],\;\;\;\alpha \sim \Gamma(0.01,100),$$ where $\Gamma(a,b)$ is the Gamma distribution with shape $a$ and scale $b$. We use the Gibbs sampler of Algorithm 2 in neal2000markov along with updates for $\alpha$ described in escobar1995bayesian.\footnote{In writing our initial prototype, we followed the code in the Particles.jl package kleinschmidt2024particles.} The reason we can directly use the Gibbs sampler of neal2000markov is that we can explicitly conduct the required marginalization and posterior updates in closed form, see Supplement (ref). We use 2,000 burn-in iterations and 10,000 post burn-in iterations (with iterations defined as full sweeps through the data) to compute the posterior mean of $\boldmu$.

For the NPMLE, we use the by-now-standard approach of koenker2014convex using discretization and the MOSEK mosek interior point convex optimization solver

Simulation study

table[table omitted — 970 chars of source]

table[table omitted — 902 chars of source]

We consider a standard sparse simulation setup johnstone2004needles, jiang2009general, koenker2014convex under (ref). We take $n=1,000$ throughout. In each simulation setting, we have that, $$ \mu_i =

cases\mu, for i=1,\ldots, n_1\\ 0, for i=n_1+1,\ldots,n,

$$ where $n_1 \in \cb{5,50,500}$ and $\mu \in \cb{3,4,5,7}$ are simulation parameters.

We consider three methods: the proposed BNP estimator, the NPMLE (as described in Section (ref)), and the separable oracle from Theorem (ref). Following standard practice, we report results in terms of unnormalized MSE, $\EEInline[\boldmu]{\NormInline{\widehat{\boldmu}-\boldmu}^2}$ which here we estimate by averaging over $100$ Monte Carlo replicates. Results are shown in Table (ref). Moreover since the BNP approach to EB has been criticized for being computationally slow, we include timings in Table (ref).

We observe that the BNP estimator performs competitively with the NPMLE across all configurations, and both methods track the oracle well. For well-separated signals ($\mu=7$), the BNP estimator is closer to the oracle than the NPMLE in every sparsity regime. The computational cost is higher, roughly a five-to tenfold increase in wall time, but remains moderate in absolute terms for $n=1{,}000$, suggesting that the BNP approach is a practical alternative to the NPMLE

\paragraph{Acknowledgements.} We thank Asaf Weinstein for helpful discussions on Bayes empirical Bayes. This work was completed in part with resources provided by the University of Chicago’s Research Computing Center. N.I. gratefully acknowledges support from the U.S. National Science Foundation (DMS-2443410).