We study function estimation in the empirical Bayes setting for Poisson and normal means. Specifically, given observations $Y_i\sim f(\cdot; θ_i)$ with latent parameters $θ_i\sim π$, the goal is to estimate $\mathbb{E}_π[\ell(θ)|X = x]$. This task lies between classical deconvolution (recovering the full prior $π$), and standard empirical Bayes mean estimation. While the minimax risk for estimating $π$ in the Wasserstein distance is known to decay only logarithmically, we show that estimating certain smooth functions admits dramatically faster rates. In particular, for polynomial functions of degree $k$ in the Poisson model, we establish a tight bound of $Θ(\frac{1}{n}(\frac{\log n}{\log \log n})^{k+1})$ and $Θ(\frac{1}{n}(\log n)^{2k+1})$ for bounded and subexponential priors, respectively, attainable by estimators mimicking those that achieve optimal regret for the mean estimation problem (Robbins, mininum distance, ERM). Our analysis identifies the approximation-theoretic origin of this improvement: smooth functions can be well-approximated by low-degree polynomials, whereas Lipschitz functions require dense polynomial approximations, incurring a $\frac{1}{k}$ loss for degree $k$ polynomial approximation. The results reveal a sharp hierarchy in the difficulty of empirical Bayes problems: ranging from slow, logarithmic deconvolution to near-parametric convergence for smooth posterior functionals, and establish new connections between nonparametric empirical Bayes theory, polynomial approximation, and statistical inverse problems. Finally, we complement our analysis with a lower bound of $Ω(\frac 1n (\frac{\log n}{\log \log n})^{k+1})$ (bounded priors) and $Ω(\frac 1n (\log n)^{k + 1})$ (subgaussian priors) for the normal means model.