--- title: Zeroth-Order Optimization on High-Res Images (ZO-high) url: https://www.emergentmind.com/topics/zeroth-order-optimization-on-high-resolution-images-zo-high type: topic --- # Zeroth-Order Optimization on High-Res Images (ZO-high) Searching arXiv for the cited papers and topic-specific context. Zeroth-order optimization on high-resolution images, often abbreviated ZO-high, denotes derivative-free optimization in image spaces whose dimensionality is determined by large spatial resolutions, typically with variables $x \in [0,1]^{H\times W\times C}$ and ambient dimension $d=HWC$. The setting is motivated by black-box objectives, non-differentiable pipelines, and regimes in which backpropagation is inaccessible or unaffordable. Its central difficulty is not merely the absence of gradients, but the interaction between image resolution, network depth, estimator variance, and query budget: coordinate-wise finite differences scale with spatial extent, random-direction estimators incur variance that grows with $d$, and conventional layer-by-layer zeroth-order training of deep networks exhibits quadratic dependence on depth [2602.10607][2502.01014][2510.19165]. ## 1. Problem class and scaling bottlenecks The canonical zeroth-order problem is $$ \min_{x\in \mathcal{X}\subset \mathbb{R}^d} f(x), $$ where only function values are observable [2102.10707]. In the image setting, $d=HWC$ can be very large even before any transform-domain lifting, and the effective dimension may become larger still in wavelet or continuous-wavelet representations; one reported audio-wavelet configuration reaches $d=1{,}776{,}000$, while an ImageNet discrete-wavelet configuration is reported with $d\approx 676{,}353$ [2102.10707]. Two scaling bottlenecks dominate ZO-high. First, high-resolution CNNs produce large feature maps $A\in \mathbb{R}^{C\times H\times W}$, so estimating gradients for all spatial positions in standard zeroth-order schemes would require $O(HW)$ forward queries per layer [2602.10607]. Second, in random-direction estimators the variance grows linearly in the ambient dimension. For the one-point multi-direction estimator used in R-AdaZO, the coordinate-wise variance bound is $$ \Sigma^2 \triangleq \frac{8(\sigma^2 + C^2)d}{K\mu^2}, $$ which decreases with the number of directions $K$ and increases as $\mu$ becomes small [2502.01014]. This creates the standard bias-variance-query trade-off: larger $K$ reduces estimator variance but increases query cost, while smaller $\mu$ reduces smoothing bias but amplifies numerical sensitivity and stochastic noise. Depth introduces an additional, distinct obstacle. In non-hierarchical activation-space ZO for deep networks of depth $L$ and average width $M$, repeated long forward paths induce $O(ML^2)$ query complexity [2602.10607]. In high-resolution models, the depth and resolution effects compound: the optimization variable is large, the intermediate activations are spatially extensive, and the cost of each function evaluation itself scales with the input size. ## 2. Gradient estimators and objective formulations A standard reference point in ZO-high is the two-point estimator $$ g(x) \;=\; \frac{f(x+\mu u) - f(x-\mu u)}{2\mu}\, u,\quad u \sim \mathcal{N}(0,I)\ \text{or Rademacher}, $$ which estimates the gradient of a smoothed objective and is widely used in black-box optimization [2602.10607]. R-AdaZO instead develops its theory around a one-point, multi-direction estimator with directions sampled uniformly from the sphere: $$ \hat{\nabla} f(\boldsymbol{\theta}, \xi) \triangleq \frac{d}{K}\sum_{k=1}^K \frac{f(\boldsymbol{\theta} + \mu \boldsymbol{u}_k;\, \xi) - f(\boldsymbol{\theta};\, \xi)}{\mu}\, \boldsymbol{u}_k, \quad \boldsymbol{u}_k \sim \mathrm{Unif}(S^{d-1}), $$ for which $$ \mathbb{E}\left[\hat{\nabla} f(\boldsymbol{\theta}, \xi)\right] = \nabla F_{\mu}(\boldsymbol{\theta}), \quad \mathbb{E}\left[\|\nabla F(\boldsymbol{\theta}) - \nabla F_{\mu}(\boldsymbol{\theta})\|\right] \le \mu L \sqrt{d}. $$ The latter inequality makes explicit that smoothing bias scales with both $\mu$ and $\sqrt{d}$ [2502.01014]. Coordinate and block estimators form the main alternative to random-direction methods. In the constrained block framework of ZOB-GDA and ZOB-SGDA, for a block $B$ of size $b$ and coordinate unit vector $e_i$, one-sided finite differences are defined by $$ \widehat{\partial_i s}(x) = \frac{s(x + \mu e_i) - s(x)}{\mu}, $$ with the block estimator $$ \widehat{\nabla}_B s(x) := P_B \sum_{i \in B} \widehat{\partial_i s}(x) e_i. $$ Under this model, a single gradient step requires $b$ shifted queries plus, at most, a single base evaluation, so the per-step query cost is $O(b)$ [2510.19165]. ZO-high objectives are often written either directly in image space or in a transform domain. A standard adversarial formulation is $$ \min_{\|\boldsymbol{\delta}\|_p \le \varepsilon} \; f(\boldsymbol{x}+\boldsymbol{\delta};\, \xi), $$ with projection onto an $\ell_p$ ball and clipping to valid pixel bounds [2502.01014]. In the wavelet-domain attack studied with ZO-BCD, the objective is written as $$ \min_{\delta}\; f(\mathrm{IWT}(\mathrm{WT}(x)+\delta)) + \lambda\|\delta\|_0, $$ optionally replacing the $\ell_0$ term by an $\ell_2$ or $\ell_\infty$ constraint [2102.10707]. These formulations are structurally different but share the same ZO constraint: the optimizer receives values, not gradients. ## 3. Hierarchical zeroth-order training for deep image models Hierarchical Zeroth-Order optimization (HZO) addresses the depth bottleneck by decomposing the composite network along depth rather than perturbing the full parameter vector or propagating finite differences layer by layer. Writing $$ F(x) = f_L \circ f_{L-1} \circ \cdots \circ f_1(x), $$ and, for a subnetwork, $$ \mathcal{N}_{i:j}(a_{i-1}) = f_j\big(f_{j-1}(\dots f_i(a_{i-1})\dots)\big), $$ HZO bisects the interval $[i,j]$ at $k=\lfloor(i+j)/2\rfloor$, computes the intermediate activation $a_k=\mathcal{N}_{i:k}(a_{i-1})$, estimates the Jacobian of the right subnetwork by symmetric coordinate perturbations, $$ J_{mn} \approx \frac{\big[\mathcal{N}_{k+1:j}(a_k + \epsilon e_n)\big]_m - \big[\mathcal{N}_{k+1:j}(a_k - \epsilon e_n)\big]_m}{2\epsilon}, $$ and propagates the target signal via $$ T_k = J^\top T_j. $$ At the leaves, it applies the local delta rule $$ \Delta W_i = -\eta\, T_i\, a_{i-1}^\top. $$ The recursion reconstructs the chain rule over depth using Jacobian-target products rather than backpropagated gradients [2602.10607]. The main complexity result is a reduction from $$ O(ML^2)\quad\text{to}\quad O(ML \log L), $$ where $L$ is depth and $M$ is average width. The argument is based on a balanced binary recursion over depth: at recursion level $k$, the right-side path length shrinks by $2^{k+1}$ while the number of subproblems grows as $2^k$, leaving per-level work constant at $ML$ and total work $ML\log_2 L$ [2602.10607]. This is a structural reduction in query complexity, not a mere constant-factor engineering improvement. The stability analysis is tied to the layer-wise Lipschitz factor $L_{\mathrm{lip}}$. Under $L_{\mathrm{lip}}$-Lipschitz continuity and $\beta$-smoothness, the recursive HZO error satisfies $$ \|\mathbf{e}_{\text{HZO}}\| = O\!\left(\beta\, \delta^2\, \frac{L_{\text{lip}^L} - 1}{L_{\text{lip} - 1}}\right). $$ Near the unitary limit, $L_{\mathrm{lip}}\approx 1$, $$ \|\mathbf{e}_{\text{HZO}}\| = O(L\, \beta\, \delta^2), $$ so cumulative error scales linearly rather than exponentially in depth [2602.10607]. The paper associates this regime with residual connections or orthogonal or normalized initialization, and reports that HZO maintains high gradient cosine similarity, approximately $0.95$ to $1.0$, in deep ResNets. For CNNs on large images, HZO adds Spatial Parallel Perturbation (SPP). If $R$ is the effective receptive field in the next layer and $\mathcal{S}$ is chosen so that $$ \forall (h,w), (h',w') \in \mathcal{S}, \quad \|(h,w)-(h',w')\|_\infty \ge R, $$ then perturbations at those positions are independent with respect to the next layer. The reported consequence is a reduction of per-layer spatial query complexity from $O(HW)$ to $O(R^2)$, making ZO-high effectively resolution-independent in the number of queries, although the cost of each forward pass still scales with $H$ and $W$ [2602.10607]. The empirical profile reported for HZO is training on CIFAR-10 to $74.2\%$ test accuracy in $8.14$ hours on a single NVIDIA A100 GPU, compared with a DeepZero baseline at $64.1\%$ in $28.15$ hours on $4\times$V100, together with $65.0\%$ top-1 accuracy on ImageNet-10 and gradient cosine similarity greater than $0.95$ throughout training [2602.10607]. The implementation notes further state that float32 is the minimum recommended precision for deep networks, float64 improves stability at increased cost, float16 leads to instability due to subtractive cancellation in finite differences, and GELU yields higher cosine similarity than ReLU. ## 4. Huge-scale block-coordinate and compressed-sensing methods ZO-BCD targets the “huge-scale” regime, explicitly described as $d>10^7$, where even basic vector operations on the decision variable are infeasible [2102.10707]. Its starting point is that full-dimensional zeroth-order estimators are impractical: coordinatewise finite differences require $O(d)$ queries per gradient and $O(d)$ memory per iteration, while compressed-sensing full-gradient schemes still incur at least $\Omega(s d \log d)$ memory and FLOPs in the sparse regime. The algorithm partitions $x$ into $J$ blocks and updates one randomly selected block per iteration. For a chosen block $j$, it uses one-sided finite differences with $m$ sensing directions $z_i$, $$ y_i \leftarrow \frac{E_f\!\big(x+\delta\,U^{(j)}z_i\big)-E_f(x)}{\sqrt{m}\,\delta}, $$ stacks them as $y\approx Z g^{(j)}$, and recovers a sparse block gradient through $$ \hat g^{(j)} = \arg\min_{\|v\|_0\le s}\|Z v - y\|_2, $$ implemented via CoSaMP [2102.10707]. The paper distinguishes ZO-BCD-R, which uses i.i.d. Rademacher directions, from ZO-BCD-RC, which uses rows of a circulant matrix generated by a single Rademacher seed vector. Under the stated assumptions, the total query complexity to attain an $\varepsilon$-optimal solution is $\tilde{\mathcal{O}}(s/\varepsilon)$ for both variants. The block strategy is paired with an “almost equisparse” effect: if blocks are assigned by a random permutation and are equal-sized, then with high probability each block has sparsity approximately $$ s_{\mathrm{block}} \approx 1.1\,\frac{s_{\mathrm{exact}}}{J}, $$ which reduces per-iteration query count by approximately $J$ [2102.10707]. The per-iteration computational profile is then reported as $\tilde{\mathcal{O}}(s\,d/J^2)$ FLOPs and memory for ZO-BCD-R, versus $\tilde{\mathcal{O}}(d/J)$ FLOPs and $\mathcal{O}(d/J)$ memory for ZO-BCD-RC. The circulant construction is central for ZO-high because it compresses both storage and matrix-vector multiply costs. Rather than storing a dense sensing matrix of size $m\times(d/J)$, ZO-BCD-RC stores only the seed vector and the row index set, reducing storage from $O(d s/J^2)$ to $O(d/J)$, and it exploits FFT-based circular convolutions to accelerate the CoSaMP inner loop [2102.10707]. The image-specific results are reported in transform space. On ImageNet with Inception-v3 in a three-level Daubechies “db45” DWT domain, ZO-BCD-R achieves attack success rate $92\%$, mean $\ell_2$ distortion $14.1$, and mean queries to first success $2131$; ZO-BCD-RC achieves $92\%$, $14.2$, and $2090$; and a large-coefficients-only variant with approximately $67\%$ dimension reduction achieves $96\%$, $13.7$, and $1662$ [2102.10707]. The same paper reports untargeted baselines ZO-SCD at $78\%$ ASR with $\ell_2$ distortion $57.5$ and $2400$ queries, ZO-SGD at $78\%$, $37.9$, and $1590$, ZO-AdaMM at $81\%$, $28.2$, and $1720$, and ZORO at $90\%$, $21.1$, and $2950$. In targeted audio attacks using a Morse CWT representation with $111$ frequencies, ZO-BCD is reported to achieve attack success rate $97.93\%$, mean relative loudness $-6.32$ dB, and mean queries approximately $7073$ [2102.10707]. ## 5. Adaptive variance-aware optimization in image-scale dimensions R-AdaZO addresses a different failure mode of ZO-high: the underutilization of moment information in adaptive zeroth-order methods [2502.01014]. Its first update is the momentum recursion $$ \boldsymbol{m}_t \gets \beta_1 \boldsymbol{m}_{t-1} + (1-\beta_1)\, \hat{\boldsymbol{g}}_t, $$ with no bias correction, specifically to preserve the variance-reduction effect of the first moment. The paper’s variance bound is $$ \mathbb{E}\left[\left|m_{t,i} - \nabla_i F_{\mu}(\boldsymbol{\theta}_{t-1})\right|^2\right] \le \frac{1-\beta_1}{1+\beta_1}\, \Sigma^2 + \text{squared bias terms}, $$ so the first moment reduces variance by the factor $(1-\beta_1)/(1+\beta_1)$ relative to the raw ZO estimator [2502.01014]. The second update departs from ZO-AdaMM by replacing $\hat g_{t,i}^2$ with $m_{t,i}^2$: $$ v_{t,i} \gets \beta_2 v_{t-1,i} + (1-\beta_2)\, m_{t,i}^2. $$ The update is then $$ \boldsymbol{\theta}_t = \boldsymbol{\theta}_{t-1} - \eta \frac{\boldsymbol{m}_t}{\sqrt{\boldsymbol{v}_t + \zeta}}. $$ The paper’s interpretation is that the refined second moment better approximates the variance-free ideal than the conventional choice based on $\hat g_t^2$, because the variance contribution to $v_t$ is reduced from $(1+\beta_1)\Sigma^2$ to $(1-\beta_1)\Sigma^2$ [2502.01014]. This is the basis for its variance-aware convergence framework. The convergence statement is given for non-convex objectives under bounded value, coordinate-wise Lipschitz smoothness, and bounded evaluation variance. The final bound takes the form $$ \frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}\big[\|\nabla F(\boldsymbol{\theta}_t)\|\big] \le \frac{(1+\beta_1)\sqrt{d}\sqrt{\beta_1(1-\beta_2)}\, \epsilon^2}{} + \left(\sqrt[4]{\zeta} + \sqrt{V\,d}\right) + \mu L \sqrt{d}, $$ with explicit dependence on $d$, the smoothed-gradient bias, and the moment statistics [2502.01014]. The paper states that this is the first variance-aware convergence framework for adaptive ZO methods and that R-AdaZO converges faster than ZO-AdaMM, with a speedup proportional to $$ \sqrt[4]{\frac{1+\beta_1}{1-\beta_1}} $$ in the variance-dominant regime. For high-resolution images, the method retains the query profile of standard direction-based ZO: $K+1$ function evaluations per iteration for the one-point estimator, or $2q$ for a two-point symmetric estimator [2502.01014]. The reported practical guidance is to use $\beta_1=0.9$, $\beta_2=0.99$, small $\zeta$, and task-dependent $\mu$; for high-dimensional images, the paper recommends $\beta_1$ in $[0.9,0.95]$, $\mu$ in $[5\mathrm{e}{-3}, 1\mathrm{e}{-2}]$, small $K$ or $q$ such as $1$ to $4$, and structured directions such as patch-wise perturbations or block-coordinate sampling. On MNIST adversarial attacks, the reported mean iterations to success are $15.6$k $\pm 3.2$k for ZO-AdaMM, $15.5$k $\pm 4.1$k for AdaZO, and $2.9$k $\pm 0.8$k for R-AdaZO, corresponding to an approximately $5.4\times$ speedup [2502.01014]. Although MNIST is not high-resolution, the paper explicitly frames the method as suitable for large-$d$ image settings because the variance term scales with $d$. ## 6. Constrained block ZO and the structure of query efficiency Query-efficient constrained ZO on high-resolution images is represented by ZOB-GDA and ZOB-SGDA, which solve $$ \min_{x \in \mathbb{R}^d} \; \max_{\lambda \in \mathcal{Y}} \; L(x,\lambda), \quad L(x,\lambda)=h(x)+\lambda^\top c(x), $$ with compact dual set $$ \mathcal{Y} = \{ \lambda \in \mathbb{R}^m \mid 0 \le \lambda \le \bar{\lambda} \} $$ [2510.19165]. The stationarity measure is built from the primal gradient and projected dual ascent: $$ \mathfrak{g}(x,\lambda) = \begin{pmatrix} \nabla_x L(x,\lambda) \\ \frac{1}{\beta}\big(\lambda - \Pi_{\mathcal{Y}}\big[\lambda + \beta \nabla_\lambda L(x,\lambda)\big]\big) \end{pmatrix}, $$ and an $\varepsilon$-stationary point satisfies $$ M(x,\lambda) := \min\big\{ \|\mathfrak{g}(x,\lambda)\|, \|\nabla \Phi_{1/2L}(x)\| \big\} \le \epsilon. $$ The two algorithms differ in their smoothing mechanism. ZOB-GDA applies block coordinate descent-ascent directly to $L$, while ZOB-SGDA introduces $$ K(x,\lambda; z) := L(x,\lambda) + \frac{p}{2}\|x-z\|^2, $$ updates $x$ using the block estimator of $\nabla_x K$, updates $\lambda$ by projected ascent on $c(x)$, and evolves the auxiliary variable by $$ z_{t+1} = \gamma x_{t+1} + (1-\gamma) z_t. $$ The reported rates are $$ \min_{t