ConceptioArchivearXiv CS
arXiv CSopen access

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

Unknown · 2026 · arxiv_cs
arXiv CS · Papers · License: Open Access · 2026
Open Source ↗Direct PDF ↓
machine learning, deep learning, neural networks

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

Kukyoung Jang 1 Taehyun Cho 2 Junrui Zhang 2 Ping Xu * 3 Kyungjae Lee * 1

arXiv:2605.27316v1 [cs.LG] 26 May 2026

Abstract

widely used in vision, inverse problems, and adversarial robustness (Blake & Zisserman, 1987; Brox & Malik, 2011; Szegedy et al., 2013; Mobahi & Fisher III, 2015; Hazan et al., 2016). Classical Gaussian homotopy relies on an explicit continuation schedule implemented through a multiloop procedure, leading to substantial computational overhead and sensitivity to hyperparameter tuning. To mitigate this cost, single-loop Gaussian smoothing methods update the optimization variable and smoothing scale simultaneously (Iwakiri et al., 2022). While these methods reduce computational burden, their theoretical guarantees remain local, characterizing convergence only to stationary points of the smoothed objective rather than the global optimum.

Probabilistic smoothing is a standard tool for global optimization, but existing methods rely on Gaussian kernels and specific transforms, often resulting in strong hyperparameter sensitivity and limited robustness. We propose a general smoothing framework that combines flexible symmetric unimodal kernels with monotonic ratio-based transformations. Under mild conditions, we show that the smoothed objective preserves the global maximizer and that all stationary points concentrate near the true optimum for sufficiently large amplification, without requiring a decreasing smoothing schedule. We further provide explicit complexity bounds for stochastic gradient ascent and show that a leave-one-out baseline provably reduces variance. Experiments on high-dimensional benchmarks and black-box adversarial attacks demonstrate improved robustness and competitive performance.

A major advance within Gaussian smoothing was introduced by Xu (2025), which replaces explicit continuation by applying power or exponential transformations before smoothing. This amplification forces stationary points of the smoothed objective to concentrate near the global maximizer, yielding the first approximate global optimality guarantee in a singleloop Gaussian smoothing setting. However, strong amplification also sharpens curvature and inflates the variance of Monte Carlo gradient estimators, resulting in pronounced sensitivity to the amplification parameter and a clear tradeoff between localization strength and optimization stability.

1. Introduction Global optimization over a compact domain is a longstanding challenge in machine learning and engineering, especially when the objective function is highly nonconvex, multimodal, or accessible only through black-box evaluations (Naser et al., 2025). In such settings, gradient-based methods are sensitive to initialization and often converge to suboptimal local extrema.

In parallel, zeroth-order optimization has shown that Gaussian smoothing enables gradient-free optimization via Monte Carlo score-function estimators (Ghadimi & Lan, 2013; Nesterov & Spokoiny, 2017; Chen et al., 2019). While subsequent work generalized the perturbation distribution beyond Gaussian kernels to reduce estimator variance (Gao & Sener, 2022), these approaches provide guarantees only for local stationarity and do not address global localization.

A classical approach to global optimization is homotopybased optimization, also known as continuation or graduated non-convexity methods, which replace the original objective with a smoothed surrogate and gradually recover the original landscape. Among these, Gaussian homotopy, which smooths the objective via Gaussian convolution, has been

This paper revisits probabilistic smoothing from a unified perspective. We view Gaussian homotopy, single-loop smoothing, and amplification-based methods as instances of a broader design space defined by (i) the smoothing distribution and (ii) the transformation applied to objective values prior to smoothing. From this viewpoint, instability arises not from smoothing itself, but from restrictive kernel choices combined with overly aggressive transformations.

1

Department of Statistics, Seoul, Korea University 2 Department of Electrical and Computer Engineering, Seoul, Seoul National University 3 Shandong University at Weihai, Weihai, China. Correspondence to: Ping Xu <[email protected]>, Kyungjae Lee <kyungjae [email protected]>.

Motivated by this observation, we propose Probabilistic Preprint. May 27, 2026.

1

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization Table 1. Comparison of representative zeroth-order optimization methods. The table summarizes the evaluation setting (stochastic or deterministic), the type of optimality guarantee (local stationarity or approximate global localization), the required assumptions on the objective f , and the iteration complexity to achieve ∥∇f ∥2 ≤ ε. Method

Sto. / Det.

Optimality

Assumptions on f

Iteration Complexity

RSGF (Ghadimi & Lan, 2013) ZO-SGD (Nesterov & Spokoiny, 2017) ZO-AdaMM (Chen et al., 2019) ZO-SLGH (Iwakiri et al., 2022) ZO-SLGH (Iwakiri et al., 2022) BeS (Gao & Sener, 2022) EPGS (Xu, 2025)

Stochastic Stochastic Stochastic Deterministic Stochastic Stochastic Deterministic

Local Local Local Local Local Local Approx. global

f smooth, Lipschitz gradient f smooth or weakly smooth f smooth, Lipschitz gradient f smooth, Lipschitz gradient f smooth, Lipschitz gradient f smooth, bounded; biased gradients f bounded on compact domain

Ours Ours with Variance Reduction

Deterministic Deterministic

Approx. global Approx. global

f bounded on compact domain f bounded on compact domain

O(d/ε2 ) Oa (M/ε2 ) O(d2 /ε2 ) Ob (d/ε) Ob (d/ε2 ) Oa (M/ε2 ) O(d4 /ε2 )  O d2 /ε2  c O (1 − C/2) d2 /ε2

a

The dependence on the dimension d is not explicit. Here, M denotes an upper bound on the second moment of the stochastic gradient estimator, which may implicitly depend on d. b The original results in Iwakiri et al. (2022) are stated for ∥∇f ∥ ≤ ε; all bounds are converted to ∥∇f ∥2 ≤ ε for consistency. c The constant C depends on the objective function f and algorithmic hyperparameters, but is independent of ε. While this does not change the asymptotic order of the iteration complexity, the exact bound is strictly smaller since (1 − C/2) < 1.

2. Related Work

Smoothing with Ratio-Monotone Transforms (ProMoT), a single-loop framework that generalizes both components. We allow a broad class of symmetric unimodal smoothing kernels, including heavy-tailed distributions, and introduce ratio-monotone transformations that subsume power and exponential forms. Under mild conditions, we show that the global maximizer of the original objective is preserved and that all stationary points of the smoothed objective concentrate near the true optimum without requiring a decreasing smoothing schedule.

Gaussian homotopy methods. Smoothing-based optimization replaces the original objective with a smoothed surrogate, most commonly via Gaussian convolution. Classical Gaussian homotopy or graduated optimization methods rely on a continuation schedule that gradually decreases the smoothing scale to recover the original landscape (Blake & Zisserman, 1987; Mobahi & Fisher III, 2015; Hazan et al., 2016). While effective in practice, such multi-loop procedures incur additional computational overhead and are sensitive to the choice of the smoothing schedule.

Finally, since probabilistic smoothing relies on Monte Carlo gradient estimation, its variance is a fundamental bottleneck. We introduce ProMoT-loo, a variance-reduced variant based on a leave-one-out baseline. We show that this estimator is unbiased and yields a strictly smaller second-moment bound, and we explicitly quantify how this variance reduction improves the iteration-complexity. We summarize our contributions as follows:

To mitigate this issue, single-loop Gaussian homotopy (SLGH) methods update the optimization variable and the smoothing scale simultaneously (Iwakiri et al., 2022). Under smoothness and Lipschitz-gradient assumptions on f , deterministic and stochastic variants of SLGH guarantee convergence to stationary points with iteration complexity O(d/ε) and O(d/ε2 ), respectively. As reported in Table 1, these guarantees are purely local, characterizing convergence to stationary points of the smoothed objective rather than to the global maximizer of f .

• We identify explicit conditions on both smoothing distributions and transformations, covering a broad class of symmetric unimodal kernels, including heavy-tailed distributions beyond Gaussian, and introducing ratiomonotone transformations that strictly generalize the power and exponential transforms of Xu (2025).

A further line of work removes explicit continuation schedules by modifying the objective prior to Gaussian smoothing. In particular, exponential and power-based Gaussian smoothing (EPGS) was the first to provide a global optimality guarantee in a Gaussian smoothing framework (Xu, 2025), with iteration complexity O(d4 /ε2 ).

• Under the proposed distributional and transformation conditions, ProMoT achieves the same ε-approximate global localization guarantees as Xu (2025). • We first introduce a leave-one-out variance reduction scheme for gradient estimation and demonstrate that it strictly improves the iteration-complexity constant.

Non-Gaussian kernels and zeroth-order smoothing. Zeroth-order optimization enables gradient-free learning by estimating gradients through randomized perturbations. Early work established convergence to stationary points under Gaussian smoothing for nonconvex objectives (Ghadimi & Lan, 2013; Nesterov & Spokoiny, 2017), with iteration

• We show that ProMoT is significantly more robust to hyperparameter choices than existing methods, with strong performance on high-dimensional benchmarks. 2

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

complexity on the order of O(d/ε2 ). Subsequent studies generalized the perturbation distribution beyond Gaussians to reduce estimator variance and improve sample efficiency (Chen et al., 2019; Gao & Sener, 2022). These approaches retain local stationarity guarantees and typically require smoothness or weak smoothness assumptions on f .

the product kernel p(z) :=



(4)

Second, the power/exponential transform is replaced by a general transform class that preserves order and expands level ratios. For an amplification parameter θ > 0 and a strictly increasing transform g(θ, ·), the transformed probabilistic smoothing is defined as    Gθ,σ (µ) := EX∼pµ,σ g θ, f (X) ; S .

(5)

Choosing p as the standard Gaussian and g(θ, y) = y recovers Gaussian homotopy (Mobahi & Fisher III, 2015); choosing p as Gaussian and g(θ, y) = y θ or eθy recovers power or exponential variants (Xu, 2025); the formulation in (5) subsumes both by allowing general p and a general ratio-monotone transform g. Our goal is to generalize the kernel p and transform g(θ, ·) so that maximizing Gθ,σ (µ) is equivalent to maximizing f (x) and admits a unique global maximizer reachable by gradient ascent. Optimization of the Smoothed Objective We next describe how the generalized smoothed objective Gθ,σ can be optimized in practice. Under mild regularity conditions on the smoothing kernel, the gradient of Gθ,σ with respect to the location parameter µ admits a score-function representation,

(1)

and x∗ is a global maximizer (i.e., x∗ ∈ arg maxx∈S f (x)). Direct gradient ascent on f can converge to a local maximizer depending on the initialization, which motivates smoothing-based approaches.

   ∇µ Gθ,σ (µ) = EX∼pµ,σ g θ, f (X) ∇µ log pµ,σ (X) ; S , (6)

which enables gradient estimation without requiring derivatives of f . Throughout this work, we adopt this formulation to derive stochastic optimization algorithms for Gθ,σ .

Gaussian homotopy (Mobahi & Fisher III, 2015) optimizes a smoothed surrogate (2)

Given a mini-batch {X(k) }B k=1 ∼ pµ,σ , define S(x) = ∇µ log pµ,σ (x) and h(x) = g(θ, f (x)) 1{x∈S} . A Monte Carlo estimator of the gradient is then given by

and gradually decreases the smoothing scale σ > 0 so that maximizers of Gσ track the global maximizer x∗ . While effective, such continuation schemes require carefully tuned multi-loop schedules and are computationally sensitive (Iwakiri et al., 2022). To remove explicit continuation, Xu (2025) proposed amplifying the objective prior to smoothing, enabling single-loop optimization with fixed σ. Specifically, power and exponential transforms, θ EPGS GPGS (µ) = E[eθf (X) ; S], θ,σ (µ) = E[f (X) ; S], Gθ,σ

p (xi − µi )/σ

denote the same kernel shifted to center µ and rescaled by σ (i.e., the original kernel moved to µ with bandwidth σ).

Consider a compact set S ⊂ Rd and a continuous function f : S → R. The global optimization problem is

  Gσ (µ) := EX∼N (µ,σI) f (X); S ,

d Y i=1

3. Probabilistic Smoothing with Ratio-Monotone Transforms

x∈S

i=1 p(zi ), and let

pµ,σ (x) := σ −d

Relation to the present work. The proposed approach builds on amplification-based probabilistic smoothing and extends prior work in two key directions. First, we generalize both the smoothing kernel and the transformation while preserving approximate global localization on compact domains. Second, we derive explicit iteration-complexity bounds that quantify the impact of variance reduction in generalized smoothing, showing that it improves the complexity constant without changing the asymptotic rate, as summarized in Table 1.

max f (x),

Qd

bµ Gθ,σ (µ) = 1 ∇ B

B X

h(X(k) ) S(X(k) ).

(7)

k=1

The parameter µ is updated by stochastic gradient ascent, bµ Gθ,σ (µ ), µt+1 = µt + ηt ∇ t

(3)

(8)

where ηt > 0 is a step size. Notably, the global localization guarantees established in this work do not rely on decreasing the smoothing scale σ over time.

concentrate stationary points near x for large θ. However, stronger amplification increases curvature and gradient variance, leading to a stability–localization trade-off and increased sensitivity to parameter tuning, a phenomenon that will be explicitly illustrated through a motivating example in a later section.

Coordinate-wise smoothing. The formulation readily extends to anisotropic smoothing. Instead of a single isotropic scale σ, we may use a diagonal matrix Σ = Qd diag(σ12 , . . . , σd2 ) and define pµ,Σ (x) = i=1 σi−1 p (xi − µi )/σi . This allows different smoothing resolutions across coordinates and is useful when the objective exhibits heterogeneous curvature or sensitivity along different dimensions.

The present formulation generalizes both lines along two axes. First, the Gaussian kernel is replaced by a general product density that can include even long-tail distribution. Let p be a one–dimensional probability density on R, define 3

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

(a) σ = 2.0 (b) σ = 2.5 (c) σ = 3.0 (d) σ = 3.5 Figure 1. Mean squared error (MSE) between the true global maximizer x∗ of the original objective f (x) and the solution obtained by maximizing the corresponding smoothed objective, plotted as a function of the amplification parameter θ (log scale). For each value of θ, we run stochastic gradient ascent on the smoothed objective starting from a fixed initialization and record the final estimate x̂(θ); the reported MSE is ∥x̂(θ) − x∗ ∥2 , averaged over 10 runs.

(a) θ = 10, σ = 2.0 (b) θ = 10, σ = 2.5 (c) θ = 10, σ = 3.0 (d) θ = 10, σ = 3.5 Figure 2. Effect of probabilistic smoothing on a multimodal one-dimensional objective. The dashed gray curve denotes the original objective f (x). Each colored curve represents a smoothed objective obtained by taking expectations under a shifted smoothing distribution with a fixed amplification parameter θ = 10 and varying the smoothing scale σ. The value θ = 10 corresponds to the gray dotted vertical reference line in Figure 1, and this figure visualizes how changing σ at that fixed θ reshapes the optimization landscape.

Variance reduction via leave-one-out baselines. While amplification improves localization, it can increase the variance of the stochastic gradient estimator. To mitigate this effect, we introduce a leave-one-out baseline that centers each sample contribution using information from the remaining samples in the same mini-batch, while preserving unbiasedness. PB For the baseline, define U := j=1 h(X(j) ) ∥S(X(j) )∥2 PB and V := j=1 ∥S(X(j) )∥2 . For each index k, the leaveone-out baseline is given by

bk = (U −h(X(k) ) ∥S(X(k) )∥2 )

3.1. A Motivating Example: Parameter Sensitivity in Probabilistic Smoothing We illustrate the effect of transformation and kernel generalization using a one-dimensional nonconvex objective f (x) with a unique global maximizer x∗ . For each method and each value of the amplification parameter θ, we construct the corresponding probabilistically smoothed surrogate objective by applying the specified transformation and smoothing operator. Starting from a fixed initialization, we run stochastic gradient ascent on each smoothed objective for a fixed number of iterations and record the final estimate x̂(θ). The reported error in Figure 1 is the mean squared error ∥x̂(θ) − x∗ ∥2 , averaged over multiple independent runs to account for stochasticity.

(V −∥S(X(k) )∥2 +λ), (9)

The key advantage of generalizing both the transformation and the smoothing kernel is improved optimization stability. As shown in Figure 1, exponential-based smoothing (EPGS) achieves low error only within a narrowly tuned range of the amplification parameter θ, and its performance degrades rapidly under mild misspecification. In contrast, ProMoT and ProMoT-loo maintain consistently low error over several orders of magnitude of θ, indicating substantially reduced sensitivity to hyperparameter choices. In particular, ProMoT-loo exhibits pronounced robustness with respect to the choice of θ.

where λ > 0 is a small ridge parameter for numerical stability. The resulting variance-reduced estimator is bµloo Gθ,σ (µ) = 1 ∇ B

B X

 h(X(k) ) − bk S(X(k) ).

(10)

k=1

This estimator remains unbiased and admits a strictly smaller second-moment bound than the baseline-free estimator. As shown in the theoretical analysis, this variance reduction leads to an explicit improvement in the iteration complexity constant, without altering the asymptotic convergence rate.

To further understand the origin of this stability, Figure 2 visualizes the smoothed optimization landscapes produced 4

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization Table 2. Examples of smoothing distributions satisfying Assumptions 4.1 and 4.2. R ′2 R Kernel p(z) I = (pp) K = |p′′ |

by different methods. Here, we directly plot the smoothed surrogate objectives obtained from the same base function f (x) at a fixed amplification level θ = 10 while varying the smoothing scale σ, without performing any optimization. The dashed gray curve denotes the original objective, and each colored curve corresponds to a different choice of σ.

2

Gaussian √12π e−z /2 Logistic e−z /(1 + e−z )2 Student-tν

The figure reveals a qualitative difference in the resulting smoothed landscapes. Under exponential amplification, multiple local extrema persist after smoothing, and increasing θ amplifies irregularities and leads to unstable gradient behavior. By contrast, combining order-preserving, moderated transformations with more flexible smoothing kernels suppresses spurious extrema and yields an effectively unimodal surrogate. This enhanced landscape regularity explains why ProMoT and ProMoT-loo converge reliably across a wide range of hyperparameters, without delicate tuning.

ν = 1 (Cauchy) ν=3 ν = 10  Hyper. Secant 12 sech πz/2 β Gen. Gaussian βe−|z| /(2Γ(1/β))

We establish δ-approximate global optimality and iteration complexity guarantees under a set of structural assumptions on the smoothing kernel and the transformation. These assumptions are introduced solely for analysis and enable control of smoothness, variance, and global optimality properties of the smoothed objective.

inf

(11)

p(x + δ) − p(x − δ) ≥ cδ (c, r).

(12)

|x−c|≤r

(a+c)θ (b+c)θ

>0 log a+c b+c (a − b) eθ(a−b) > 0 α α (aα − bα ) eθ(a −b ) > 0 β (a − b) (a+c) eθ(a−b) > 0 (b+c)β

log(1 + eθy )

a eθa (1+eθb )−b eθb (1+eθa ) >0 (1+eθa )(1+eθb ) sinh(θa) (a coth(θa) − b coth(θb)) > 0 sinh(θb)

Assumptions on Transformation Our framework also extends to a broad class of transformations. Specifically, we impose two structural assumptions on the transformation. Assumption 4.3 (Monotonicity). For any fixed θ > 0, y 7→ g(θ, y) is strictly increasing with respect to y ∈ R. Assumption 4.4 (Ratio-Monotonicity). For any fixed a > b > 0, the map θ 7→ g(θ, a)/g(θ, b) is strictly increasing. Concrete examples of transformations satisfying Assumption 4.4 are listed in Table 3, illustrating that the condition accommodates classical power and exponential forms as well as more general constructions.

Assumption 4.2 (Symmetry and Unimodality). Assume p is symmetric and unimodal so that p(z) = p(−z) and p′ (z) ≤ 0 for z > 0 to ensure limz→∞ p(z) = 0. Theorem 4.1 (Uniform non-degeneracy on a compact window). Fix δ > 0. For every r > 0 and every window center c such that |c| > r + δ (i.e., [c − r − δ, c + r + δ] does not intersect 0), there exist constants Cδ (c, r) > 0 and cδ (c, r) > 0 such that  F (x + δ) − F (x − δ) ≥ Cδ (c, r),

0.82691 0.87870 0.92883 π/2 3.36400

summarized in Table 2, together with the corresponding constants I and K that govern variance and curvature bounds in the subsequent analysis.

Assumption 4.1 (Fisher Information and Regularity). As2 I := Rsume′ p ∈2 C . Then, there exist two R constants, ′′ (p (z)) /p(z) dz < ∞ and K := |p (z)| dz < ∞ R R where I is the Fisher information of p and K reflects a regularity condition controlling the total variation of p′ .



(numeric)

(y + c)θ eθy α eθy (y + c)β eθy sinh(θy)

Assumptions on Smoothing Probability To extend the Gaussian smoothing framework to a more general class of smoothing distributions, we identify two fundamental distributional conditions.

inf

0.96749 0.38496

Table 3. Examples of ratio-monotone transformations g(θ, y).   g(θ,a) d g(θ, y) for a > b dθ g(θ,b)

4. Complexity Analysis

|x−c|≤r

1 1/3 ν+1 ν+3 0.50000 0.66667 0.84615 π 2 /8 4.05587

4.1. Invariance of Global-Optimum Theorem 4.2. Suppose the smoothing kernel p(x) and transformation g(θ, y) satisfy Assumption 4.1–4.4. Then, for any M > 0 and δ > 0 such that cube(x∗ ; δ) := {x|∀i ∈ [d], |xi − x∗i | ≤ δ} and cube(x∗ ; δ) ⊂ S, there exists θδ,σ,M > 0, such that whenever θ > θδ,σ,M , for any ∥µ∥∞ < M and any i ∈ [d], we have ∂Gθ,σ (µ)/∂µi > 0 if µi < x∗i − δ, and ∂Gθ,σ (µ)/∂µi < 0 if µi > x∗i + δ. Here, µi and x∗i denote the ith dimension of µ and x∗ .

The proof can be found in Appendix B.1. This theorem provides the sign structure and tail behavior needed for localization arguments, ensuring that derivative signs are consistent and that mass vanishes at infinity. Representative smoothing kernels satisfying Assumptions 4.1 and 4.2 are

Remark. See Appendix B.2 for the proof. The inequalities in Theorem 4.2 imply that any stationary point of Gθ,σ within the region {∥µ∥∞ < M } must lie inside the cube(x∗ ; δ). In other words, the gradient cannot vanish outside this cube, 5

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

since each coordinate derivative is strictly positive to the left of x∗i − δ and strictly negative to the right of x∗i + δ. Thus, all stationary points (if any exist) are confined to a neighborhood of the true maximizer x∗ . While Xu (2025) analyzes a similar property under an ℓ2 -ball with Gaussian kernels, such analysis becomes intractable for general product kernels. We instead use an ℓ∞ -cube, which better suits the product kernel and coordinate-wise control.

+

(13)

The proof can be found in Appendix B.6. Theorem 4.6 establishes a generic inequality for the score-based stochastic ascent applied to Gθ,σ . In contrast to Gaussian homotopy analyses (Iwakiri et al., 2022; Xu, 2025), this argument does not rely on Gaussian distribution or gradients of the original objective f (x), but applies directly to a general kernel through the score-function formulation.

Theorem 4.3. Fix σ > 0 and θ > 0. Under Assumption 4.1–4.4, if there exists a unique global maximizer, then the objective Gθ,σ (µ) admits the global maximizer.

Corollary 4.7. Under Assumption 4.1–4.4, fix σ > 0, M > 0, and δ > 0. Then there exists θδ,σ,M > 0 such that for all θ ≥ θδ,σ,M , ProMoT converges into a local maximizer of Gθ,σ (µ) within {∥µ∥∞ < M } lies in cube(x∗ ; δ). Moreover, consider the step-size ηt = 1 σ 2 /d · (t + 1)−( 2 +γ) for γ ∈ (0, 1/2). Then, for any  ε > 0, ProMoT achieves E ∥∇µ Gθ,σ (µT )∥2 < ε, when2  1−2γ 3 ever T > Cγ d · (g∗,θ + I max(K, I) g∗,θ )/(σ 2 ε)

The proof can be found in Appendix B.3. Theorem 4.3 guarantees that the smoothed objective Gθ,σ actually attains its maximum. This existence result does not depend on Gaussian convolution structure and therefore applies to the broader kernel class considered in this work. 4.2. Iteration Complexity of ProMoT This subsection derives an iteration-complexity guarantee for ProMoT by reducing the analysis to two reusable ingredients: a global smoothness constant for Gθ,σ and a uniform second-moment bound for the score estimator. We first show that the Lipschitz constant of ∇µ Gθ,σ is controlled by K (Lemma 4.4). Next, we bound the second moment of the score estimator via I (Lemma 4.5). We apply a standard one-step smoothness inequality and telescope to control the cumulative squared gradient norm (Theorem 4.6), which yields explicit rates under polynomial steps (Corollary 4.7) and their anisotropic variants (Corollaries 4.8). Finally, we incorporate variance reduction via a leave-one-out baseline, proving unbiasedness and deriving an improved secondmoment bound that translates into a strictly better complexity (Lemmas 4.9–4.10, Theorem 4.11, Corollaries 4.12).

1

where Cγ indicates ( 21 − γ)/(2 2 −γ − 1). The proof can be found in Appendix B.7. Corollary 4.7 combines Theorem 4.2 with the iteration bound of Theorem 4.6. As a result, ProMoT is guaranteed not only to approach a stationary point of Gθ,σ , but also to approach one that lies within a prescribed δ-neighborhood of the true global maximizer x∗ . Moreover, as γ → 0, the iteration complexity bound approaches the optimal rate shown in Table 1. Corollary 4.8. Under Assumption 4.1–4.4, assume the single base kernel setting with product density pµ,Σ (x) =  Qd −1 −2 i=1 σi p (xi − µi )/σi . Define S2 (Σ) := d maxi σi . Then, for any fixed (θ, Σ): 2 Lθ (Σ) = g∗,θ max(K, I) S2 (Σ), Qθ (Σ) = g∗,θ I S2 (Σ). (14)

The following two lemmas characterize the global Lipschitz continuity of the gradient and the second moment of its Monte–Carlo gradient estimator.

1

Moreover, consider the step-size ηt = (t + 1)−( 2 +γ) · 1 S2 (Σ)−1 for Then, ProMoT  any γ ∈ 2 (0, 2 ). achieves E ∥∇µ Gθ,σ (µT )∥ < ε, whenever T > 2  1−2γ 3 Cγ · S2 (Σ) · (g∗,θ + I max(K, I) g∗,θ )/ε .

Lemma 4.4. Under Assumption 4.1–4.4, the gradient d ∇µ Gθ,σ (µ) is globally Lipschitz on  2R with Lipschitz constant Lθ = g∗,θ d max(K, I) σ , where g∗,θ := g(θ, f (x∗ )). Lemma 4.5. Under its second moment  Assumption24.1–4.4,  ˆ µ Gθ,σ (µ)∥ ≤ Qθ := d I g 2 /σ 2 . is bounded as E ∥∇ ∗,θ

The proof can be found in Appendix B.8. Corollary 4.8 shows that ProMoT allows different smoothing scales to be assigned to different coordinates. Each coordinate i is smoothed with its own scale σi , and both the smoothness constant and the second-moment P −2bound depend on the aggregate quantity S2 (Σ) = i σi . Coordinates with smaller σi contribute more strongly to the curvature and variance terms, making explicit how dimension-wise sensitivity affects the convergence rate. In contrast, standard Gaussian smoothing enforces a single isotropic scale across all coordinates.

The proofs can be found in Appendix B.4 and Appendix B.5. Under these lemmas, we can derive the following inequality. Theorem 4.6. Under Assumption 4.1–4.4, ∇µ Gθ,σ is Lθ Lipschitz and the score estimator has a uniform secondmoment bound with Qθ . Consequently, for any horizon −1 T ≥ 1 and nonnegative steps {ηt }Tt=0 , T −1 X

T −1 g∗,θ d max(K, I) d I 2 X 2 · 2 g∗,θ ηt . 2 σ 2σ t=0

    ηt E ∥∇µ Gθ,σ (µt )∥2 ≤ g∗,θ − E Gθ,σ (µ0 )

t=0

6

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

4.3. Iteration Complexity of ProMoT-loo

See Appendix B.11 for the proof. Theorem 4.11 selects an explicit ridge parameter that balances baseline estimation bias and variance. The key implication of Theorem 4.11 is that introducing a leave-one-out baseline strictly reduces the variance of the score estimator compared to the baseline-free estimator. The leading second-moment term is multiplied by (1 − Cθ,Σ /2), which is always no larger than one and is strictly smaller whenever the baseline is informative. Corollary 4.12. Under the same conditions of Theorem 4.11. Let the step size be ηt = S2 (Σ)−1 (t + 1 1)−( 2 +γ) with any γ ∈ (0, 21 ). Then, ProMoT achieves E[∥∇µ Gθ,Σ (µT )∥2 ] < ε, whenever

We now turn to the theoretical analysis of ProMoT-loo, which incorporates a leave-one-out baseline into the scorefunction estimator. Lemma 4.9. Under Assumption 4.1–4.4, consider the leaveone-out estimator defined in (10). Then, conditional on µt , h i bµloo Gθ,Σ (µ ) µ = ∇µ Gθ,Σ (µ ). E ∇ t t t

(15)

See Appendix B.9 for the proof. Lemma 4.9 shows that introducing a leave-one-out baseline does not change the expected update direction. This allows variance reduction to be incorporated without biasing the optimization. Lemma 4.10. Under Assumption 4.1–4.4, define 2 Rθ,µ,Σ :=

2 E[h(X)∥S(X)∥2 | µ] . E[h(X)2 ∥S(X)∥2 |µ] E[∥S(X)∥2 |µ]

 T >

(16)

The proof can be found in Appendix B.12. Corollary 4.12 shows that the leave-one-out baseline yields a provable reduction in the iteration complexity of ProMoT by strictly improving the second-moment bound of the score estimator. Moreover, as γ → 0, the iteration complexity bound approaches the rate shown in Table 1. To the best of our knowledge, this is the first analysis that makes the benefit of variance reduction explicit at the level of iteration complexity, rather than variance stabilization alone. This theoretical improvement is consistent with the empirical results in Figure 1, where ProMoT-loo achieves the best convergence behavior among all compared methods.

Then, conditionally on µt ,   2 2 2 loo b E ∇µ Gθ,Σ (µt ) µt ≤ (1 − Rθ,µ ) g∗,θ I S2 (Σ) t ,Σ + Cloo (B, λ) I S2 (Σ),

(17)

where 2  ⋆ 2 Std(U ) 2µU Std(V ) b λ √ √ + + Cloo (B, λ) := µV + λ λ s λ2 s    1  2 Var(U ) 4µ2 Var(V ) 1 +O + + U 4 + O 4 2 , (18) 2 2 λ s λ s λ s λ s

with s := B − 1 and µU = E[h(X)∥S(X)∥2 |µ], µV = E[∥S(X)∥2 |µ], b⋆ = µU /µV . Here, O(·) ignores constant depending only on U and V (but not on B or λ).

5. Experiments In this section, we evaluate ProMoT and ProMoT-loo on (i) canonical high-dimensional non-convex optimization benchmarks and (ii) real-world black-box targeted adversarial attack tasks. We compare against the following baselines: EPGS (Xu, 2025), RSGF (Ghadimi & Lan, 2013), ZOSGD (Nesterov & Spokoiny, 2017), ZO-AdaMM (Chen et al., 2019), ZO-SLGHd/r (Iwakiri et al., 2022), and CMAES (Hansen & Ostermeier, 2001; Hansen et al., 2019).

See Appendix B.10 for the proof. Lemma 4.10 characterizes how variance reduction through a leave-one-out baseline affects the second moment of the score estimator. The leading term shows that the variance (or the second moment) is 2 2 reduced proportionally to (1 − Rθ,µ,Σ ), where Rθ,µ,Σ measures the alignment between the score magnitude and the transformed objective. The additional term Cloo (B, λ) captures the bias induced by using a regularized finite-sample baseline defined in (9), instead of the variance-optimal one b⋆ , and is explicitly controlled by the ridge parameter λ. Theorem 4.11. With the setting of Lemma 4.9 and 4.10, choose λ = (B − 1)−1/8 . Then Cloo (B, λ) = O(B −1/4 ) and  2 bµloo Gθ,Σ (µ ) µ ∇ t t   2 2 ≤ (1 − Rθ,µ ) g∗,θ + O(B −1/4 ) IS2 (Σ). t ,Σ

Evaluation metrics. For canonical optimization benchmarks, we report: (i) Mean Squared Error (MSE), defined as the minimum squared distance between the optimization trajectory and the true global optimum; (ii) hitting time, defined as the iteration at which this minimum MSE is first achieved; and (iii) best value, defined as the objective value corresponding to the minimum MSE. All metrics are averaged over independent runs, with standard deviations reported in parentheses.

 E

Define Cθ,Σ

h i 2 := mint∈[T ] Eµt Rθ,µ and B t ,Σ

(19)

:=

For adversarial attack experiments, we report: (i) success rate (SR), defined as the fraction of inputs for which a successful attack is achieved; (ii) mean R2 (x, x + µ), measuring similarity between the original input and the perturbed

2 Ω(16/(Cθ,Σ g∗,θ )4 ), then,

  2 2 bµloo Gθ,Σ (µ ) E ∇ ≤ (1 − Cθ,Σ /2) g∗,θ IS2 (Σ). t

 2 o 1−2γ Cγ S2 (Σ) n C 3 g∗,θ + I max(K, I) g∗,θ 1 − θ,Σ . 2 ε (21)

(20)

7

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization Table 7. CIFAR-10 SR mean R2 1 0.97(0.02) 1 0.98(0.01) 1 0.97(0.02) 0.34 -0.50(1.50) 1 0.98(0.02) 1 0.96(0.03) 0.99 0.98(0.01) 1 0.95(0.04) 1 0.94(0.04) Table 8. VitalDB Method SR mean R2 ProMoT 1 0.97(0.01) ProMoT-loo 1 0.98(0.01) EPGS 1 0.97(0.01) RSGF 1 0.91(0.04) ZO-SGD 1 0.83(0.09) ZO-AdaMM 1 0.96(0.02) ZO-SLGHd 1 0.87(0.04) ZO-SLGHr 1 0.95(0.02) CMA-ES 1 0.99(0.01)

input; and (iii) mean L∞ norm of the perturbation, capturing the maximum per-coordinate distortion. Higher R2 and lower L∞ indicate less perceptible adversarial attacks. Method ProMoT ProMoT-loo EPGS RSGF ZO-SGD ZO-AdaMM ZO-SLGHd ZO-SLGHr CMA-ES

Table 4. Ackley MSE Hitting Time 0.49(0.04) 393.60(7.70) 0.04(0.00) 392.60(9.06) 5.98(0.93) 399.50(1.07) 24.90(0.13) 22.65(38.56) 0.28(0.05) 397.75(3.30) 12.25(0.70) 395.40(5.65) 0.58(0.15) 398.70(1.62) 0.57(0.11) 399.25(1.22) 0.09(0.01) 397.85(1.35)

Method ProMoT ProMoT-loo EPGS RSGF ZO-SGD ZO-AdaMM ZO-SLGHd ZO-SLGHr CMA-ES

Table 5. Rosenbrock MSE Hitting Time 0.13(0.01) 398.95(1.50) 0.02(0.01) 396.40(4.66) 0.37(0.01) 399.85(0.36) 3.99(0.01) 0.85(1.24) 0.06(0.01) 400.00(0.00) 0.05(0.01) 400.00(0.00) 0.44(0.02) 400.00(0.00) 0.40(0.02) 400.00(0.00) 0.65(0.03) 96.90(8.04)

Method ProMoT ProMoT-loo EPGS RSGF ZO-SGD ZO-AdaMM ZO-SLGHd ZO-SLGHr CMA-ES

Table 6. Griewank MSE Hitting Time 1.95(0.12) 399.05(1.02) 0.23(0.01) 394.60(7.53) 3.32(0.17) 399.30(1.58) 24.86(0.17) 24.00(29.72) 3.62(0.06) 400.00(0.00) 0.24(0.05) 397.10(5.10) 2.77(0.16) 228.00(8.88) 2.81(0.17) 227.35(9.98) 1.16(0.11) 395.70(2.59)

Method ProMoT ProMoT-loo EPGS RSGF ZO-SGD ZO-AdaMM ZO-SLGHd ZO-SLGHr CMA-ES

Best Value -4.33(0.11) -1.71(0.05) -9.41(0.49) -13.42(0.69) -3.37(0.16) -10.76(0.20) -4.32(0.31) -4.33(0.23) -2.61(0.10)

Best Value -38.63(2.30) -3.11(0.21) -113.43(4.32) -3622.45(29.14) -3.81(0.44) -5.95(0.45) -120.81(4.36) -111.46(3.61) -449.42(56.74)

mean L∞ 0.28(0.07) 0.20(0.06) 0.28(0.07) 1.49(0.54) 0.25(0.07) 0.31(0.08) 0.22(0.06) 0.31(0.09) 0.32(0.04) mean L∞ 0.58(0.17) 0.57(0.17) 0.65(0.18) 0.92(0.23) 1.24(0.33) 0.82(0.27) 1.24(0.25) 0.80(0.19) 0.60(0.27)

transformations, together with the leave-one-out estimator, yields stable and consistently superior performance. 5.2. Black-Box Targeted Adversarial Attacks We evaluate ProMoT and ProMoT-loo on black-box targeted adversarial attacks against models trained on CIFAR-10 and VitalDB, with results summarized in Tables 7 and 8. For CIFAR-10, the input dimension is d = 3,072, corresponding to 32 × 32 × 3 images. For VitalDB, the final input dimension after preprocessing is d = 42. Across both domains, ProMoT-loo consistently achieves the smallest L∞ perturbations while maintaining high attack success rates and competitive mean R2 , indicating less perceptible adversarial examples. On CIFAR-10, ProMoT-loo outperforms all baselines in terms of L∞ without sacrificing reconstruction quality, while on VitalDB—a non-smooth tree-based classifier with severe class imbalance—it again yields the smallest perturbations and remains competitive in R2 . Notably, although ProMoT-loo attains a slightly lower mean R2 than CMA-ES (by approximately 1%), it achieves about a 5% reduction in mean L∞ , demonstrating improved imperceptibility with comparable reconstruction quality. Overall, these results highlight the robustness of the proposed framework across heterogeneous black-box models and data modalities.

Best Value -1.03(0.35) 74.71(4.41) -1.42(0.02) -4.11(0.02) -1.45(0.01) 105.38(1.79) -0.75(0.31) -0.89(0.30) 38.55(4.15)

5.1. Canonical Non-Convex Benchmarks All canonical benchmark experiments are conducted in a high-dimensional setting with dimension d = 500. Across the Ackley, Rosenbrock, and Griewank benchmarks (Tables 4–6), ProMoT-loo consistently achieves the lowest MSE, indicating reliable convergence toward the global optimum across diverse non-convex landscapes. On Ackley and Rosenbrock, it attains both the lowest MSE and the best fitted values, demonstrating robustness to strong multimodality and curvature anisotropy, while on Griewank it maintains the lowest MSE with competitive objective values. One contributing factor to this advantage in high dimensions is the choice of the smoothing kernel. As the dimension increases, Gaussian smoothing concentrates most samples within a narrow neighborhood (e.g., within 3σ), which limits effective exploration, whereas the heavier-tailed kernels used in ProMoT maintain a higher probability of sampling distant regions, enabling more effective global search and improved performance. Overall, these results highlight that combining heavier-tailed smoothing kernels with ratio-monotone

6. Conclusion We introduced ProMoT, a single-loop probabilistic smoothing framework that achieves approximate global optimization using general smoothing distributions and ratiomonotone transformations. We established global localization and iteration-complexity guarantees with provable variance reduction, and validated robustness on highdimensional benchmarks and black-box adversarial attacks. 8

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

Impact Statement

Lee, H.-C., Park, Y., Yoon, S. B., Yang, S. M., Park, D., and Jung, C.-W. Vitaldb, a high-fidelity multi-parameter vital signs database in surgical patients. Scientific Data, 9(1): 279, 2022.

This paper presents work whose goal is to advance the field of Machine Learning. There are many potential societal consequences of our work, none which we feel must be specifically highlighted here.

Mobahi, H. and Fisher III, J. W. A theoretical analysis of optimization by gaussian continuation. In Bonet, B. and Koenig, S. (eds.), Proceedings of the Twenty-Ninth AAAI Conference on Artificial Intelligence, January 2530, 2015, Austin, Texas, USA, pp. 1205–1211. AAAI Press, 2015.

References Blake, A. and Zisserman, A. Visual reconstruction. MIT press, 1987. Brox, T. and Malik, J. Large displacement optical flow: Descriptor matching in variational motion estimation. IEEE Trans. Pattern Anal. Mach. Intell., 33(3):500–513, 2011.

Naser, M. Z., Al-Bashiti, M. K., Tapeh, A. T. G., Naser, A., Kodur, V., Hawileh, R., Abdalla, J., Khodadadi, N., Gandomi, A. H., and Eslamlou, A. D. A review of benchmark and test functions for global optimization algorithms and metaheuristics. WIREs Comput. Stat., 17(2), May 2025.

Carlini, N. and Wagner, D. Towards evaluating the robustness of neural networks. In 2017 ieee symposium on security and privacy (sp), pp. 39–57. Ieee, 2017.

Nesterov, Y. E. and Spokoiny, V. G. Random gradient-free minimization of convex functions. Found. Comput. Math., 17(2):527–566, 2017.

Chen, X., Liu, S., Xu, K., Li, X., Lin, X., Hong, M., and Cox, D. Zo-adamm: Zeroth-order adaptive momentum method for black-box optimization. Advances in neural information processing systems, 32, 2019.

Papernot, N., McDaniel, P., Wu, X., Jha, S., and Swami, A. Distillation as a defense to adversarial perturbations against deep neural networks. In 2016 IEEE symposium on security and privacy (SP), pp. 582–597. IEEE, 2016.

Gao, K. and Sener, O. Generalizing gaussian smoothing for random search. In International Conference on Machine Learning, 17-23 July 2022, Baltimore, Maryland, USA, volume 162 of Proceedings of Machine Learning Research, pp. 7077–7101. PMLR, 2022.

Szegedy, C., Zaremba, W., Sutskever, I., Bruna, J., Erhan, D., Goodfellow, I. J., and Fergus, R. Intriguing properties of neural networks. In Bengio, Y. and LeCun, Y. (eds.), 2nd International Conference on Learning Representations, Banff, AB, Canada, April 14-16, 2014, Conference Track Proceedings, 2013.

Ghadimi, S. and Lan, G. Stochastic first-and zeroth-order methods for nonconvex stochastic programming. SIAM journal on optimization, 23(4):2341–2368, 2013.

Xu, C. Global optimization with a power-transformed objective and gaussian smoothing. In Forty-second International Conference on Machine Learning, Vancouver, BC, Canada, July 13-19, 2025. OpenReview.net, 2025.

Hansen, N. and Ostermeier, A. Completely derandomized self-adaptation in evolution strategies. Evolutionary computation, 9(2):159–195, 2001. Hansen, N., Akimoto, Y., and Baudis, P. CMA-ES/pycma on Github, February 2019. Hazan, E., Levy, K. Y., and Shalev-Shwartz, S. On graduated optimization for stochastic non-convex problems. In Balcan, M. and Weinberger, K. Q. (eds.), Proceedings of the 33nd International Conference on Machine Learning, ICML, New York City, NY, USA, June 19-24, 2016, volume 48, pp. 1833–1841. JMLR.org, 2016. Iwakiri, H., Wang, Y., Ito, S., and Takeda, A. Single loop gaussian homotopy method for non-convex optimization. In Annual Conference on Neural Information Processing Systems, New Orleans, LA, USA, November 28 - December 9, 2022, 2022. Krizhevsky, A., Hinton, G., et al. Learning multiple layers of features from tiny images. 2009. 9

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

A. Verification of Transformation Assumptions In this section, we verify that the transformation families listed in Table 3 satisfy Assumption 4.3 (measurability, nonnegativity, and monotonicity) and Assumption 4.4 (ratio monotonicity). When applicable, we also explicitly state a boundedness constant Gθ . Throughout, we fix θ > 0 and assume that y lies in the domain where each transformation is well-defined (e.g., y > −c when (y + c) appears). Power: g(θ, y) = (y + c)θ , c ≥ 0, y > −c Assumption 4.3. We have ∂y g(θ, y) = θ(y + c)θ−1 ≥ 0

for θ > 0,

so g is nondecreasing in y and nonnegative. Measurability is immediate. Assumption 4.4. For a > b > −c, g(θ, a)  a + c θ , = g(θ, b) b+c and since (a + c)/(b + c) > 1, the ratio is strictly increasing in θ. Exponential: g(θ, y) = eθy Assumption 4.3. ∂y g(θ, y) = θeθy ≥ 0, so g is nondecreasing, nonnegative, and measurable. Assumption 4.4. For a > b, g(θ, a) = eθ(a−b) , g(θ, b) which is strictly increasing in θ. Softplus (log-exp): g(θ, y) = log 1 + eθy



Assumption 4.3. We have ∂y g(θ, y) = θ σ(θy) ≥ 0,

σ(t) =

1 , 1 + e−t

so g is nondecreasing in y for θ > 0. Moreover, g(θ, y) ≥ 0 and is measurable, hence Assumption 4.3 is satisfied. Assumption 4.4. For a > b, consider R(θ) := log g(θ, a) − log g(θ, b) = log

 log(1 + eθa )  log(1 + eθb )

.

It suffices to show that R′ (θ) ≥ 0. Define ϕ(t) := t σ(t) − log(1 + et ). A direct calculation yields R′ (θ) =

ϕ(θa) − ϕ(θb)  . log(1 + eθa ) log(1 + eθb )

Since  ϕ′ (t) = t σ(t) 1 − σ(t) ≥ 0

for t ≥ 0,

the function ϕ is increasing on [0, ∞). Thus, for a > b > 0 and θ > 0, we have ϕ(θa) ≥ ϕ(θb) and hence R′ (θ) ≥ 0. Therefore, θ 7→ g(θ, a)/g(θ, b) is nondecreasing, establishing Assumption 4.4. Remark. Compared to pure exponential transforms, the softplus form grows exponentially only for large θy while remaining nearly linear for small θy, which moderates amplification and improves gradient conditioning in practice. 10

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

Sigmoid–power (bounded): g(θ, y) = σ(αy)θ , α > 0 Assumption 4.3. Since σ(αy) is increasing in y and takes values in (0, 1), g(θ, y) is nondecreasing, nonnegative, and measurable. Moreover, 0 < g(θ, y) ≤ 1, so the boundedness requirement in Assumption 4.3 holds with Gθ = 1. Assumption 4.4. For a > b, g(θ, a)  σ(αa) θ = , g(θ, b) σ(αb) which is strictly increasing in θ. α Fractional exponential: g(θ, y) = exp{θ y+ }, α > 0, y+ = max{y, 0} α Assumption 4.3. The map y 7→ y+ is nondecreasing and nonnegative, hence g is nondecreasing, nonnegative, and measurable.

Assumption 4.4. For a > b, g(θ, a) α = exp{θ(aα + − b+ )}, g(θ, b) which is strictly increasing in θ. Power–exponential hybrid: g(θ, y) = (y + c)β eθy , β ≥ 0, c ≥ 0, y > −c Assumption 4.3. ∂y g(θ, y) = β(y + c)β−1 eθy + θ(y + c)β eθy ≥ 0, so g is nondecreasing and nonnegative. Assumption 4.4. For a > b > −c, g(θ, a)  a + c β θ(a−b) = e , g(θ, b) b+c which is strictly increasing in θ. Hyperbolic variant: g(θ, y) = sinh(θ(y + c))+ Restricting to y ≥ −c, g is nondecreasing and nonnegative, satisfying Assumption 4.3. For a > b, the ratio sinh(θ(a + c))/ sinh(θ(b + c)) is strictly increasing in θ, establishing Assumption 4.4. Summary All transformation families above satisfy Assumptions 4.3 and 4.4 under the stated parameter and domain conditions. Exponential and fractional-exponential forms provide strong contrast at the cost of increased variance, while softplus and sigmoid-based transformations offer moderated amplification and improved stability through effective boundedness.

B. Proofs of Theoretical Results B.1. Proof of Theorem 4.1 The proof can be done by following two lemmas. Lemma B.1 (CDF gap on a compact window). Fix δ > 0. Let p be a density that is continuous and strictly positive on an interval [a, b]. Let F be its CDF. If [x − δ, x + δ] ⊂ [a, b], then Z x+δ F (x + δ) − F (x − δ) = p(t) dt ≥ 2δ · inf p(t). t∈[a,b]

x−δ

Proof. By definition, F (x + δ) − F (x − δ) = in the integral range. Hence

R x+δ x−δ

p(t) dt. Since [x − δ, x + δ] ⊂ [a, b], we have p(t) ≥ inf [a,b] p for all t

Z x+δ

Z x+δ p(t) dt ≥

x−δ

inf p dt = 2δ inf p. x−δ

11

[a,b]

[a,b]

(22)

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

Lemma B.2 (PDF difference gap away from the mode). Assume p ∈ C 1 (R) is symmetric and unimodal, and strictly decreasing on (0, ∞). Fix δ > 0 and assume |c| > r + δ. Then there exists cδ (c, r) > 0 such that inf

|x−c|≤r

p(x + δ) − p(x − δ) ≥ cδ (c, r).

(23)

Proof. By symmetry it suffices to consider c > 0. The condition c > r + δ implies x − δ > 0 and x + δ > 0 for all |x − c| ≤ r; hence x ± δ lie in the compact interval J := [c − r − δ, c + r + δ] ⊂ (0, ∞).

(24)

Define ϕ(x) := p(x − δ) − p(x + δ) for x ∈ [c − r, c + r]. Since p is strictly decreasing on (0, ∞) and x − δ < x + δ, we have ϕ(x) > 0 for every x in the domain. Moreover, ϕ is continuous because p is continuous. A continuous positive function on a compact set attains a positive minimum. Therefore, m :=

min

ϕ(x) > 0.

(25)

x∈[c−r,c+r]

Finally, |p(x + δ) − p(x − δ)| = p(x − δ) − p(x + δ) = ϕ(x) ≥ m > 0 ∀|x − c| ≤ r.

(26) (27)

Setting cδ (c, r) := m completes the proof. B.2. Proof of Theorem 4.2 Proof. Consider the partial derivative ∂Gθ,σ (µ)/∂µi . We first decompose this derivative into two terms as ∂Gθ,σ (µ) = ∂µi

Z

∂pµ,σ (x) dx ∂µi x∈S Z Z ∂pµ,σ (x) ∂pµ,σ (x) = g(θ, f (x)) dx + g(θ, f (x)) dx ∂µi ∂µi x∈cube(x∗ ;δ) x∈S\cube(x∗ ;δ) g(θ, f (x))

:=

∂H(µ) ∂R(µ) + , ∂µi ∂µi

(28) (29) (30)

where cube(x∗ ; δ) := {x|∀i ∈ [d], |xi − x∗i | ≤ δ} is a d dimensional cube centered at x∗ with 2δ length. Then, the main ∂R(µ) strategy of the proof is to show that the first term ∂H(µ) ∂µi dominates the second term ∂µi for sufficiently large θ, and then, the sign of ∂H(µ) ∂µi behaves like the statement. First, we bound ∂R(µ) ∗ ;δ) f (u). / ∂µi . Let Vδ := supu∈cube(x ∂R(µ) = ∂µi

Z g(θ, f (x)) x∈S\cube(x∗ ;δ)

∂pµ,σ (x) dx ≤ ∂µi

Z g(θ, f (x)) ∗ ;δ) x∈cube(x /

∂pµ,σ (x) dx ∂µi

Z ∂pµ,σ (x) ∂pµ,σ (x) dx ≤ g(θ, Vδ ) dx ∂µ ∂µi ∗ d i x∈cube(x / ;δ) x∈R Z d  x − µ  dx Z Y 1 ′  xi − µi  \i j j = g(θ, Vδ ) p p dxi d−1 2 σ σ σ x\i ∈Rd−1 xi ∈R σ

(31)

Z

≤ g(θ, Vδ )

(32) (33)

j̸=i

= g(θ, Vδ )

d Z Y j̸=i

Z = g(θ, Vδ )

p xj ∈R

 x − µ  dx Z j

j

σ

j

σ

1 ′  xi − µi  p dxi 2 σ xi ∈R σ

2g(θ, Vδ )p(0) 1 ′ p (z) dz ≤ , σ σ z∈R

where x\i indicates d − 1 dimensional variables except for xi . 12

(34) (35)

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization ∗ ∗ ′ Second, we bound ∂H(µ) ∂µi . Let Dδ := (f (x ) + Vδ )/2. Let ϵδ := f (x ) − Dδ . Then, there exists δ such that, for all ∗ ′ x ∈ cube(x , δ ),

f (x) > f (x∗ ) − ϵδ = Dδ > Vδ .

(36)

Then, under the condition |µi − x∗i | > δ in the statement, the sign of the partial derivative does not change in the cube, cube(x∗ ; δ). In fact, if µi < x∗i − δ, then for all x ∈ cube(x∗ ; δ), xi − µi = xi − x∗i + x∗i − µi > −δ + δ = 0 and ∂pµ,σ (x) ∂p (x) ≥ 0. If µi > x∗i + δ, then for all x ∈ cube(x∗ ; δ), xi − µi = xi − x∗i + x∗i − µi < δ − δ = 0 and µ,σ ≤ 0. ∂µi ∂µi Hence, for ∂H(µ) ∂µi , we have, ∂H(µ) = ∂µi

Z g(θ, f (x)) x∈cube(x∗ ;δ)

∂pµ,σ (x) dx = ∂µi

Z g(θ, f (x)) x∈cube(x∗ ;δ)

∂pµ,σ (x) dx ∂µi

Z ∂pµ,σ (x) ∂pµ,σ (x) dx ≥ g(θ, Dδ ) dx ∂µ ∂µi ∗ ′ ∗ ′ i x∈cube(x ;δ ) x∈cube(x ;δ )   d Z  x − µ  dx Z Y 1 ′  xi − µ i  j j j = g(θ, Dδ )  p dxi p 2 σ σ σ ∗ ′ ∗ ′ ′ ∗ ′ σ xi ∈[x∗ j̸=i xj ∈[xj −δ ,xj +δ ] i −δ ,xi +δ ]   Z d Z g(θ, Dδ ) Y ′ " ∗ # p (z) dz   ∗  p (z) dz =   x −µj −δ ′ x∗ −µj +δ ′ x −µi −δ ′ x∗ −µi +δ ′ j σ z∈ z∈ i , j , i

(37)

Z

g(θ, f (x))

j̸=i

=

σ

σ

σ

(38)

(39)

(40)

σ

  ∗   ∗   ∗  d   ∗ xj − µj − δ ′ xj − µj + δ ′ g(θ, Dδ ) Y xi − µ i − δ ′ xi − µ i + δ ′ −F −p F p σ σ σ σ σ

(41)

         d g(θ, Dδ ) Y z + δ′ z + δ′ z − δ′ z − δ′ F p − F inf − p inf ∗ ∗ z∈[x∗ z∈[x∗ σ σ σ σ σ i +M,xi −M ] j +M,xj −M ]

(42)

j̸=i

j̸=i

d

g(θ, Dδ ) Y · Cδ′ /σ (x∗j /σ, M/σ)cδ′ /σ (x∗i /σ, M/σ) σ

∵ Theorem 4.1

(43)

j̸=i

Hence, there exists the positive number θ such that the following inequality holds, d

g(θ, Dδ ) Y g(θ, Vδ ) · 2p(0) · Cδ′ /σ (x∗j /σ, M/σ)cδ′ /σ (x∗i /σ, M/σ) > . σ σ

(44)

j̸=i

Since Dδ > Vδ holds and g(θ, Dδ )/g(θ, Vδ ) is increasing due to the Assumption 4.4, we can always pick a sufficiently large θ that makes the above inequality hold. Choose θ large enough so that ∂H(µ) ∂R(µ) > . ∂µi ∂µi Then, by the decomposition

∂Gθ,σ (µ) ∂R(µ) = ∂H(µ) ∂µi ∂µi + ∂µi , we have

    ∂Gθ,σ (µ) ∂H(µ) sign = sign , ∂µi ∂µi because the perturbation term cannot flip the sign: indeed, ∂H(µ) ∂Gθ,σ (µ) ∂H(µ) ∂R(µ) >0 ⇒ ≥ − > 0, ∂µi ∂µi ∂µi ∂µi and similarly, ∂Gθ,σ (µ) ∂H(µ) ∂R(µ) ∂H(µ) <0 ⇒ ≤ + < 0. ∂µi ∂µi ∂µi ∂µi 13

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization ∂p

(x)

µ,σ ∗ It remains to identify the sign of ∂H(µ) ∂µi . On cube(x ; δ), the sign of ∂µi More precisely, since   ∂pµ,σ (x) 1 xi − µ i =− s pµ,σ (x), ∂µi σ σ

is constant under the condition |µi − x∗i | > δ. s(z) =

p′ (z) , p(z)

and s(z) < 0 for z > 0 while s(z) > 0 for z < 0 (because p is strictly decreasing on (0, ∞) and symmetric), we obtain: • If µi < x∗i − δ, then xi − µi > 0 for all x ∈ cube(x∗ ; δ), hence follows that ∂H(µ) ∂µi ≥ 0 (in fact > 0), and therefore

∂Gθ,σ (µ) > 0. ∂µi

• If µi > x∗i + δ, then xi − µi < 0 for all x ∈ cube(x∗ ; δ), hence < 0), and therefore

∂pµ,σ (x) ≥ 0 on the cube. Since g(θ, f (x)) > 0, it ∂µi

∂pµ,σ (x) ≤ 0 on the cube. Thus ∂H(µ) ∂µi ∂µi ≤ 0 (in fact

∂Gθ,σ (µ) < 0. ∂µi

This matches the claimed invariance of the sign outside cube(x∗ ; δ) and completes the proof. B.3. Proof of Theorem 4.3 Proof. Under Assumptions 4.1–4.1 and 4.3–4.4, we have the following two properties: (1) Lipschitz continuity. Gθ,σ is Lipschitz (hence continuous) on Rd . |Gθ,σ (µ1 ) − Gθ,σ (µ2 )| = |∇µ Gθ,σ (ν)⊺ (µ1 − µ2 )|

√ g∗,θ dI ≤ ∥∇µ Gθ,σ (ν)∥2 ∥µ1 − µ2 ∥2 ≤ ∥µ1 − µ2 ∥2 σ

(45) (46)

(2) Vanishing at infinity. lim∥µ∥→∞ Gθ,σ (µ) = 0 and there exists µ0 with Gθ,σ (µ0 ) > 0. Since S ⊂ Rd is compact, it is bounded. Assume that S ⊂ cube(0, M ) for some M > 0. Then, g(θ, f (x)) = 0 if ∥x∥∞ > M . Hence, Gθ,σ (µ) ≤ g∗,θ

d Z Y j=1

= g∗,θ

p xj ∈[−M,M ]

 x − µ  dx j

j

σ

j

(47)

σ

   d   Y M − µj −M − µj F −F σ σ j=1

→ 0 as ∥µ∥∞ → ∞,

∵ For any fixed δ,

lim F (x) − F (x − δ) = 0

x→±∞

(48) (49)

Given (1) and (2), the standard compactness-continuity argument (as in Chen Xu’s case) applies verbatim: choose a closed ball where the interior value exceeds the uniform tail bound, invoke the extreme value theorem to attain a maximizer on the ball, and compare with the exterior to conclude global maximality. We omit further details. B.4. Proof of Lemma 4.4 Proof. The gradient can be expressed in score form as   ∇µ Gθ,σ (µ) = E g(θ, f (X)) ∇µ log pµ,σ (X); S .

(50)

∂ 1 log pµ,σ (x) = − s(zi ), ∂µi σ

(51)

For the ith coordinate,

(z) i where zi = xi −µ and s(z) = pp(z) . Differentiating once more with respect to µi gives σ

 ∂2 1  p′′ (zi ) log pµ,σ (x) = 2 − s(zi )2 . 2 ∂ µi σ p(zi )

(52)

Hence the diagonal entries of the Hessian H(µ) = ∇2µ Gθ,σ (µ) are bounded in absolute value by g∗,θ K/σ 2     ∇2µ Gθ,σ (µ) =E g(θ, f (X)) ∇µ log pµ,σ (X) (∇µ log pµ,σ (X))⊺ ; S + E g(θ, f (X)) ∇2µ log pµ,σ (X); S .

14

(53)

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

Furthermore, we have,   1  p′′ (zi ) 1 [∇2µ Gθ,σ (µ)]ii = E g(θ, f (X)) · 2 s(zi )2 + g(θ, f (X)) · 2 − s(zi )2 σ σ p(zi )   Z p′′ (zi ) g∗θ g∗θ K = E g(θ, f (X)) · 2 . ≤ 2 |p′′ (z)| dz = σ p(zi ) σ σ2

(54) (55)

For the off-diagonal entries i ̸= j,

g∗,θ [∇2µ Gθ,σ (µ)]ij ≤ 2 σ

since

R R

|p′ (z)|dz

2

R R

Z

Z

2

g∗,θ I, σ2

(56)

  g∗,θ g∗,θ K + (d − 1)I ≤ 2 d max K, I . 2 σ σ

(57)

g∗,θ |s(zi )s(zj )|pµ,σ (x)dx = 2 σ d R

|p (z)|dz

R

R (p′ (z))2 p(z) dz R p(z)dz = I by using Cauchy–Schwarz.

Since H(µ) is symmetric, we use

λmax (H) ≤ ∥H∥∞ = max i

d X

|Hij | ≤

j=1

By the fundamental theorem of calculus, Z 1 ∥∇µ Gθ,σ (µ) − ∇µ Gθ,σ (ν)∥2 =

H(ν + t(µ − ν))(µ − ν) dt 0

(58) 2

≤ ∥H(ν + t(µ − ν))∥op ∥(µ − ν)∥2 2

≤ g∗,θ d max(K, I)/σ ∥(µ − ν)∥2 ,

(59) (60)

which proves the claim. B.5. Proof of Lemma 4.5 Proof. Define Y = g(θ, f (X)) ∇µ log p(X) with mean m = E[Y ] = ∇µ Gθ,σ (µ). For B i.i.d. copies Y1 , . . . , YB , independence gives E

1 B

B X

2

Yk

∥m∥2 . = B1 E∥Y ∥2 + B−1 B

(61)

k=1

Pd To bound E∥Y ∥2 , note that ∥∇µ log p(X)∥2 = σ −2 i=1 s(Zi )2 with Zi = (Xi − µi )/σ and s(z) = p′ (z)/p(z). Taking 2 expectations yields E∥∇µ log p(X)∥2 = (dI)/σ 2 . Since g(θ, f (X)) ≤ g∗,θ , it follows that E∥Y ∥2 ≤ g∗,θ (dI)/σ 2 . Moreover, ∥m∥2 ≤ E∥Y ∥2 by Cauchy-Schwarz, so the same bound holds for ∥m∥2 . Substituting back into the decomposition gives the uniform bound   2 dI ˆ µ Gθ,σ (µ)∥2 ≤ g∗,θ E ∥∇ . σ2

(62)

B.6. Proof of Theorem 4.6 Proof. By Lemma 4.4, Gθ,σ is Lθ -smooth, so for all x, y ∈ Rd , Gθ,σ (y) ≥ Gθ,σ (x) + ∇µ Gθ,σ (x)⊤ (y − x) −

Lθ ∥y − x∥2 . 2

(63)

b µ Gθ,σ (µt ) to obtain Apply (63) with x = µt and the update y = µt+1 = µt + ηt ∇ 2

b µ Gθ,σ (µ ) − Lθ ηt Gθ,σ (µt+1 ) ≥ Gθ,σ (µt ) + ηt ∇µ Gθ,σ (µt )⊤ ∇ t 2

15

b µ Gθ,σ (µ ) 2 . ∇ t

(64)

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

b µ Gθ,σ (µt ) | µt ] = ∇µ Gθ,σ (µt ) converts Taking conditional expectation given µt and using the score unbiasedness E[∇ 2 the mixed inner product into ηt ∥∇µ Gθ,σ (µt )∥ . Taking total expectation and invoking the uniform bound from Lemma 4.5 yields       Lθ ηt2 E Gθ,σ (µt+1 ) −E Gθ,σ (µt ) ≥ ηt E ∥∇µ Gθ,σ (µt )∥2 − Qθ . 2

(65)

Summing (65) over t = 0, . . . , T − 1 telescopes the left-hand side to E[Gθ,σ (µT )] − E[Gθ,σ (µ0 )]. Hence −1       Lθ Qθ TX ηt2 . ηt E ∥∇µ Gθ,σ (µt )∥2 ≤ E Gθ,σ (µT ) − E Gθ,σ (µ0 ) + 2 t=0 t=0

T −1 X

(66)

B.7. Proof of Corollary 4.7 Proof. The localization claim follows from Theorem 4.2: for each δ > 0 and M > 0 there exists θδ,σ,M such that the coordinatewise derivatives have fixed signs outside cube(x∗ ; δ) within {∥µ∥∞ < M }, hence no stationary point (and therefore no local maximizer) can lie outside cube(x∗ ; δ). For the rate, start from Theorem 4.6. With the polynomial steps, ∞ X

1 −1  TX σ 2 (2 2 −γ − 1) 21 −γ σ 4 1  T 1 + , η ≥ . t 2γ d2 d 12 − γ t=0

ηt2 ≤

t=0

(67)

2 From Gθ ≤ (dI/σ 2 ) g∗,θ and Lθ ≤ g∗,θ (d max(K, I))/σ 2 , we have



min

0≤t≤T −1

2

E ∥∇µ Gθ,σ (µt )∥

g∗,θ + 21 Lθ Gθ P ≤ t ηt

2 t ηt

P

3 g∗,θ + I max(K, I) g∗,θ 1

σ 2 T 2 −γ

· d · Cγ ,

(68)

1

where Cγ := ( 21 − γ)/(2 2 −γ − 1) which implies the stated T bound after inversion. Finally, we get the iteration complexity 2  1−2γ 3 g∗,θ + I max(K,I) g∗,θ T > Cγ d · . 2 σ ε B.8. Proof of Corollary 4.8  Qd Proof. We work in the single-kernel, anisotropic setting pµ,Σ (x) = i=1 σi−1 p (xi − µi )/σi . It is enough to derive the Lipschitz constant and the second moment bound. For the Lipschitz constant, let H(µ) := ∇2µ Gθ,Σ (µ) and use the induced 2-operator norm. For the diagonal Hessian entries, ∂µ2 i µi pµ,Σ (x) =

1 p′′ (zi ) pµ,Σ (x), σi2 p(zi )

(69)

and therefore ∂µ2 i µi Gθ,Σ (µ) ≤

g∗,θ σi2

Z

|p′′ (z)| dz =

R

g∗,θ K. σi2

(70)

For the off-diagonal entries i ̸= j, ∂µ2 i µj pµ,Σ (x) =

1 s(zi )s(zj )pµ,Σ (x), σi σj

(71)

so that ∂µ2 i µj Gθ,Σ (µ) ≤

g∗,θ σi σj

Z

|p′ (z)| dz R

16

2 ≤

g∗,θ I. σi σj

(72)

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

Thus,  X I K . ∥H(µ)∥op ≤ ∥H(µ)∥∞ = max |Hij (µ)| ≤ g∗,θ max  2 + 1≤i≤d 1≤i≤d σi σi σj j=1

(73)

 X I K . Lθ (Σ) ≤ g∗,θ max  2 + 1≤i≤d σi σi σj

(74)

S2 (Σ) := d max σi−2 .

(75)

1 ≤ max σℓ−2 , 1≤ℓ≤d σi σj

(76)

X I K + ≤ d max(K, I) max σℓ−2 . 1≤ℓ≤d σi2 σi σj

(77)

Lθ (Σ) ≤ g∗,θ max(K, I) S2 (Σ).

(78)

h i   2 ˆ µ Gθ,Σ (µ)∥2 ≤ E ∥g(θ, f (X))∇µ log pµ,Σ (X)∥2 ≤ g∗,θ E ∥∇ I S2 (Σ).

(79)

d X

j̸=i

Hence 

j̸=i

Now define 1≤i≤d

Since K ≤ max(K, I), I ≤ max(K, I), and

we obtain

j̸=i

Therefore

For the second moment bound,

2 Hence Qθ (Σ) = g∗,θ I S2 (Σ). The remaining steps follow exactly as in Corollary 4.7.

B.9. Proof of Lemma 4.9 Proof. Given µt , bk is a function of {(X(j) )}j̸=k and thus independent of X(k) . Using E[S(X(k) ) | µt ] = 0, h i h i h i  E h(X(k) ) − bk S(X(k) ) | µt = E h(X(k) ) S(X(k) ) | µt − E bk S(X(k) ) | µt h i h i = E h(X(k) ) S(X(k) ) | µt − E[bk | µt ] E S(X(k) ) | µt {z } | =0 i h = E h(X(k) ) S(X(k) ) | µt .

Averaging over k yields the claim. B.10. Proof of Lemma 4.10 Proof. Let Yk = (h(X(k) ) − bk ) S(X(k) ). Then, B h i    P 1 X  2 E B1 k Yk µt ≤ E ∥Yk ∥2 µt = E ∥Y1 ∥2 µt . B k=1

Since ∥Y1 ∥2 = (h(X(1) ) − b1 )2 ∥S(X(1) )∥2 , we introduce an auxiliary constant baseline b⋆ and write   h(X(1) ) − b1 = h(X(1) ) − b⋆ + b⋆ − b1 . 17

(80) (81)

(82)

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

Expanding and taking conditional expectation yields h i h i   E ∥Y1 ∥2 µt = E (h(X(1) ) − b⋆ )2 ∥S(X(1) )∥2 µt + E (b⋆ − b1 )2 ∥S(X(1) )∥2 µt {z } | {z } | =:A =:B h i + 2E (b⋆ − b1 )(h(X(1) ) − b⋆ )∥S(X(1) )∥2 µt .

By construction, b1 is independent of X(1) , and hence independent of (h(X(1) ) − b⋆ )∥S(X(1) )∥2 . Moreover, b⋆ is defined as the optimal constant baseline minimizing E[(h(X) − b)2 ∥S(X)∥2 |µt ], which implies the first-order optimality condition h i E (h(X(1) ) − b⋆ )∥S(X(1) )∥2 |µt = 0. Therefore, h i E (b⋆ − b1 )(h(X(1) ) − b⋆ )∥S(X(1) )∥2 |µt = 0, and the final term vanishes after taking conditional expectation. Consequently, it suffices to analyze the two terms A and B 2 and derive explicit upper bounds for each. For the first term A, by the definition of Rθ,µ,Σ and g(θ, f ) ≤ g∗,θ , 2 2 2 A = (1 − Rθ,µ ) E[h(X)2 ∥S(X)∥2 ] ≤ (1 − Rθ,µ ) g∗,θ I S2 (Σ). t ,Σ t ,Σ

(83)

For the second term B, note that bk is independent of X(k) conditional on µt , so E[(b⋆ − bk )2 ∥S(X)∥2 | µt ] = E[(b⋆ − bk )2 | µt ] E[∥S(X)∥2 | µt ] ≤ E[(b⋆ − bk )2 ] I S2 (Σ).

(84)

We now bound E[(b⋆ − bk )2 ] using the mean value theorem for the two-variable map g(u, v) = u/(v + λ). Decompose b⋆ − b k =



µU µV

   U =: δ0 − ∆, − µVµU+λ + µVµU+λ − V \k +λ

(85)

\k

λ so that δ0 = µbV +λ captures the deterministic ridge bias and ∆ is the random fluctuation around (µU , µV ).

Let Ut = µU + t∆u and Vt = µV + t∆v with ∆u = U\k − µU and ∆v = V\k − µV . By the (integral) mean value representation, Z 1 ∆ = g(U\k , V\k ) − g(µU , µV ) =

∇g(Ut , Vt ) · (∆u , ∆v ) dt.

(86)

0

 Since ∇g(u, v) = (v + λ)−1 , −u(v + λ)−2 , it follows that |∆| ≤

Z 1 0

|∆u | |Ut | |∆v |  |∆u | (µU + |∆u |) |∆v | + dt ≤ + , Vt + λ (Vt + λ)2 λ λ2

(87)

where we used Vt + λ ≥ λ and |Ut | ≤ µU + |∆u |. Squaring and using (x + y)2 ≤ 2x2 + 2y 2 gives ∆2 ≤

2(µU + |∆u |)2 ∆2v 2 ∆2u 4µ2 ∆2 4 ∆2u ∆2v 2 ∆2u + ≤ + U4 v + . 2 4 2 λ λ λ λ λ4

(88)

Taking expectations and using the sample-mean scalings E[∆2u ] = Var(U )/s, E[∆2v ] = Var(V )/s, and E[∆2u ∆2v ] ≤ (E∆4u )1/2 (E∆4v )1/2 = O(s−2 ), we obtain E[∆2 ] ≤

 1  2 Var(U ) 4µ2 Var(V ) · + 4U · +O 4 2 . 2 λ s λ s λ s

(89)

Finally, E[(b⋆ − bk )2 ] = δ02 − 2δ0 E[∆] + E[∆2 ] ≤



b⋆ λ µV + λ

2

2 Std(U ) 2µU Std(V ) √ √ + +O λ s λ2 s  1  4µ2U Var(V ) 2 Var(U ) + + O + , λ2 s λ4 s λ4 s2 +



1 λ2 s



which is a mean-value-theorem alternative to the Taylor-based bound and has the same rate profile in s = B − 1. 18

(90)

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

B.11. Proof of Theorem 4.11 Proof. The claim follows directly from Lemma 4.10 with an appropriate choice of the ridge parameter. From (17), for any λ > 0 we have     2 loo 2 2 b E ∇µ Gθ,Σ (µt ) µt ≤ (1 − Rθ,µ ) g∗,θ + Cloo (B, λ) I S2 (Σ). t ,Σ

(91)

We choose λ = (B − 1)−1/8 . Substituting this choice into the expression of Cloo (B, λ) in (18) shows that each term in Cloo (B, λ) decays at rate O((B − 1)−1/4 ), and therefore Cloo (B, λ) = O(B −1/4 ). This yields the first inequality stated in the corollary. Taking expectation of (91) with respect to µt gives 

bµloo Gθ,Σ (µ ) ∇ t

E

2



   2 2 ≤ g∗,θ 1 − E[Rθ,µ ] + C (B, λ) I S2 (Σ). loo ,Σ t

(92)

2 By definition of Cθ,Σ , we have E[Rθ,µ ] ≥ Cθ,Σ for all t, and hence t ,Σ

 E

bµloo Gθ,Σ (µ ) ∇ t

2



  2 ≤ g∗,θ (1 − Cθ,Σ ) + Cloo (B, λ) I S2 (Σ).

(93)

Since Cloo (B, λ) = O(B −1/4 ), there > 0 such that Cloo (B, λ) ≤ c B −1/4 . Under the  exists a numerical constant c −1/4 2 4 2 assumption B = Ω 16/(Cθ,Σ g∗,θ ) , the remainder term satisfies c B ≤ 12 Cθ,Σ g∗,θ . Substituting this bound into the previous inequality yields     2 Cθ,Σ 2 bµloo Gθ,Σ (µ ) E ∇ ≤ 1 − g∗,θ I S2 (Σ), t 2

(94)

which proves the claim. B.12. Proof of Corollary 4.12 Proof. The one-step expected improvement inequality in Theorem 4.6 holds with the same smoothness constant Lθ (Σ) = 2 g∗,θ max(K, I) S2 (Σ) from Lemma 4.4. Replacing the plain score second moment Qθ (Σ) = g∗,θ IS2 (Σ) by the LOO constant from Theorem 4.11, we have   (loo) C −1/4 2 Qθ (Σ) = (1 − θ,Σ ) + O(B ) g∗,θ I S2 (Σ). 2 (loo)

Plugging Lθ (Σ) and Qθ (Σ) into the telescoping argument of Theorem 4.6, and then using the polynomial step schedule PT −1 P∞ 1 1 1 ηt = S2 (Σ)−1 (t + 1)−( 2 +γ) (whose sums satisfy t=0 ηt ≥ Cγ S2 (Σ)−1 T 2 −γ and t=0 ηt2 ≤ (1 + 2γ )S2 (Σ)−2 ).

C. Experimental Detail This appendix provides full experimental details for the canonical optimization benchmarks and black-box targeted adversarial attack experiments reported in the main text, including hyperparameter selection procedures, optimization settings, and dataset-specific configurations. C.1. Canonical Benchmark Objective Functions Ackley. v  u D u1 X Ackley(x) = 20 exp −0.2t x2d  + exp D 

d=1

which attains its global optimum value of 0 at x = 0. 19

D

1 X cos(2πxd ) D d=1

! − 20 − e,

(95)

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization

Rosenbrock. D−1

Rosenbrock(x) = −

 1 X 100(xd+1 − x2d )2 + (1 − xd )2 , D−1

(96)

d=1

which attains its global optimum value of 0 at x = 1. Griewank. Griewank(x) = − 1 −

D

D

d=1

d=1

1 X 2 Y xd + 1.050.2 cos 4000



x √d d

 ,

(97)

which attains its global optimum value of 1.050.2D − 1 at x = 0. C.2. Hyperparameters of Canonical Benchmarks Hyperparameters for the benchmark test functions were determined by averaging results from 20 trajectories initialized at different points. For each configuration, we calculated the average of the best MSE values obtained across these 20 runs and selected the one with the lowest overall mean. These selected values are reported in the table 9, 10, and 11. The experimental environment was fixed with dimension D = 500, the total number of optimization steps per run was set to T = 400, and the number of Monte Carlo samples for gradient estimation was set to B = 50. Table 9. Hyperparameters for Optimizing Ackley. The initial points for the Ackley function were sampled from a multivariate normal distribution N (51, 0.012 I). The candidate set for initial learning rate (η0 ) is L := {0.1, 0.5}, and for initial smoothing scale (σ) is S := {0.1, 0.5}. The candidate set for the amplification parameter (θ) is A := {1, 3, 5}. For both ProMoT and ProMoT-loo, a Logistic kernel was employed as the smoothing kernel, with the transformation function defined as g(θ, y) = (y + c)β eθy (where c = 600 and β = 10). In RSGF and ZO-SLGHd/r, γ denotes the decreasing factor for the smoothing scale σ, while in the baseline ZO-SLGHd, α represents the step size used to update the smoothing scale.

Method ProMoT ProMoT-loo EPGS RSGF ZO-SGD ZO-AdaMM ZO-SLGHd ZO-SLGHr CMA-ES

Selected Values η0 = 0.5, σ = 0.5, θ = 5 η0 = 0.5, σ = 0.1, θ = 1 η0 = 0.5, σ = 0.5, θ = 5 η0 = 0.1, σ = 0.1, γ = 0.9 η0 = 10, σ = 0.1 η0 = 0.5, σ = 0.1, β1 = 0.5, β2 = 0.5 η0 = 5, σ = 0.5, γ = 0.99, α = 0.001 η0 = 5, σ = 0.5, γ = 0.99 σ = 0.5

Candidates η0 ∈ L, σ ∈ S, θ ∈ A η0 ∈ L, σ ∈ S, θ ∈ A η0 ∈ L, σ ∈ S, θ ∈ A η0 ∈ L, σ ∈ S, γ ∈ {0.8, 0.9, 0.99} η0 ∈ {1, 5, 10}, σ ∈ S η0 ∈ L, σ ∈ S, β1 ∈ {0.5, 0.7, 0.9}, β2 ∈ {0.1, 0.3, 0.5} η0 ∈ {1, 5, 10}, σ ∈ S, γ ∈ {0.95, 0.99}, α ∈ {0.1, 0.01, 0.001} η0 ∈ {1, 5, 10}, σ ∈ S, γ ∈ {0.95, 0.99} σ∈S

Table 10. Hyperparameters for Optimizing Rosenbrock. The initial points for the Rosenbrock function were sampled from a multivariate normal distribution N (31, 0.012 I). The candidate set for initial learning rate (η0 ) is L := {0.1, 0.5}, and for initial smoothing scale (σ) is S := {0.1, 0.5}. The candidate set for the amplification parameter (θ) is A := {0.1, 0.01, 0.001}. For both ProMoT and ProMoT-loo, a Logistic kernel was employed as the smoothing kernel, with the transformation function defined as g(θ, y) = (y + c)β eθy (where c = 6, 000 and β = 10). In RSGF and ZO-SLGHd/r, γ denotes the decreasing factor for the smoothing scale σ, while in the baseline ZO-SLGHd, α represents the step size used to update the smoothing scale.

Method ProMoT ProMoT-loo EPGS RSGF ZO-SGD ZO-AdaMM ZO-SLGHd ZO-SLGHr CMA-ES

Selected Values η0 = 0.1, σ = 0.1, θ = 0.1 η0 = 0.1, σ = 0.1, θ = 0.001 η0 = 0.1, σ = 0.1, θ = 0.1 η0 = 0.1, σ = 0.5, γ = 0.99 η0 = 0.01, σ = 0.1 η0 = 0.5, σ = 0.1, β1 = 0.3, β2 = 0.9 η0 = 0.001, σ = 0.5, γ = 0.95, α = 0.1 η0 = 0.001, σ = 0.5, γ = 0.95 σ = 0.1

Candidates η0 ∈ L, σ ∈ S, θ ∈ A η0 ∈ L, σ ∈ S, θ ∈ A η0 ∈ L, σ ∈ S, θ ∈ A η0 ∈ L, σ ∈ S, γ ∈ {0.9, 0.95, 0.99} η0 ∈ {0.01, 0.001}, σ ∈ S η0 ∈ L, σ ∈ S, β1 ∈ {0.5, 0.7, 0.9}, β2 ∈ {0.1, 0.3, 0.5} η0 ∈ {0.01, 0.001}, σ ∈ S, γ ∈ {0.95, 0.99}, α ∈ {0.1, 0.01, 0.001} η0 ∈ {0.01, 0.001}, σ ∈ S, γ ∈ {0.95, 0.99} σ∈S

20

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization Table 11. Hyperparameters for Optimizing Griewank. The initial points for the Griewank function were sampled from a multivariate normal distribution N (51, 0.012 I). The candidate set for initial learning rate (η0 ) is L := {0.1, 0.5}, and for initial smoothing scale (σ) is S := {1, 2}. The candidate set for the amplification parameter (θ) is A := {1, 3, 5}. For both ProMoT and ProMoT-loo, a Logistic kernel was employed as the smoothing kernel, with the transformation function defined as g(θ, y) = (y + c)β eθy (where c = 1, 000 and β = 10). In RSGF and ZO-SLGHd/r, γ denotes the decreasing factor for the smoothing scale σ, while in the baseline ZO-SLGHd, α represents the step size used to update the smoothing scale.

Method ProMoT ProMoT-loo EPGS RSGF ZO-SGD ZO-AdaMM ZO-SLGHd ZO-SLGHr CMA-ES

Selected Values η0 = 0.1, σ = 2, θ = 5 η0 = 0.1, σ = 1, θ = 1 η0 = 0.1, σ = 2, θ = 5 η0 = 0.1, σ = 1, γ = 0.9 η0 = 10, σ = 1 η0 = 0.5, σ = 1, β1 = 0.1, β2 = 0.7 η0 = 10, σ = 1, γ = 0.99, α = 0.1 η0 = 10, σ = 2, γ = 0.95 σ=2

Candidates η0 ∈ L, σ ∈ S, θ ∈ A η0 ∈ L, σ ∈ S, θ ∈ A η0 ∈ L, σ ∈ S, θ ∈ A η0 ∈ L, σ ∈ S, γ ∈ {0.8, 0.9, 0.99} η0 ∈ {1, 5, 10}, σ ∈ S η0 ∈ L, σ ∈ S, β1 ∈ {0.5, 0.7, 0.9}, β2 ∈ {0.1, 0.3, 0.5} η0 ∈ {1, 5, 10}, σ ∈ S, γ ∈ {0.95, 0.99}, α ∈ {0.1, 0.01, 0.001} η0 ∈ {1, 5, 10}, σ ∈ S, γ ∈ {0.95, 0.99} σ∈S

C.3. Experimental Settings for Black-Box Targeted Adversarial Attacks We first describe the black-box attack objective shared across all adversarial experiments, and then provide dataset-specific model architectures, preprocessing procedures, and optimization hyperparameters. Attack objective. Let C denote a black-box classifier and x an input sample. A targeted adversarial attack seeks a perturbation µ such that arg max C(x + µ)y = ytgt , y

where the target label is chosen as the most unlikely class ytgt = arg miny C(x)y . We optimize the C&W-style black-box objective (Carlini & Wagner, 2017; Xu, 2025)   L(µ) = − max max C(x + µ)y − C(x + µ)ytgt , κ − λ∥µ∥2 , (98) y̸=ytgt

which encourages the target logit to exceed all non-target logits by a margin κ while penalizing the perturbation magnitude. Once the margin constraint is satisfied, optimization focuses on minimizing ∥µ∥2 , leading to less perceptible adversarial examples. CIFAR-10. We evaluated the proposed methods and baselines against a CNN-based classifier trained on the CIFAR-10 dataset (Krizhevsky et al., 2009). Specifically, the target model’s architecture, following the configurations of the model in Tables 1 and 2 of (Carlini & Wagner, 2017), was implemented using the PyTorch framework. To evaluate the effectiveness of our methods and baselines against a robust target, defensive distillation (Papernot et al., 2016) was applied during the training process. After training, the model achieved a test accuracy of 75.7%, serving as a representative target for our adversarial attack evaluations. Detailed hyperparameter configurations for these experiments are listed in Table 12. The input dimension was set to D = 3, 072 (corresponding to 32 × 32 × 3 pixels), and we limited the total number of optimization steps per run to T = 500. For gradient estimation, we employed B = 30 Monte Carlo samples. VitalDB. We used the case-level clinical table (cases) from the VitalDB dataset (Lee et al., 2022) and defined death inhosp as the label. We binarized sex and removed categorical fields as well as metadata identifiers and timestamps. Missing values were imputed as follows: (i) intraop ebl, intraop uo, and intraop crystalloid were set to 0; (ii) arterial blood gas variables were filled with normal reference values; and (iii) all remaining features were imputed using the feature-wise median. We applied a log(1 + x) transform to skewed lab and medication variables and applied a Yeo–Johnson transform to preop be. Finally, min–max scaling was used to map all features to [−2, 2]. This resulted in a final input dimension of D = 42. All preprocessing transforms were fitted on the training split only and then applied to the test split. The XGBoost classifier was trained on a highly imbalanced dataset consisting of 4,791 training samples and 1,597 test samples, with an event rate of only 0.9%. Despite this extreme class imbalance, the model achieved an AUROC of 0.843 and 21

Probabilistic Smoothing with Ratio-Monotone Transforms for Global Optimization Table 12. Hyperparameters for CIFAR-10 Attack. In this experiment, η0 , σ, and θ denote the initial learning rate, the initial smoothing scale, and the amplification parameter, respectively. For the smoothing kernel, ProMoT employed a Gaussian kernel, while ProMoT-loo utilized a Generalized Gaussian kernel. Both algorithms shared the same transformation function, defined as g(θ, y) = (y + c)β eθy (where c = 10, 000 and β = 10). Following the baseline configurations, γ represents the decreasing factor for the smoothing scale σ in RSGF and ZO-SLGHd/r, and α denotes the step size used for updating the smoothing scale in ZO-SLGHd.

Method ProMoT ProMoT-loo EPGS RSGF ZO-SGD ZO-AdaMM ZO-SLGHd ZO-SLGHr CMA-ES

Hyperparameters η0 = 0.005, σ = 0.1, θ = 0.03 η0 = 0.005, σ = 0.01, θ = 0.03 η0 = 0.005, σ = 0.1, θ = 0.03 η0 = 0.03, σ = 1.0, γ = 0.8 η0 = 0.00005, σ = 0.1 η0 = 0.03, σ = 0.1 ,β1 = 0.9, β2 = 0.1 η0 = 0.00003, σ = 0.1, γ = 0.999, α = 0.1/3072 η0 = 0.00003, σ = 0.1, γ = 0.995 σ = 0.05

an AUPRC of 0.535, with an overall accuracy of 91.1%. At an operating threshold of 0.011, the model yielded an F1-score of 0.113, corresponding to 1,446 TN, 137 FP, 5 FN, and 9 TP. For adversarial optimization, we limited the total number of attack iterations per run to T = 500 and used B = 30 Monte Carlo samples for gradient estimation. Detailed hyperparameter settings for the adversarial optimization are provided in Table 13. Table 13. Hyperparameters for VitalDB Attack. In this experiment, η0 , σ, and θ denote the initial learning rate, the initial smoothing scale, and the amplification parameter, respectively. Both ProMoT and ProMoT-loo used the same Gaussian smoothing kernel and shared an identical transformation function defined as g(θ, y) = (y + c)β eθy , where c = 1, 000 and β = 10. Following the baseline configurations, γ represents the decreasing factor for the smoothing scale σ in RSGF and ZO-SLGHd/r, and α denotes the step size used for updating the smoothing scale in ZO-SLGHd.

Method ProMoT ProMoT-loo EPGS RSGF ZO-SGD ZO-AdaMM ZO-SLGHd ZO-SLGHr CMA-ES

Hyperparameters η0 = 0.01, σ = 0.3, θ = 5 η0 = 0.01, σ = 0.3, θ = 5 η0 = 0.03, σ = 0.3, θ = 5 η0 = 0.03, σ = 3.0, γ = 0.95 η0 = 0.3, σ = 0.3 η0 = 0.5, σ = 0.3 ,β1 = 0.7, β2 = 0.3 η0 = 0.1, σ = 1.0, γ = 0.999, α = 0.001 η0 = 0.1, σ = 0.3, γ = 0.995 σ = 0.5

22

Record · ID 229482 · SHA-256 3b3f4bc337370bf4
Retrieved via Conceptio — every document is proof-bundled with source, license, and retrieval metadata.