Quasi-Bayes empirical Bayes estimation of sums of random variables

Technical Summary: Quasi-Bayes Empirical Bayes Estimation of Sums of Random Variables

Problem Setting
The paper addresses the long-standing statistical problem of estimating sums of functions of observable and unobservable variables, denoted as $S_n = \sum_{i=1}^n u(X_i, \theta_i)$. Here, $X_i$ are observable random variables, while $\theta_i$ are unobservable parameters (or constants) and $u$ is a utility function. This framework encompasses diverse applications, including treatment effects in clinical trials, species sampling (estimating unseen species or missing mass), data confidentiality (disclosure risk measures), and network analysis.

The underlying model is a mixture model where $X_i | \theta_i \sim k(\cdot | \theta_i)$ and $\theta_i \sim G$, with $G$ being an unknown mixing distribution. While Empirical Bayes (EB) methods offer a natural framework for this task, existing approaches face significant limitations:

• Parametric g-modeling: Assumes $G$ belongs to a specific parametric family. While computationally efficient and asymptotically efficient under correct specification, it is vulnerable to model misspecification.
• Nonparametric "u, v" method: A flexible nonparametric approach (Robbins, 1988) that is asymptotically efficient but is restricted to a narrow class of utility functions $u$ for which explicit solutions to integral equations exist.
• Nonparametric g-modeling: General nonparametric estimation of $G$ (e.g., via maximum likelihood) lacks established theoretical guarantees for the resulting plug-in estimates of $S_n$ and is often computationally demanding due to the infinite-dimensional nature of the problem.

Methodology: Quasi-Bayes EB
The authors propose a nonparametric g-modeling approach termed Quasi-Bayes Empirical Bayes. The core of the methodology is a recursive estimation of the mixing distribution $G$ using Newton's algorithm (Smith and Makov, 1978; Newton et al., 1998).

1. Recursive Learning: Starting from an initial guess $G_0$, the mixing distribution is updated sequentially as new observations $X_n$ arrive. The update rule is:
   $$G_{n+1}(d\theta) = (1 - \alpha_{n+1})G_n(d\theta) + \alpha_{n+1} \frac{k(X_{n+1} | \theta)G_n(d\theta)}{\int_\Theta k(X_{n+1} | \theta)G_n(d\theta)}$$
   where $(\alpha_n)$ is a learning rate sequence satisfying standard stochastic approximation conditions ($\sum \alpha_n = \infty, \sum \alpha_n^2 < \infty$).

2. Quasi-Bayesian Interpretation: This recursive scheme is interpreted as a learning model rather than a data-generating process. It approximates the Bayesian learning of latent parameters $\theta_i$ by replacing the unknown true mixing distribution in posterior expectations with the recursive estimate $G_n$. Specifically, the conditional distribution of $\theta_i$ given $X_i$ is approximated by $G_n(\cdot | X_i)$.

3. Estimation: The Quasi-Bayes EB estimate of the sum is the plug-in estimator:
   $$\hat{S}_n^{[Q-B]} = \sum_{i=1}^n \int_\Theta u(X_i, \theta) G_n(d\theta | X_i)$$
   This estimate is computationally efficient, incurring a constant per-observation cost dependent only on the grid resolution used for numerical integration.

4. Uncertainty Quantification: The authors derive asymptotic credible intervals for $S_n$ based on a conditional Central Limit Theorem (CLT) applied to the sum of independent random variables under the quasi-Bayesian learning process.

Key Theoretical Contributions
The paper establishes large-sample asymptotic guarantees for the proposed method:

• Merging of Distributions: The authors prove that as $n \to \infty$, the probability measure induced by the quasi-Bayesian learning process ($P^{(n)}$) merges with the corresponding Bayesian learning process ($P$). Specifically, the conditional distributions of the latent parameters $(\theta_1, \dots, \theta_m)$ given $X_{1:n}$ converge in Prohorov distance. This rigorously justifies interpreting $\hat{S}_n^{[Q-B]}$ as an approximate Bayesian estimate in large samples.
• Frequentist Consistency: Under the assumption that data are generated from a "true" mixing distribution $G^*$, the paper proves that the Quasi-Bayes EB estimate is consistent relative to the oracle estimate $\hat{S}_n(G^*)$. Specifically, $n^{-1}|\hat{S}_n^{[Q-B]} - \hat{S}_n(G^*)| \to 0$ almost surely.
• Convergence Rates: Under additional assumptions (specifically, discretizing the parameter space onto a finite grid), the paper establishes convergence rates for the estimation error, showing it is $o(n^\delta)$ for specific $\delta$ depending on the learning rate.

Empirical Results
The methodology is validated through synthetic and real-data analyses, primarily focusing on Poisson and Gaussian mixture models.

• Synthetic Data: Experiments compare the Quasi-Bayes EB estimate against parametric EB/Bayes EB estimates and the nonparametric "u, v" method.
  • For utility functions where the "u, v" method is applicable (e.g., estimating total intensity for specific accident counts), the Quasi-Bayes EB approach achieves accuracy comparable to the "u, v" method and superior to parametric methods when the parametric assumption is misspecified.
  • For utility functions where the "u, v" method is inapplicable (e.g., estimating the number of motorists whose accidents exceed their intensity), the Quasi-Bayes EB method outperforms parametric alternatives and provides a viable nonparametric solution.
  • The method demonstrates robustness to grid resolution choices and computational scalability.
• Real Data: The approach is applied to National Hockey League (NHL) goal-scoring data (2017–2019 seasons) to estimate future scoring contributions and the number of players with declining performance. The Quasi-Bayes EB estimates show stable out-of-sample performance, often yielding lower absolute deviations than parametric competitors, particularly for complex functionals where "u, v" estimates are unavailable.

Significance and Claims
The paper claims to bridge the gap between the flexibility of nonparametric methods and the computational tractability of parametric EB. Its primary significance lies in:

1. Generality: It provides a unified framework applicable to a broad class of utility functions $u$, overcoming the restrictive solvability conditions of the "u, v" method.
2. Computational Efficiency: Unlike full nonparametric Bayesian approaches (e.g., Dirichlet processes requiring MCMC) or nonparametric maximum likelihood (which can be computationally intensive), the recursive Newton-based update is scalable and linear in the number of observations.
3. Theoretical Rigor: It provides the first theoretical guarantees (merging with Bayes, frequentist consistency, and convergence rates) for nonparametric g-modeling estimates of sums of random variables.
4. Uncertainty Quantification: It offers a principled way to construct asymptotic credible intervals for the target sum without relying on full Bayesian posterior sampling.

The authors acknowledge limitations, noting that convergence rate proofs currently rely on discretizing the parameter space and that the frequentist coverage properties of the derived credible intervals remain an open question. They also identify future directions, including multidimensional extensions and applications to species sampling problems where "u, v" estimates are unavailable.