Simple Averaging (ZO-Avg) in Zeroth-Order Optimization

• Simple Averaging (ZO-Avg) is a zeroth-order variance-reduction method that averages independent gradient surrogates computed from function evaluations.
• It reduces estimator variance by a factor of 1/m while increasing per-iteration query costs, requiring a careful trade-off under a fixed query budget.
• The approach highlights that, for fixed resources, using a single query per iteration can be optimal, emphasizing practical and computational implications.

Simple Averaging (ZO-Avg) is a zeroth-order variance-reduction scheme in which multiple independent gradient surrogates, constructed solely from function evaluations, are averaged before or during optimization. In its most explicit Euclidean form, ZO-Avg draws random directions $u_1,\dots,u_m$ and replaces a single zeroth-order estimator by their arithmetic mean at the current iterate, thereby reducing estimator variance at the cost of additional oracle queries (Di et al., 2024, Lin et al., 19 Sep 2025). Closely related constructions also appear under other names—notably gradient sample averaging, observation averaging, consensus averaging, moving-average estimators, and transported uniform averages on manifolds—so the literature does not attach a single invariant meaning to the label across all zeroth-order settings.

1. Terminological scope and canonical definitions

In the fixed-point Euclidean setting studied most directly, ZO-Avg means averaging $m$ independent zeroth-order gradient estimates at the same point $x$. One representative definition is

$$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$

with independent Gaussian directions $u_j\sim \mathcal N(0,I)$ (Di et al., 2024). A parallel formulation writes

$$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$

where $U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$ (Lin et al., 19 Sep 2025). In that line of work, ZO-Avg is specifically a multi-query aggregation rule.

Two estimator families recur in this literature. For a smoothing radius $\mu_s>0$, the one-point finite-difference estimator is

$\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$

while the two-point version is

$$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$

The fixed-budget conclusions for ZO-Avg are stated to hold for either one- or two-point estimators, because averaging $m$0 samples yields $m$1 variance reduction while reducing the iteration count to $m$2 when the total query budget is $m$3 (Lin et al., 19 Sep 2025).

The terminology broadens in adjacent literatures. In zeroth-order online ADMM, simple averaging appears as “gradient sample averaging” and “observation averaging,” both of which average multiple independent zeroth-order samples within an iteration (Liu et al., 2017). In distributed stochastic zeroth-order optimization, the closest explicit analogue is a consensus-averaging primal scheme that averages neighboring iterates rather than multiple directions at a single point (Yi et al., 2021). In Riemannian zeroth-order optimization, a baseline called “Simple Averaging (ZO-Avg)” is not introduced explicitly; instead, the paper states that a natural variant is uniform averaging of transported past zeroth-order estimators (Li et al., 2023). Another recent line of work uses “averaging” both for averaging $m$4 Gaussian oracle calls within each iteration and for averaging the stochastic ZO trajectory in expectation, yielding a perturbed first-order recursion (Farzin et al., 28 Apr 2026).

Setting	Averaged object	Representative form
Euclidean multi-query ZO	$m$5 direction estimators at fixed $m$6	$m$7
ZOO-ADMM minibatching	directions, observations, or both	$m$8, $m$9
Distributed ZO	neighbor iterates	$x$0
Riemannian ZO	transported past estimators	$x$1

2. Statistical mechanism: unbiasedness, smoothing bias, and variance structure

In the idealized Euclidean analysis of multi-query ZO-Avg, the averaged estimator is unbiased: $x$2 and its mean-squared error satisfies

$x$3

This is the basic $x$4 variance-reduction law: averaging $x$5 independent queries divides the estimator variance by $x$6 (Lin et al., 19 Sep 2025).

A more detailed composite-optimization analysis separates sampling variance from coordinate-wise variance. In that setting,

$x$7

with each $x$8 $x$9-smooth and $$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$0-strongly convex, and $$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$1 convex and proximable (Di et al., 2024). The paper uses the one-sided Gaussian-smoothing estimator

$$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$2

and shows that its expectation differs from $$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$3 by a smoothing-bias term proportional to $$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$4 and Gaussian moments growing polynomially in $$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$5. The same analysis states that

$$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$6

so both sampling and coordinate-wise variances are reduced by a factor $$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$7 under independence (Di et al., 2024).

The coordinate-wise term is central in zeroth-order settings. The paper gives bounds of the form

$$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$8

and for a single component,

$$\hat g_{\mathrm{avg}}(x)=\frac{1}{m}\sum_{j=1}^m \hat \nabla f(x;u_j),$$9

These bounds make explicit that random-direction noise scales with dimension and with gradient magnitude (Di et al., 2024).

In composite proximal problems, the nonvanishing-gradient issue is sharper because optimality is

$u_j\sim \mathcal N(0,I)$0

so $u_j\sim \mathcal N(0,I)$1 in general. The paper therefore argues that direction-induced variance need not disappear at the solution; as a consequence, naive ZO-Avg may require large $u_j\sim \mathcal N(0,I)$2 to suppress persistent noise, which directly raises per-iteration query cost (Di et al., 2024).

The same variance-reduction logic appears in zeroth-order ADMM minibatching. There, the base one-point estimator is

$u_j\sim \mathcal N(0,I)$3

with $u_j\sim \mathcal N(0,I)$4 sampled from a distribution satisfying $u_j\sim \mathcal N(0,I)$5. Averaging $u_j\sim \mathcal N(0,I)$6 directions, $u_j\sim \mathcal N(0,I)$7 observations, or $u_j\sim \mathcal N(0,I)$8 direction-observation pairs reduces the variance term in the regret bound, replacing the dimension-dependent factor $u_j\sim \mathcal N(0,I)$9 by $$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$0 and improving the convex-rate expression to $$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$1 (Liu et al., 2017).

3. Fixed query budgets and the multi-query paradox

A major recent result is that the variance reduction of ZO-Avg does not automatically translate into better optimization efficiency once oracle cost is accounted for. Under a fixed total query budget $$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$2, using $$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$3 queries per iteration leaves only

$$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$4

iterations. The paper analyzes this trade-off for Gaussian directions across strongly convex, convex, non-convex, and stochastic convex settings and proves that, for ZO-Avg, allocating more than one query per iteration is always query-inefficient under the stated assumptions (Lin et al., 19 Sep 2025).

In the strongly convex case, with stepsize

$$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$5

the bound is

$$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$6

Substituting $$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$7 yields

$$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$8

so the exponent is best at $$\hat g_m(x)=\frac{1}{m}\sum_{i=1}^m \hat g(x;u_i) \approx \left(\frac{1}{m}UU^\top\right)\nabla f(x),$$9. The paper therefore identifies the single-query allocation as optimal under fixed $U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$0 (Lin et al., 19 Sep 2025).

The same monotonicity appears in convex and non-convex regimes. For convex $U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$1-smooth objectives, the paper gives

$U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$2

which is minimized at $U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$3. For non-convex $U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$4-smooth objectives,

$U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$5

again minimized at $U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$6 (Lin et al., 19 Sep 2025).

The stochastic convex result is formulated for a weighted-average iterate $U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$7 under a diminishing stepsize $U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$8, with $U=[u_1,\dots,u_m]\in\mathbb R^{d\times m}$9. The stated asymptotics show that both the numerator and denominator in the bound worsen as $\mu_s>0$0 increases, while the iteration count shrinks to $\mu_s>0$1, leading once more to the single-query optimum (Lin et al., 19 Sep 2025).

The scope of this “always” claim is explicitly limited. The paper states it under $\mu_s>0$2-smoothness, Gaussian or isotropic direction sampling, negligible finite-difference bias in the idealized analysis, and the theorem-specific stepsize choices. Within that framework, the conclusion is categorical: the de facto multi-query ZO-Avg rule is dominated, in query efficiency, by its $\mu_s>0$3 specialization (Lin et al., 19 Sep 2025).

4. Composite proximal optimization and the limits of naive averaging

The fixed-budget paradox does not eliminate the utility of averaging in all zeroth-order settings; rather, it clarifies that where and how averaging is inserted matters. In composite finite-sum optimization with a nonsmooth proximable regularizer, ZO-Avg reduces both sampling and coordinate-wise variances by $\mu_s>0$4, but it pays $\mu_s>0$5 stochastic zeroth-order oracle calls per iteration. Because $\mu_s>0$6 in general at a proximal optimum, the coordinate-wise variance can remain nonvanishing near the solution, so achieving high accuracy with ZO-Avg may require large $\mu_s>0$7, which the paper describes as costly (Di et al., 2024).

This observation motivates a distinct variance-reduction strategy in the same paper. Zeroth-order Proximal Double Variance Reduction (ZPDVR) combines an SVRG/SPIDER-like control variate

$\mu_s>0$8

with a temporal averaging recursion for an auxiliary estimator $\mu_s>0$9,

$\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$0

The scaling $\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$1 is tied to the Gaussian identity $\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$2, and the paper proves a contraction for $\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$3 showing progressive reduction of coordinate-wise variance over time (Di et al., 2024).

The resulting complexity statement is

$\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$4

stochastic zeroth-order oracle queries in the strongly convex smooth setting, with expected $\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$5 oracle calls per iteration when the snapshot probability is chosen as $\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$6. The point of comparison is explicit: naive ZO-Avg spends $\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$7 fresh queries every iteration, whereas ZPDVR amortizes averaging across time and suppresses both sampling and direction-induced variance without $\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$8 coordinate-wise finite differences (Di et al., 2024).

The empirical comparison in that paper reinforces the theoretical distinction. On $\hat g(x;u)=\frac{f(x+\mu_s u)-f(x)}{\mu_s}\,u,$9-regularized logistic regression for a9a, w8a, covtype, and gisette, ZPDVR exhibits linear convergence and the best performance in stochastic zeroth-order complexity plots, while ZPSVRG with random-direction gradients only stalls near the optimum, including stagnation around $$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$0 on gisette, which the paper attributes to persistent coordinate-wise variance (Di et al., 2024).

5. Minibatching, online ADMM, and decentralized averaging

In zeroth-order online ADMM, simple averaging appears not as a fixed-budget paradox but as a minibatching device embedded in a structured primal-dual algorithm. The base estimator is a one-point randomized finite difference,

$$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$1

and the paper introduces three averaging schemes: gradient sample averaging over $$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$2 directions, observation averaging over $$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$3 samples $$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$4, and hybrid averaging over $$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$5 directions and $$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$6 observations (Liu et al., 2017).

The structural role of ZO-Avg in that framework is straightforward. The ADMM $$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$7-update uses $$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$8, the $$\hat g^{\mathrm{2p}}(x;u)=\frac{f(x+\mu_s u)-f(x-\mu_s u)}{2\mu_s}\,u.$$9-update is a proximal step for the regularizer, and the dual variable is updated as usual. Averaging reduces the variance term in the regret analysis, leading to an improved convex rate

$m$00

where $m$01 denotes the dimension of the optimization variable. When $m$02, the dimension dependence is mitigated and the rate matches the first-order $m$03 form. In the strongly convex case, the single-sample ZOO-ADMM rate is $m$04 (Liu et al., 2017).

The same paper also reports application-level consequences. In black-box model blending on Yahoo! Music, both gradient sample averaging and observation averaging accelerate convergence, and with minibatch size $m$05, ZOO-ADMM reaches the best reported RMSE among the compared approaches after roughly $m$06 iterations. In sensor selection and sparse Cox regression, the method preserves the proximal and feasibility advantages of ADMM while remaining gradient-free (Liu et al., 2017).

Distributed stochastic zeroth-order optimization uses a different notion of simple averaging. The paper does not name a baseline “ZO-Avg,” but it identifies a distributed primal algorithm that performs consensus averaging over neighbors followed by a local zeroth-order descent step: $m$07 Each agent draws $m$08, evaluates the function at $m$09 and $m$10, and thus uses exactly two function queries per iteration (Yi et al., 2021).

Under connected-graph, smoothness, and state-dependent variance assumptions, this consensus-based scheme achieves

$m$11

with a known-horizon schedule, and under the Polyak–Łojasiewicz condition it achieves

$m$12

The leading terms are linear-speedup results in the number of agents $m$13 (Yi et al., 2021). Here, simple averaging is fundamentally consensus averaging of iterates, not multi-query averaging of directions.

6. Geometric and dynamical-systems perspectives

On Riemannian manifolds, the paper “Zeroth-order Riemannian Averaging Stochastic Approximation Algorithms” does not explicitly define a baseline called “Simple Averaging (ZO-Avg).” Instead, it proposes Zo-RASA, whose core recursion is an exponential moving average of transported zeroth-order gradient estimators,

$m$14

implemented either with parallel transport or with vector transport plus retractions (Li et al., 2023). The paper states that a natural uniform-average variant consistent with the framework is

$m$15

which aligns all past estimates in the current tangent space.

The distinction between EMA and uniform averaging is consequential. The paper proves its optimal $m$16 oracle-complexity results for schedules $m$17 or $m$18, with one-sample iterations after a warm-start in one regime, but it does not state an explicit theorem for the uniform-weight choice $m$19 (Li et al., 2023). Thus, in the Riemannian setting, “simple averaging” is best viewed as an admissible construction within the transport-and-averaging template rather than as the theorem-bearing algorithm of the paper.

A different generalization appears in the input-to-state stability analysis of zeroth-order methods. That work does not explicitly name an algorithm “ZO-Avg,” but it uses averaging in two senses. First, it averages $m$20 independent Gaussian oracle calls within each iteration, reducing the variance bound from $m$21 to $m$22. Second, it defines the averaged iterate as the conditional expectation of the random ZO trajectory and shows that the resulting dynamics satisfy the first-order update plus a bounded perturbation (Farzin et al., 28 Apr 2026). For gradient descent,

$m$23

with $m$24 bounded by a smoothing-bias term proportional to $m$25 and a variance term proportional to $m$26, up to the paper’s step-size-dependent constants. Analogous perturbed-recursion results are stated for heavy-ball and Nesterov acceleration (Farzin et al., 28 Apr 2026).

Under strong convexity or under $m$27-regularized non-strongly-convex dynamics, the paper applies input-to-state stability theorems to conclude that the averaged ZO recursion has the same decay rate as its first-order counterpart and converges to a neighborhood of the first-order fixed point whose radius can be made arbitrarily small by tuning $m$28, $m$29, and $m$30 (Farzin et al., 28 Apr 2026). In that sense, averaging is no longer merely a finite-sample variance reducer; it becomes the analytical bridge that expresses zeroth-order methods as perturbed first-order systems.

The broader conceptual background is classical first-order stochastic approximation. In the survey of stochastic-gradient averaging techniques, simple averaging refers to Polyak–Ruppert averaging of iterates,

$m$31

rather than averaging multiple gradient surrogates inside a single step. Under unbiased stochastic gradients, suitable smoothness and convexity assumptions, and classical stepsize conditions $m$32, $m$33, the survey states asymptotic normality with optimal covariance and reports non-asymptotic $m$34 rates for smooth convex problems (Lakshmanan, 10 Mar 2026). Tail averaging, window averaging, exponential moving averages, and stochastic weight averaging are then presented as refinements of the same bias-variance logic. The zeroth-order literature inherits that logic but redistributes it across queries, time, agents, or tangent spaces depending on the information model.

A recurrent misconception is therefore that “more averaging” is uniformly better. The literature shows a sharper picture. When ZO-Avg means averaging more random directions at a fixed point under a fixed query budget, the optimal choice can be the single-query case $m$35 (Lin et al., 19 Sep 2025). When averaging is amortized across time, embedded in control variates, aligned by transport, or used to form an averaged iterate with bounded perturbations, it can recover strong stability or complexity properties that naive multi-query averaging does not attain (Di et al., 2024, Li et al., 2023, Farzin et al., 28 Apr 2026).