Sharpness-Aware Black-Box Optimization

## Abstract

Abstract Black-box optimization algorithms have been widely used in various machine learning problems, including reinforcement learning and prompt fine-tuning. However, directly optimizing the training loss value, as commonly done in existing black-box optimization methods, could lead to suboptimal model quality and generalization performance.
To address those problems in black-box optimization, we propose a novel Sharpness-Aware Black-box Optimization (SABO) algorithm, which applies a sharpness-aware minimization strategy to improve the model generalization. Specifically, the proposed SABO method first reparameterizes the objective function by its expectation over a Gaussian distribution.
Then it iteratively updates the parameterized distribution by approximated stochastic gradients of the maximum objective value within a small neighborhood around the current solution in the Gaussian distribution space.
Theoretically, we prove the convergence rate and generalization bound of the proposed SABO algorithm.
Empirically, extensive experiments on the black-box prompt fine-tuning tasks demonstrate the effectiveness of the proposed SABO method in improving model generalization performance.

## 1 Introduction

Black-box optimization involves optimizing one objective function by using function queries only. In this work, we study the black-box optimization problem [], which is formulated as

$$ $\min_{\bm{x}}F(\bm{x}),\ ~{}\mathrm{s.t.}\ ~{}\bm{x}\in\mathcal{X},$ (1) $$

where $\mathcal{X}\subseteq\mathbb{R}^{d}$, and $d$ represents the parameter dimension. The objective function $F:\mathbb{R}^{d}\to\mathbb{R}$, which satisfies $F(\bm{x})\geq-\infty$ (lower bounded), can only be queried to obtain function values and we cannot get the gradient of $F$ w.r.t. $\bm{x}$. In this work, we focus on the online setting for black-box optimization, where different from the offline setting [, ], we do not have a prior dataset containing the variable $\bm{x}$ and its corresponding objective value.

Black-box optimization has drawn intensive attention in a wide range of applications, such as deep reinforcement learning [, ], black-box adversarial attacks of deep neural networks [, ], etc. Recently, black-box optimization has shown increasing power on real-world natural language processing tasks, especially with the emergence of large language models (LLMs). Since a common practice is to release LLMs as a service and allow users to access it through their inference APIs. In such a scenario, called Languaged-Model-as-a-Service (LMaaS) [, , ], users cannot access or tune model parameters but can only tune their prompts without model backpropagation to accomplish language tasks of interest, which directly increase the demand of black-box optimization methods.

Although black-box optimization algorithms have been successfully applied to various learning tasks, most existing works directly optimize the training loss value [], which may lead to suboptimal model quality and generalization performance.
Since the training loss landscape is complex and has many local minima with different generalization abilities [], the learned model may suffer from the overfitting problem, causing poor generalization performance []. Hence, reducing the performance gap between training and testing is an important research topic in deep learning []. Recently, there have been many works exploring the close relationship between loss geometry and generalization performance, and it has been observed that flat minima often imply better generalization [, , ]. This inspires us to design a black-box optimization algorithm to improve the model generalization by finding the flat minima.

Sharpness-aware minimization (SAM) [] is a state-of-the-art method to seek flat minima in white-box cases by solving a min-max optimization problem. SAM minimizes the maximum objective value within a small neighborhood around current solution. Since SAM considers the geometry of the Euclidean parameter space, it uses the Euclidean ball to define the neighborhood. In SAM, each update consists of two forward-backward computations: one for computing the perturbation and the other for computing the actual update direction. SAM has been proven to perform better than SGD and its variants [, , , , , ] yield significant performance gains in various fields such as computer vision and natural language processing [, ]. However, SAM and its variants rely on the availability of true gradients or stochastic gradients w.r.t. the variable $\bm{x}$, so they are inapplicable to black-box optimization.

To take advantage of SAM to improve the generalization performance of black-box optimization, we propose a Sharpness-Aware Black-box Optimization (SABO) algorithm. Specifically, SABO first reparameterizes the objective function via its expectation over a Gaussian distribution, which can help to optimize the objective by only accessing the function value [, ]. Then the SABO method seeks to identify the robust minimum region over the space of Gaussian distributions, which is different from SAM that finds the flat minimum over the parameter space. To achieve that, the SABO method iteratively updates the parameterized distribution via a search direction obtained by approximated stochastic gradients for the maximum objective value within a small neighborhood around the current solution in the space of Gaussian distributions.
Theoretically, we analyze the convergence rate and provide a generalization error bound for the proposed SABO algorithm.
Empirically, we verify the convergence result of the proposed algorithm on the synthetic problems, and extensive experimental results on a black-box prompt fine-tuning problem demonstrate the effectiveness of the proposed SABO method.
Our contributions are summarized as follows.

- •
We propose the SABO algorithm for black-box optimization. To the best of our knowledge, we are the first to design a stochastic gradient approximation algorithm to improve the model generalization in black-box optimization by using the sharpness-aware minimization strategy.
- •
Theoretically, we prove that the proposed SABO algorithm possesses a convergence rate $\mathcal{O}(\frac{\log T}{T})$ in a full-batch function query setting and $\mathcal{O}(\frac{1}{\sqrt{T}})$ in a mini-batch function query setting, respectively. Moreover, we provide a generalization error analysis for the proposed SABO method.
- •
Empirically, we verify the convergence result of the SABO algorithm on the synthetic numerical problems. Moreover, extensive experiments on black-box prompt fine-tuning tasks demonstrate the effectiveness of the proposed SABO method in improving the model generalization performance.

#### Notation and Symbols.

We denote by $\|\cdot\|_{2}$ and $\|\cdot\|_{\infty}$ the $l_{2}$ norm and the $l_{\infty}$ norm for vectors, respectively. $\|\cdot\|_{\mathrm{F}}$ denotes the Frobenius norm for matrices. $\mathcal{S}^{+}$ denotes the set of positive semi-definite matrices. For a square matrix $\mathbf{X}$, $\mathrm{diag}(\mathbf{X})$ represents a vector with diagonal entries in $\mathbf{X}$, and if $\mathbf{x}$ is a vector, $\mathrm{diag}(\mathbf{x})$ represents a diagonal matrix with $\mathbf{x}$ as its diagonal entries. We define $\|X\|_{Y}:=\sqrt{\langle X,YX\rangle}$ for a positive semi-definite
matrix $Y\in\mathcal{S}^{+}$ or a non-negative vector $Y$, where $\langle\cdot,\cdot\rangle$ denotes the inner product under the Frobenius norm for matrices and inner product under the $l_{2}$ norm for vectors. $\frac{X}{Y}$ denotes the elementwise division operation when $X$ and $Y$ are vectors (for the diagonal matrix), and the elementwise division operation for diagonal elements in $X$ and $Y$ when they are diagonal matrices.

## 2 Background

#### Stochastic Gradient Approximation

The stochastic gradient approximation method [, , ] is a representative strategy for solving black-box optimization problems, which instead of maintaining a population of searching points, iteratively updates a search distribution by stochastic gradient approximation. The general procedure of stochastic gradient approximation methods is to first generate a batch of sample points by a parameterized search distribution. Then the sample points allow the algorithm to capture the local structure of the fitness function and appropriately estimate the stochastic gradient to update the distribution.

Specifically, the stochastic gradient approximation method reparameterizes $F(\bm{x})$ as

$$ $\displaystyle J(\bm{\theta})=\mathbb{E}_{p_{\bm{\theta}}(\bm{x})}[F(\bm{x})]=\int F(\bm{x})p(\bm{x};\bm{\theta})\mathrm{d}\bm{x},$ (2) $$

where $\bm{\theta}$ denotes the parameters of density $p(\bm{x};\bm{\theta})$ or $p_{\bm{\theta}}$ and $F(\bm{x})$ is also referred to as the fitness function for $\bm{x}$.
Based on this definition, we can obtain the Monte Carlo estimation of the search gradient as

$$ $\displaystyle\bar{\nabla}_{\bm{\theta}}J(\bm{\theta})=\frac{1}{N}\sum\nolimits_{j=1}^{N}F(\bm{x}_{j})\nabla_{\bm{\theta}}\log p_{\bm{\theta}}(\bm{x}_{j}),$ (3) $$

where $\bm{x}_{j}$ denotes the $j$-th sample and $N$ denotes the number of samples. Therefore, the stochastic gradient $\bar{\nabla}_{\bm{\theta}}J(\bm{\theta})$ provides a search direction in the space of search distributions.

#### Sharpness-Aware Minimization

SAM [] attempts to improve generalization by finding flat minima. This is achieved by minimizing the worst-case loss within some perturbation radius. Mathematically, it is formulated as the following minimax optimization problem:

$$ $\displaystyle\min_{\bm{x}\in\mathcal{X}}\max_{\|\bm{\epsilon}\|_{2}\leq\rho^{2}}F(\bm{x}+\bm{\epsilon}),$ (4) $$

where $F:\mathbb{R}^{d}\to\mathbb{R}$ is the objective function, $\bm{x}$ denotes variables that can represent model parameters, $\rho>0$ is a positive constant, and $\bm{\epsilon}$ is the perturbation whose magnitude is bounded by $\rho^{2}$. By taking the first-order approximation of $F(\bm{x}+\bm{\epsilon})$ over $F(\bm{x})$, a solution of $\bm{\epsilon}$ for the maximization subproblem can be obtained as

$$ $\displaystyle\bm{\epsilon}(\bm{x})=\frac{\rho^{2}\nabla F(\bm{x})}{\|\nabla F(\bm{x})\|_{2}}.$ (5) $$

Then problem ([4](#S2.E4)) can be solved by performing the gradient descent method for the minimization subproblem as $\bm{x}_{t+1}=\bm{x}_{t}-\beta_{t}\nabla_{\bm{x}}F(\bm{x}_{t}+\bm{\epsilon}(\bm{x}_{t})),$ where $\beta_{t}$ represents the step size in the $t$-th iteration. Note that this gradient implicitly depends on the Hessian of $F(\bm{x})$ because $\bm{\epsilon}(\bm{x})$ is a function of $\bm{x}$. To accelerate the computation, a common approach in SAM-based methods [, , ] is to apply a first-order gradient approximation, so we obtain the update rule for SAM as

$$ $\displaystyle\bm{x}_{t+1}=\bm{x}_{t}-\beta_{t}\nabla_{\bm{x}}F(\bm{x}_{t}+\bm{\epsilon}_{t}),$ (6) $$

where $\bm{\epsilon}_{t}=\bm{\epsilon}(\bm{x}_{t})$ is viewed as a constant w.r.t. $\bm{x}$.

## 3 Methodology

In this section, we introduce the proposed SABO algorithm. Firstly, we formulate the sharpness-aware black-box optimization as a min-max optimization problem and solve it in Section [3.1](#S3.SS1), and in Section [3.2](#S3.SS2), we derive the update formula of parameters in the search distribution. The detailed derivations in this section are put in Appendix [A](#A1).

### 3.1 Sharpness-Aware Black-box Optimization

Suppose we are given a training set $\mathcal{D}$ with i.i.d. samples $\{(X_{i},y_{i})\}$. The main objective is defined as

$$ $F(\bm{x};\mathcal{D})=\frac{1}{|\mathcal{D}|}\sum\nolimits_{(X_{i},y_{i})\in\mathcal{D}}l(\bm{x};(X_{i},y_{i})),$ (7) $$

where $\bm{x}$ denotes model parameters, $|\mathcal{D}|$ denotes the number of data in the dataset $\mathcal{D}$, and $l(\bm{x};(X,y))$ denotes the loss function (e.g., the cross-entropy loss for classification). To simplify the notation, we define $F(\bm{x}):=F(\bm{x};\mathcal{D})$.

In black-box optimization, we aim at minimizing the objective function $F(\bm{x}),$ with only function queries. Due to the lack of gradient information, we first apply the stochastic gradient approximation method [, ]. We denote by $\bm{\theta}$ the parameters of the search distribution $p_{\bm{\theta}}$ and define the expected fitness of $F(\bm{x})$ under the parametric search distribution $p_{\bm{\theta}}(\bm{x})$ as $J(\bm{\theta})=\mathbb{E}_{p_{\bm{\theta}}(\bm{x})}[F(\bm{x})]$. Then the optimal parameter $\bm{\theta}$ can be found by minimizing the reparameterized objective $J(\bm{\theta})$.

Inspired by SAM [], we attempt to improve generalization by finding flat minima of $\bm{\theta}$. However, for the reparameterized objective $J(\bm{\theta})$, the geometry of the corresponding distribution space is not Euclidean but a statistical manifold, where the distance between two probability distributions is defined by some statistical distance, e.g., Kullback-Leibler (KL) divergence. Therefore, instead of restricting the perturbation in an Euclidean ball, we restrict the perturbation distribution to be inside a small neighborhood of the unperturbed distribution w.r.t. the KL divergence []. The proposed optimization problem for black-box optimization is formulated as

$$ $\min_{\bm{\theta}}\max_{\bm{\delta}\in\mathcal{C}(\bm{\theta})}J(\bm{\theta}+\bm{\delta}),$ (8) $$

where $\mathcal{C}(\bm{\theta})=\{\bm{\delta}\mid\mathrm{KL}(p_{\bm{\theta}+\bm{\delta}}\|p_{\bm{\theta}})\leq\rho^{2}\}$, $\rho$ is a positive constant, and $J(\bm{\theta}+\bm{\delta})=\mathbb{E}_{x\sim p_{\bm{\theta}+\bm{\delta}}}[F(\bm{x})]$. Note that $\mathcal{C}(\bm{\theta})$ defines the neighborhood around a given distribution $p_{\bm{\theta}}$ in the distribution space, which is different from the neighborhood of SAM that is defined in the parameter space. Here $\rho$ defines the size of the neighborhood.

Generally, problem ([8](#S3.E8)) is applicable to any family of distribution $p_{\bm{\theta}}$. For computational consideration, the search distribution is assumed to be a Gaussian distribution, i.e., $p_{\bm{\theta}}(\bm{x})=\mathcal{N}(\bm{x}\mid\bm{\mu},\bm{\Sigma})$ where $\bm{\mu}$ denotes the mean and $\bm{\Sigma}$ denotes the covariance matrix, and correspondingly $\bm{\theta}$ includes $\bm{\mu}$ and $\bm{\Sigma}$, i.e., $\bm{\theta}=\{\bm{\mu},\bm{\Sigma}\}$. In this work, we assume that the covariance matrix $\bm{\Sigma}$ is a diagonal matrix. For a perturbation $\bm{\delta}=\{\bm{\delta}_{\bm{\mu}},\bm{\delta}_{\bm{\Sigma}}\}$ where $\bm{\delta}_{\bm{\Sigma}}$ is a digaonal matrix, the perturbed distribution is a Gaussian distribution $p_{\bm{\theta}+\bm{\delta}}(\bm{x})=\mathcal{N}(\bm{x}\mid\bm{\mu}+\bm{\delta}_{\bm{\mu}},\bm{\Sigma}+\bm{\delta}_{\bm{\Sigma}})$.

We need to solve problem ([8](#S3.E8)) to derive the update formulation for $\bm{\theta}$. By using the first-order Taylor expansion, the maximization subproblem can be approximated as a quadratically constrained linear programming problem:

$$ $\max_{\bm{\delta}\in\mathcal{C}(\bm{\theta})}J(\bm{\theta}+\bm{\delta})\approx\max_{\bm{\delta}\in\mathcal{C}(\bm{\theta})}\left<\nabla_{\bm{\theta}}J(\bm{\theta}),\bm{\delta}\right>,$ (9) $$

The corresponding Lagrangian of problem ([9](#S3.E9)) is

$$ $\mathcal{L}(\bm{\delta},\lambda)=-\left<\nabla_{\bm{\theta}}J(\bm{\theta}),\bm{\delta}\right>+\lambda(\mathrm{KL}(p_{\bm{\theta}+\bm{\delta}}\|p_{\bm{\theta}})-\rho^{2}),$ (10) $$

where $\lambda$ is the Lagrange multiplier. We can see that problem ([10](#S3.E10)) is convex with respect to $\bm{\delta}$. Therefore, by setting the derivatives w.r.t. $\bm{\delta}_{\bm{\mu}}$ and $\bm{\delta}_{\bm{\Sigma}}$ to zero, we can obtain $\bm{\delta}$ as

$$ $\displaystyle\bm{\delta}_{\bm{\mu}}(\bm{\theta})=\frac{1}{\lambda}\bm{\Sigma}\nabla_{\bm{\mu}}J(\bm{\theta}),\ ~{}\bm{\delta}_{\bm{\Sigma}}(\bm{\theta})=\frac{2\bm{\Sigma}\nabla_{\bm{\Sigma}}J(\bm{\theta})}{\lambda\bm{\Sigma}^{-1}-2\nabla_{\bm{\Sigma}}J(\bm{\theta})}.$ (11) $$

As shown in Eq. ([11](#S3.E11)), to calculate the perturbation, we need to calculate the inverse covariance matrix, i.e., $\bm{\Sigma}^{-1}$, which is computationally expensive for high-dimensional problems. Therefore, assuming that $\bm{\Sigma}$ is a diagonal matrix can significantly reduce the computation cost.
Then by plugging $\bm{\delta}_{\bm{\mu}}(\bm{\theta})$ and $\bm{\delta}_{\bm{\Sigma}}(\bm{\theta})$ into the neighborhood constraint, i.e., $\bm{\delta}\in\mathcal{C}(\bm{\theta})$, we can determine the optimal $\lambda$ as

$$ $\displaystyle\lambda=\frac{1}{\rho}\sqrt{\|\bm{\Sigma}\nabla_{\bm{\Sigma}}J(\bm{\theta})\|_{\mathrm{F}}^{2}+0.5\|\bm{\Sigma}^{\frac{1}{2}}\nabla_{\bm{\mu}}J(\bm{\theta})\|^{2}_{2}}.$ (12) $$

Based on Eqs. ([11](#S3.E11)) and ([12](#S3.E12)), we obtain the approximated closed-form solution $\bm{\delta}(\bm{\theta})$ for a given $\bm{\theta}$. With $\bm{\delta}(\bm{\theta})$, the minimization subproblem of problem ([8](#S3.E8)) can be reformulated as

$$ $\displaystyle\min_{\bm{\theta}}J(\bm{\theta}+\bm{\delta}(\bm{\theta})).$ (13) $$

To solve problem ([13](#S3.E13)), in the $t$-th iteration, we add a regularization term $\frac{1}{\beta_{t}}\mathrm{KL}(p_{\bm{\theta}}\|p_{\bm{\theta}_{t}})$ to problem ([13](#S3.E13)) to enforce $\bm{\theta}$ to be close to $\bm{\theta}_{t}$, and obtain $\bm{\theta}_{t+1}$ by solving the following problem as

$$ $\bm{\theta}_{t+1}=\arg\min_{\bm{\theta}}J(\bm{\theta}+\bm{\delta}(\bm{\theta_{t}}))-J(\bm{\theta}_{t})+\frac{1}{\beta_{t}}\mathrm{KL}(p_{\bm{\theta}}\|p_{\bm{\theta}_{t}}).$ (14) $$

Following standard SAM-based methods [], we treat $\bm{\delta}(\bm{\theta}_{t})$ as a constant $\bm{\delta}_{t}$ instead of a function of $\bm{\theta}_{t}$ to accelerate the computation. Then by using the first-order Taylor expansion, problem ([14](#S3.E14)) can be approximated as

$$ $\displaystyle\bm{\theta}_{t+1}$ $\displaystyle=\arg\min_{\bm{\theta}}\langle\bm{\theta}-\bm{\theta}_{t},\nabla_{\bm{\theta}}J(\bm{\theta}_{t}+\bm{\delta}_{t})\rangle+\frac{1}{\beta_{t}}\mathrm{KL}(p_{\bm{\theta}}\|p_{\bm{\theta}_{t}}).$ (15) $$

By solving problem ([15](#S3.E15)), we can obtain the update formulations for $\bm{\mu}$ and $\bm{\Sigma}$ in the $t$-th iteration as

$$ $\displaystyle\bm{\mu}_{t+1}=\bm{\mu}_{t}-\beta_{t}\bm{\Sigma}_{t}\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t}),\ ~{}\bm{\Sigma}_{t+1}^{-1}=\bm{\Sigma}_{t}^{-1}+2\beta_{t}\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t}+\bm{\delta}_{t}),$ (16) $$

where $\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$ and $\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$ denote the derivative of $J(\bm{\theta})$ w.r.t. $\bm{\mu}$ and $\bm{\Sigma}$ at $\bm{\mu}=\bm{\theta}_{t}+\bm{\delta}_{\bm{\mu}_{t}}$ and $\bm{\Sigma}=\bm{\Sigma}_{t}+\bm{\delta}_{\bm{\Sigma}_{t}}$, respectively.

### 3.2 Update Formulations for SABO

The gradients of the reparameterized objective $J(\bm{\theta})$ w.r.t. $\bm{\mu}$ and $\bm{\Sigma}$ rely on the expectations of the black-box function and can be obtained with only function queries [] (see Theorem [A.1](#A1.Thmtheorem1) in Appendix [A.4](#A1.SS4)).
Hence, we estimate them by Monte Carlo sampling.
Specifically, the stochastic approximation of the gradients $\nabla_{\bm{\mu}}J(\bm{\theta}_{t})$ and $\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})$ are given as

$$ $\displaystyle\bm{g}_{t}^{\prime}$ $\displaystyle=\frac{1}{N}\sum\nolimits_{j=1}^{N}\bm{\Sigma}_{t}^{-1}(\bm{x}_{j}^{\prime}-\bm{\mu}_{t})\big{(}F(\bm{x}_{j}^{\prime})-F(\bm{\mu}_{t})\big{)},$ (17) $\displaystyle\bm{G}_{t}^{\prime}$ $\displaystyle=\frac{1}{2N}\sum\nolimits_{j=1}^{N}\mathrm{diag}\Big{[}\bm{\Sigma}_{t}^{-1}\big{[}\mathrm{diag}\big{(}(\bm{x}_{j}^{\prime}-\bm{\mu}_{t})(\bm{x}_{j}^{\prime}-\bm{\mu}_{t})^{\top}\bm{\Sigma}_{t}^{-1}-\bm{I}\big{)}(F(\bm{x}_{j}^{\prime})-F(\bm{\mu}_{t}))\big{]}\Big{]},$ (18) $$

where $\bm{x}_{j}^{\prime}$ denotes the $j$-th sample sampled from the distribution $\mathcal{N}(\bm{x}\mid\bm{\mu}_{t},\bm{\Sigma}_{t})$. Note that $\bm{g}_{t}^{\prime}$ is an unbiased estimator for the gradient $\nabla_{\bm{\mu}}J(\bm{\theta})$ as proved in Lemma [C.5](#A3.Thmtheorem5), and inspired by [], we subtract $F(\bm{\mu}_{t})$ to improve the computational stability while keeping them as unbiased estimations.

Then according to Eq. ([11](#S3.E11)), we obtain the perturbation $\bm{\delta}_{t}$ in the $t$-th iteration as

$$ $\displaystyle\bm{\delta}_{t}=\left\{\frac{1}{\lambda}\bm{\Sigma}_{t}\bm{g}_{t}^{\prime},\frac{2\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}}{\lambda\bm{\Sigma}_{t}^{-1}-2\bm{G}_{t}^{\prime}}\right\},$ (19) $$

where $\lambda$ is approximated by $\frac{1}{\rho}\sqrt{\|\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}\|_{\mathrm{F}}^{2}+0.5\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}^{\prime}\|^{2}_{2}}$.
Similarly, the gradients $\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$ and $\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$ can be approximated as follows:

$$ $\displaystyle\bm{g}_{t}$ $\displaystyle=\frac{1}{N}\sum\nolimits_{j=1}^{N}\widehat{\bm{\Sigma}}_{t}^{-1}(\bm{x}_{j}-\widehat{\bm{\mu}}_{t})\big{(}F(\bm{x}_{j})-F(\widehat{\bm{\mu}}_{t})\big{)},$ (20) $\displaystyle\bm{G}_{t}$ $\displaystyle=\frac{1}{2N}\sum\nolimits_{j=1}^{N}\mathrm{diag}\Big{[}\widehat{\bm{\Sigma}}_{t}^{-1}\big{[}\mathrm{diag}\big{(}(\bm{x}_{j}-\widehat{\bm{\mu}}_{t})(\bm{x}_{j}-\widehat{\bm{\mu}}_{t})^{\top}\widehat{\bm{\Sigma}}_{t}^{-1}-\bm{I}\big{)}(F(\bm{x}_{j})-F(\widehat{\bm{\mu}}_{t}))\big{]}\Big{]},$ (21) $$

where $\widehat{\bm{\Sigma}}_{t}=\bm{\Sigma}_{t}+\bm{\delta}_{\bm{\Sigma}_{t}}$, $\widehat{\bm{\mu}}_{t}=\bm{\mu}_{t}+\bm{\delta}_{\bm{\mu}_{t}}$, and $\bm{x}_{j}$ denotes the $j$-th sample sampled from the distribution $\mathcal{N}(\bm{x}\mid\widehat{\bm{\mu}}_{t},\widehat{\bm{\Sigma}}_{t})$.
Then the updated formulations for $\bm{\mu}$ and $\bm{\Sigma}$ are rewritten as

$$ $\displaystyle\bm{\mu}_{t+1}=\bm{\mu}_{t}-\beta_{t}\bm{\Sigma}_{t}\bm{g}_{t},\ ~{}~{}~{}\bm{\Sigma}_{t+1}^{-1}=\bm{\Sigma}_{t}^{-1}+2\beta_{t}\bm{G}_{t}.$ (22) $$

The entire algorithm is shown in Algorithm [1](#alg1). For the proposed SABO method, the computation cost per iteration is of order $\mathcal{O}(Nd)$.

#### Mini-batch SABO.

In Algorithm [1](#alg1), we assume full access to the objective function $F(\bm{x};\mathcal{D})$, while in practice, a full-batch function query might be costly. Therefore, we can perform a mini-batch function query. Specifically, in each iteration, we query the expected fitness by a mini-batch of data and approximate $F(\bm{x};\mathcal{D})$ by

$$ $F(\bm{x};\mathcal{B})=\frac{1}{|\mathcal{B}|}\sum\nolimits_{(X,y)\in\mathcal{B}}l(\bm{x};(X,y)),$ (23) $$

where $|\mathcal{B}|$ denotes the number of data in the mini-batch $\mathcal{B}$. The corresponding expected fitness of $F(\bm{x};\mathcal{B})$ under the distribution $p_{\bm{\theta}}(\bm{x})$ is formulated as $J(\bm{\theta};\mathcal{B})=\mathbb{E}_{p_{\bm{\theta}}(\bm{x})}[F(\bm{x};\mathcal{B})]$. Then similar to the full-batch function query setting, we can approximate the gradients $\nabla_{\bm{\mu}}J(\bm{\theta}_{t},\mathcal{B})$ and $\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t},\mathcal{B})$, where the detailed formulations can be found in Appendix [B](#A2). The entire SABO algorithm with mini-batch function queries is shown in Algorithm [2](#alg2) in Appendix [B](#A2).

Figure: Algorithm 1 SABO

## 4 Analysis

In this section, we provide comprehensive theoretical analyses for the proposed SABO method with all the detailed proofs in Appendix [D](#A4).

### 4.1 Convergence Analysis of SABO

Firstly, we make an assumption for the reparameterized objective function.

###### Assumption 4.1 .

The function $J(\bm{\theta})$ is $H$-Lipschitz and $L$-smoothness w.r.t. $\bm{\theta}=\{\bm{\mu},\bm{\Sigma}\}\in\Theta$, where $\Theta:=\{\bm{\mu},\bm{\Sigma}\mid\bm{\mu}\in\mathbb{R}^{d},\bm{\Sigma}\in\mathcal{S}^{+}\}$.

The smoothness assumption in Assumption [4.1](#S4.Thmtheorem1) has been widely adopted in optimization literature []. Note that the proposed SABO algorithm approximates the gradients of the reparameterized objective function.
It is necessary to study the relation between the optimal solutions of the original objective functions $F(\bm{x})$ and $J(\bm{\theta})$, and we put the results in the following proposition.

###### Proposition 4.2 .

[]
Suppose $p_{\bm{\theta}}(\bm{x})$ is a Gaussian distribution with $\bm{\theta}=\{\bm{\mu},{\bm{\Sigma}}\}$ and $F(\bm{x})$ is a convex function. Let $J(\bm{\theta})=\mathbb{E}_{p_{\bm{\theta}}}[F(\bm{x})]$, and $J(\bm{\mu}^{*},\bm{0}):=F(\bm{\mu}^{*})$. Then we have

$$ $\displaystyle F(\bm{\mu})-F(\bm{\mu}^{*})\leq J(\bm{\mu},\bm{\Sigma})-J(\bm{\mu}^{*},\bm{0}),$ (24) $$

where $\bm{0}$ denotes a zero matrix with appropriate size.

The convexity assumption of the objective function in Theorem [4.3](#S4.Thmtheorem3) has been widely adopted in the area of stochastic gradient approximation black-box optimization [, , , ].
Since the Gaussian-smooth approximation function is always an upper bound of the true target function in convex cases, i.e., $F(\mu)\leq\mathbb{E}_{\mathcal{N}(\mu,\Sigma)}[F(x)]$.
When $\bm{\mu}^{*}$ is an optimal solution of minimization problem $\min_{\bm{x}}F(\bm{x})$, Proposition [4.2](#S4.Thmtheorem2) implies that the difference between the objective value at $\bm{\mu}$ and the optimal objective value of the original problem is upper-bounded by that of the expected objective function. Then for Algorithm [1](#alg1), the following theorem captures the convergence of $\bm{\mu}$ for a convex objective function.

###### Theorem 4.3 .

Suppose that $F(\bm{x})$ is a convex function, $J(\bm{\theta})$ is $c$-strongly convex w.r.t. $\bm{\mu}$, the gradient estimator $\bm{G}_{t}$ (w.r.t. the covariance matrix) is positive semi-definite matrix such that $\xi\bm{I}\preceq\bm{G}_{t}\preceq\frac{c\bm{I}}{4}$ with $\xi\geq 0$, $\bm{\Sigma}_{0}\in\mathcal{S}^{+}$, and $\bm{\Sigma}_{0}\preceq R\bm{I}$ where $R>0$. Suppose the sequence $\{\bm{\mu}_{t}\}$ generated by Algorithm [1](#alg1) satisfies that the distance between the sequence $\{\bm{\mu}_{t}\}$ and the optimal solution of $F(\bm{x})$ is bounded, i.e., $\|\bm{\mu}_{t}-\bm{\mu}^{*}\|\leq D$, $\beta_{t}=\mathcal{O}(1)$, and $\rho<\frac{\sqrt{d}}{2}$ satisfies $\rho=\mathcal{O}(\frac{1}{\sqrt{T}})$, then with Assumption [4.1](#S4.Thmtheorem1), we have

$$ $\frac{1}{T}\sum\nolimits_{t=0}^{T-1}\mathbb{E}\left[J(\bm{\mu}_{t+1},\bm{\Sigma}_{t})-J(\bm{\mu}^{*},0)\right]=\mathcal{O}\left(\frac{\log T}{T}\right).$ (25) $$

Based on Theorem [4.3](#S4.Thmtheorem3) and Proposition [4.2](#S4.Thmtheorem2), when $\beta_{t}=\mathcal{O}(1)$ and $\rho=\mathcal{O}(\frac{1}{\sqrt{T}})$, we have

$$ $\displaystyle\frac{1}{T}\sum\nolimits_{t=0}^{T-1}\mathbb{E}\left[F(\bm{\mu}_{t+1})-F(\bm{\mu}^{*})\right]=\mathcal{O}(\frac{\log T}{T}).$ (26) $$

Therefore, the proposed SABO algorithm with full-batch function query possesses a convergence rate $\mathcal{O}(\frac{\log T}{T})$ for the convex objective function. Additionally, in Theorem [4.3](#S4.Thmtheorem3), when $\beta=\mathcal{O}(1)$ and $\rho=\mathcal{O}(1)$, the proposed SABO algorithm still maintains a convergence rate $\mathcal{O}(T^{-\frac{1}{2}})$. The detailed discussion is provided in Remark [D.1](#A4.Thmtheorem1).

For the mini-batch setting, we make an additional assumption for the objective function $F(\bm{x};\mathcal{D})$.

###### Assumption 4.4 .

It is assumed that the datasets $\mathcal{D}$ and $\mathcal{B}$ are i.i.d. sampled from a data distribution $P(X,y)$, and the variance of the mini-batch estimation of the objective function is bounded, i.e., $\|F(\bm{x};\mathcal{B})-\mathbb{E}F(\bm{x};\mathcal{B})\|_{2}^{2}\leq\varepsilon_{\mathrm{B}}^{2}$ and $\|F(\bm{x};\mathcal{D})-\mathbb{E}F(\bm{x};\mathcal{D})\|_{2}^{2}\leq\varepsilon_{\mathrm{D}}^{2}$.

Note that for the standard stochastic gradient descent (SGD) method [], the unbiased estimation and bounded variance assumptions were made for the approximated gradient. However, in black-box optimization, the gradient of the objective function $F(\bm{x})$ w.r.t. $\bm{x}$ is unavailable. Hence we can only make assumptions for the batch estimations $F(\bm{x};\mathcal{B})$ and $F(\bm{x};\mathcal{D})$.

Then we have the following result for the mini-batch estimation.

###### Proposition 4.5 .

Suppose Assumption [4.4](#S4.Thmtheorem4) holds, then we have $\mathbb{E}F(\bm{x};\mathcal{B})=\mathbb{E}F(\bm{x};\mathcal{D})$, and $\|F(\bm{x};\mathcal{B})-F(\bm{x};\mathcal{D})\|_{2}^{2}\leq\varepsilon^{2}$, where $\varepsilon^{2}=2(\varepsilon_{\mathrm{B}}^{2}+\varepsilon_{\mathrm{D}}^{2})$.

Then with Assumption [4.4](#S4.Thmtheorem4), the following theorem shows the convergence of $\bm{\mu}$ for Algorithm [2](#alg2).

###### Theorem 4.6 .

Suppose that $F(\bm{x})$ is a convex function, $J(\bm{\theta})$ is $c$-strongly convex w.r.t. $\bm{\mu}$, the gradient estimator $\bm{G}_{t}$ (w.r.t. the covariance matrix) is positive semi-definite matrix such that $\xi\bm{I}\preceq\bm{G}_{t}\preceq\frac{c\bm{I}}{4}$ with $\xi\geq 0$, $\bm{\Sigma}_{0}\in\mathcal{S}^{+}$, and $\bm{\Sigma}_{0}\preceq R\bm{I}$ where $R>0$. Suppose the sequence $\{\bm{\mu}_{t}\}$ generated by Algorithm [2](#alg2) satisfies that the distance between the sequence $\{\bm{\mu}_{t}\}$ and the optimal solution of $F(\bm{x})$ is bounded, i.e., $\|\bm{\mu}_{t}-\bm{\mu}^{*}\|\leq D$, $\beta_{t}=\mathcal{O}(1)$, and $\rho<\frac{\sqrt{d}}{2}$ satisfies $\rho=\mathcal{O}(1)$. Then with Assumptions [4.1](#S4.Thmtheorem1) and [4.4](#S4.Thmtheorem4), we have

$$ $\begin{split}\frac{1}{T}\sum\nolimits_{t=0}^{T-1}\mathbb{E}\left[J(\bm{\mu}_{t+1},\bm{\Sigma}_{t})-J(\bm{\mu}^{*},0)\right]=\mathcal{O}\left(\frac{1}{\sqrt{T}}\right).\end{split}$ (27) $$

Based on Theorem [4.6](#S4.Thmtheorem6) and Proposition [4.2](#S4.Thmtheorem2), the proposed SABO algorithm with mini-batch function query possesses a convergence rate $\mathcal{O}(T^{-\frac{1}{2}})$.

###### Remark 4.7 .

Note that in Theorem [4.3](#S4.Thmtheorem3) and Theorem [4.6](#S4.Thmtheorem6), the convergence results do not require the objective function $F(\bm{x})$ to be strongly convex or differentiable. Hence, the convergence holds for a non-smooth convex function $F(\bm{x})$ (as long as $J(\bm{\theta})$ being a $c$-strongly convex function w.r.t. $\bm{\mu}$).

###### Assumption 4.1 .

###### Proposition 4.2 .

###### Theorem 4.3 .

###### Assumption 4.4 .

###### Proposition 4.5 .

###### Theorem 4.6 .

###### Remark 4.7 .

### 4.2 Generalization Error Analysis

In this subsection, we analyze the generalization bound of the proposed SABO algorithm. Specifically, we bound the expectation of the objective function over the Gaussian perturbation.

We denote by $(X,y)$ a data pair drawn from a data distribution $P(X,y)$ and by $F(\bm{x};(X,y))$ the corresponding loss of parameter $\bm{x}$ on $(X,y)$. So we have
$F(\bm{x};(X,y))=l(\bm{x};(X,y))$.
We define the population loss over the data distribution $P(X,y)$ as $\mathbb{E}_{{P}(X,y)}[F(\bm{x};(X,y))]$, and the empirical loss over a dataset $\mathcal{S}$, which consists of $M$ i.i.d. samples drawn from ${P}(X,y)$, as $F(\bm{x};\mathcal{S})=\frac{1}{M}\sum_{i=1}^{M}l(\bm{x};(X_{i},y_{i}))$. Then we have following result.

###### Theorem 4.8 .

Let the loss function $F(\bm{x};(X,y))$ be a convex function w.r.t. $\bm{x}$, then for any $\bm{\mu}\in\mathbb{R}^{d}$, with probability at least $1-\kappa$, we have

$$ $\mathbb{E}_{P(X,y)}[F(\bm{\mu};(X,y))]\leq\max_{\bm{\delta}\in\mathcal{C}(\bm{\theta})}\mathbb{E}_{p_{\bm{\theta}+\bm{\delta}}}\big{[}F(\bm{x};S)\big{]}+\sqrt{\frac{\rho^{2}+\log(\frac{M}{\kappa})}{2(M-1)}},$ (28) $$

where $p_{\bm{\theta}}:=\mathcal{N}(\bm{\mu},\bm{\Sigma})$, $\mathcal{C}(\bm{\theta})=\{\bm{\delta}\mid\mathrm{KL}(p_{\bm{\theta}+\bm{\delta}}\|p_{\bm{\theta}})\leq\rho^{2}\}$, and $\mathcal{S}$ denotes the training set that consists of $M$ i.i.d. samples drawn from data distribution ${P}(X,y)$.

Theorem [4.8](#S4.Thmtheorem8) provides a generalization bound for the proposed SABO algorithm. Compared with the generalization bound of SAM presented in Appendix A.1 of [], Theorem [4.8](#S4.Thmtheorem8) has an asymptotically identical order in the complexity term. However, the expected generalization loss on the right-hand side of Eq. ([28](#S4.E28)) is different in that we reparameterize the objective function and have perturbation of $\bm{\theta}$ in its neighborhood of a statistical manifold, i.e., $\bm{\delta}\in\mathcal{C}(\bm{\theta})$, while SAM bounds the generalization loss averaged over a spherical Gaussian perturbation on parameters.

###### Theorem 4.8 .

## 5 Related Works

Black-Box Optimization.
Many methods have been proposed for black-box optimization, including Bayesian optimization (BO) methods [, , ], stochastic optimization methods such as evolution strategies (ES) [, , ], and genetic algorithms (GA) [, ]. Among those methods, BO achieves good performance for low-dimensional problems, but it often fails to handle high-dimensional problems through a global surrogate model, as shown in [] and []. As a result, TuRBO [] and GIBO [] try to address this problem with the local BO approach. [] further showed that decomposition is important for alleviating the high-dimensional problems in BO. Although BO is not our main focus, we further compare our method with these local Bayesian optimization methods on the black-box prompt fine-tuning problem, and the corresponding exploration is shown in Appendix [F.2](#A6.SS2).
GA method is computationally expensive for machine learning problems and usually lacks convergence analysis.
The stochastic optimization methods such as CMA-ES [] and INGO [] can scale up to higher-dimensional problems compared with BO. Hence we mainly consider stochastic optimization methods as baseline methods in our experiments.

Sharpness-Aware Minimization. SAM has been widely studied for improving the model generalization. Among previous works on SAM [, , , ], the most relevant method to our approach is the FSAM method [] which also finds the worst-case objective function via a statistical manifold instead of the Euclidean space. However, the loss function of the model studied in FSAM is a predictive distribution conditional on both model parameters and data, while in our case, we consider the parameter as a Gaussian distribution. The bSAM method [] builds a connection between the SAM objective and Bayes objective by Fenchel biconjugate of the loss function. [] shares a similarity with our work in developing SAM w.r.t. the expected loss. However, The bSAM method relies on the derivation of a convex lower bound of the expected loss by the Fenchel biconjugate and the perturbation is still w.r.t. each point inside the expected loss as standard SAM. Hence FSAM and bSAM are different from the proposed SABO method. Additionally, like other variants of SAM, FSAM and bSAM are both inapplicable to black-box optimization. The STABLEOPT method [] proposes a SAM-like optimization formulation in the Bayesian optimization area. They aim to improve the robustness w.r.t. the adversarial perturbation of the return point by a GP-based optimization. Their method relies on a GP surrogate model that is expensive for training and inference. In addition, it is challenging for the proposed adversarial robust GP-based optimization to handle high-dimensional problems. In contrast, our work aims to improve the generalization property in high-dimensional black-box optimization.

## 6 Empirical Study

In this section, we empirically evaluate the proposed SABO method, and compare it with four representative black-box methods, i.e., CMA-ES [], MMES [], BES [], and INGO [].
All the experiments are conducted on a single NVIDIA GeForce RTX 3090 GPU.

### 6.1 Synthetic Problems

To verify the convergent results of the proposed SABO method in Section [4](#S4). We compare the proposed SABO method with baseline methods on minimizing four $d$-dimensional synthetic benchmark test functions, i.e., ellipsoid function, $l_{\frac{1}{2}}$-ellipsoid function, different powers function, and Levy function. All the test functions are listed in Appendix [F.1](#A6.SS1).

The results are evaluated by calculating the Euclidean distance between the solution $\bm{x}$ and the optimal solution $\bm{x}^{*}$, i.e., $\mathcal{E}=\|\bm{x}-\bm{x}^{*}\|_{2}$. We then assess the baseline methods using varying dimensions, i.e., $d\in\{200,500,1000\}$. Due to the page limitation, the implementation details and more detailed experimental results are put in Appendix [F.1](#A6.SS1).

Figure: (a) Ellipsoid.
Refer to caption: /html/2410.12457/assets/x1.png

#### Result.

Figure [1](#S6.F1) shows the results on four $d$-dimensional synthetic problems with $d=500$ and population size $N=50$. The proposed SABO method approximately achieves a linear convergence rate similar to the INGO method. This is reasonable as these two methods have the same theoretical convergent rate, i.e., $\mathcal{O}(\frac{\log T}{T})$, according to Theorem [4.3](#S4.Thmtheorem3) and Theorem 5 in []. Since the SABO method perturbs the main objective in each iteration, its practical convergence speed is slightly slower than INGO. The CMA-ES method and MMES methods can converge on the ellipsoid problem and the different powers problem, but they do not converge as fast as SABO and INGO methods. Moreover, CMA-ES and MMES methods fail on the $l_{\frac{1}{2}}$-ellipsoid problem and Levy problem. The BES method fails on all test problems. This shows that it could be challenging for BES to optimize non-smooth or high-dimensional test functions without adaptively updating mean and covariance. These results demonstrate the superiority of the SABO method in optimizing high-dimensional problems, and verify our theoretical convergence results.

### 6.2 Black-box Prompt Fine-tuning

Black-box prompt fine-tuning of large language models [, , , ] is a promising direction to achieve expertise models efficiently for downstream tasks. In such an LMaaS setting, we cannot access the model parameter and can only tune their prompts without backpropagation.
We evaluate the proposed SABO method in improving generalization performance on the black-box prompt fine-tuning task.

#### Datasets.

We conduct experiments on six language understanding benchmark datasets: SST-2 [] and Yelp polarity [] for sentiment analysis, AG’s News [] for topic classification, MRPC [] for paraphrase, RTE [] and SNLI [] for natural language inference. Each dataset contains a classification task.
The statistics of six datasets are summarized in Table 1 of []. By following [], the testing accuracy is used to measure the performance of all the methods on the SST-2, AG’s News, RTE, SNLI, and Yelp P. datasets, and the F1 score is used to measure the performance on the MRPC datasets.

#### Implementation Details.

Following [], we employ a fixed randomly initialized matrix $\bm{A}\in\mathbb{R}^{d\times D}$ to project a vector $\bm{v}\in\mathbb{R}^{d}$ onto the token embedding space $\mathbb{R}^{D}$. Then we optimize the vector $\bm{v}\in\mathbb{R}^{d}$ instead of directly optimizing the prompt $\bm{p}\in\mathbb{R}^{D}$.
The pre-trained $\text{RoBERTa}_{\text{LARGE}}$ model [] is used as the backbone model. The matrix $\bm{A}$ is sampled from the normal distribution as described in [], i.e., $\mathcal{N}(0,\frac{\sigma_{e}}{\sqrt{d}})$, where $\sigma_{e}$ is the standard deviation of word embeddings in $\text{RoBERTa}_{\text{LARGE}}$.
The templates and label words in Table 1 of [] are used to conduct the zero-shot baseline.

For CMA-ES, MMES, BES, INGO, and SABO methods, we employ the cross-entropy loss of training data as the black-box objective for six datasets and optimize the vector $\bm{v}$ with $100$ iterations. The Gaussian distributions are initialized as $\bm{\mu}_{0}=\bm{0}$ and $\bm{\Sigma}_{0}=\bm{I}$, and the population size $N$ is set to $100$.
We perform a grid search for hyperparameters of INGO, SABO, and BES methods. Specifically, we search the learning rate $\beta$ over $\{0.1,0.5,1,5\}$ for INGO, SABO, and BES, the neighborhood size $\rho$ over $\{10,50,100,500\}$ for SABO, and the spacing $c$ over $\{0.1,1,10\}$ for BES.
Additionally, we evaluate the performance of all methods on different dimensions of $\bm{v}$, specifically $d\in\{200,500,1000\}$. All the experiments are
performed in three independent runs, and the experimental results of mean objective $\pm$ std are reported.

**Table 1: Performance ($\%$) on SST-2, AG’s News, MRPC, RTE, SNLI and Yelp P. datasets. We report the mean and standard deviation over $3$ random seeds. The best result across all groups is highlighted in bold and the best result in each group is marked with underlined.**
| Methods | SST-2 | AG’s News | MRPC | RTE | SNLI | Yelp P. |
| --- | --- | --- | --- | --- | --- | --- |
| Zero-shot | $79.82$ | $76.96$ | $67.40$ | $51.62$ | $38.82$ | $89.64$ |
| Dimension $d=200$ |  |  |  |  |  |  |
| CMA-ES | $85.74_{\pm 0.35}$ | $82.09_{\pm 0.56}$ | $74.98_{\pm 2.16}$ | $51.02_{\pm 2.14}$ | $34.27_{\pm 1.18}$ | $90.57_{\pm 0.05}$ |
| MMES | $83.98_{\pm 0.78}$ | $80.52_{\pm 0.99}$ | $76.54_{\pm 4.34}$ | $48.50_{\pm 0.45}$ | $40.39_{\pm 1.83}$ | $90.94_{\pm 0.36}$ |
| BES | $83.52_{\pm 0.11}$ | $75.44_{\pm 0.31}$ | $79.23_{\pm 0.20}$ | $53.07_{\pm 0.29}$ | $38.73_{\pm 0.17}$ | $89.65_{\pm 0.01}$ |
| INGO | $83.57_{\pm 0.11}$ | $76.47_{\pm 0.03}$ | $78.87_{\pm 0.20}$ | $53.07_{\pm 0.00}$ | $38.86_{\pm 0.06}$ | $89.84_{\pm 0.04}$ |
| SABO | $\underline{87.88}_{\pm 0.53}$ | $\underline{82.22}_{\pm 0.41}$ | $\underline{79.35}_{\pm 0.12}$ | $\underline{\textbf{53.67}}_{\pm 0.17}$ | $\underline{40.72}_{\pm 0.15}$ | $\underline{91.50}_{\pm 0.13}$ |
| Dimension $d=500$ |  |  |  |  |  |  |
| CMA-ES | $86.12_{\pm 0.59}$ | $82.50_{\pm 0.23}$ | $77.10_{\pm 1.90}$ | $52.71_{\pm 0.51}$ | $41.34_{\pm 1.49}$ | $91.19_{\pm 0.44}$ |
| MMES | $85.28_{\pm 0.94}$ | $81.67_{\pm 0.80}$ | $77.31_{\pm 1.24}$ | $48.74_{\pm 0.59}$ | $42.07_{\pm 2.62}$ | $91.39_{\pm 0.24}$ |
| BES | $83.56_{\pm 0.05}$ | $75.93_{\pm 0.17}$ | $79.21_{\pm 0.09}$ | $52.95_{\pm 0.17}$ | $38.64_{\pm 0.28}$ | $89.62_{\pm 0.07}$ |
| INGO | $84.29_{\pm 0.34}$ | $76.54_{\pm 0.20}$ | $79.09_{\pm 0.15}$ | $53.19_{\pm 0.17}$ | $38.91_{\pm 0.10}$ | $89.90_{\pm 0.13}$ |
| SABO | $\underline{87.31}_{\pm 0.38}$ | $\underline{82.65}_{\pm 0.59}$ | $\underline{79.62}_{\pm 0.07}$ | $\underline{53.55}_{\pm 0.17}$ | $\underline{\textbf{42.29}}_{\pm 2.48}$ | $\underline{91.83}_{\pm 0.16}$ |
| Dimension $d=1000$ |  |  |  |  |  |  |
| CMA-ES | $86.85_{\pm 0.57}$ | $82.21_{\pm 0.36}$ | $78.98_{\pm 0.17}$ | $52.35_{\pm 0.17}$ | $38.40_{\pm 1.83}$ | $90.46_{\pm 0.62}$ |
| MMES | $84.98_{\pm 0.52}$ | $80.86_{\pm 1.95}$ | $76.43_{\pm 0.82}$ | $49.22_{\pm 1.23}$ | $39.82_{\pm 3.43}$ | $91.63_{\pm 0.20}$ |
| BES | $83.11_{\pm 0.11}$ | $75.66_{\pm 0.09}$ | $79.09_{\pm 0.08}$ | $53.19_{\pm 0.17}$ | $38.57_{\pm 0.13}$ | $89.61_{\pm 0.04}$ |
| INGO | $84.36_{\pm 0.23}$ | $76.35_{\pm 0.14}$ | $78.97_{\pm 0.08}$ | $53.07_{\pm 0.29}$ | $39.05_{\pm 0.06}$ | $89.95_{\pm 0.08}$ |
| SABO | $\underline{\textbf{87.96}}_{\pm 0.83}$ | $\underline{\textbf{82.77}}_{\pm 0.41}$ | $\underline{\textbf{79.68}}_{\pm 0.23}$ | $\underline{53.31}_{\pm 0.17}$ | $\underline{40.32}_{\pm 0.27}$ | $\underline{\mathbf{91.96}}_{\pm 0.41}$ |

#### Results.

Table [1](#S6.T1) presents experimental results on these six benchmark datasets for three different dimensions of the vector $\bm{v}$. We can see that the SABO method consistently outperforms all baselines in terms of testing classification accuracy or testing F1 scores across different settings, highlighting its effectiveness in improving generalization performance. Notably, even in the high-dimensional setting (i.e., $d=1000$), our method maintains good performance. Moreover, we can see that when $d=1000$, SABO achieves the best performance on SST-2, AG’s News, MRPC, and Yelp P. datasets. The SABO method also achieves the best performance on RTE and SNLI datasets with $d=200$ and $d=500$, respectively.

## 7 Conclusion

In this work, we have introduced SABO, a novel black-box optimization algorithm that improves generalization by utilizing a sharpness-aware minimization strategy. Theoretically, we provide a convergence guarantee for the proposed SABO algorithm in both full-batch function query and mini-batch function query settings. Additionally, we prove the generalization bound for the proposed method. Empirical studies on synthetic numerical problems verify the convergence properties of the proposed method. Moreover, extensive experimental results on black-box prompt fine-tuning problems demonstrate the effectiveness of the proposed SABO method in improving the generalization performance.

## Appendix

## Appendix A Additional Material for Section 3

### A.1 Determine the perturbation

The Lagrangian of problem ([9](#S3.E9)) is

$$ $\displaystyle\mathcal{L}(\bm{\delta},\lambda)$ $\displaystyle=-\left<\nabla_{\bm{\theta}}J(\bm{\theta}),\bm{\delta}\right>+\lambda(\mathrm{KL}(p_{\bm{\theta}+\bm{\delta}}\|p_{\bm{\theta}})-\rho^{2})$ (29) $\displaystyle\begin{split}&=-\bm{\delta}_{\bm{\mu}}^{\top}\nabla_{\bm{\mu}}J(\bm{\theta}_{t})-\text{tr}(\bm{\delta}_{\bm{\Sigma}}\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t}))\\ &\ ~{}~{}~{}~{}~{}+\frac{\lambda}{2}\left[\text{tr}(\bm{\Sigma}^{-1}(\bm{\Sigma}+\bm{\delta}_{\bm{\Sigma}}))+\bm{\delta}_{\bm{\mu}}^{\top}\bm{\Sigma}^{-1}\bm{\delta}_{\bm{\mu}}+\log\frac{|\bm{\Sigma}|}{|(\bm{\Sigma}+\bm{\delta}_{\bm{\Sigma}})|}-d\right]-\lambda\rho^{2}.\end{split}$ (30) $$

Taking the derivative $\bm{\delta}_{\bm{\mu}}$ and $\bm{\delta}_{\bm{\Sigma}}$ and setting them to zero, we can obtain that

$$ $\displaystyle\nabla_{\bm{\mu}}J(\bm{\theta})-\lambda\bm{\Sigma}^{-1}\bm{\delta}_{\bm{\mu}}=0,$ (31) $\displaystyle\nabla_{\bm{\Sigma}}J(\bm{\theta})-\frac{\lambda}{2}[(\bm{\Sigma}+\bm{\delta}_{\bm{\Sigma}})^{-1}-\bm{\Sigma}^{-1}]=0.$ (32) $$

Note that $\bm{\Sigma}$ is a diagonal matrix. Therefore, we can achieve that

$$ $\displaystyle\bm{\delta}_{\bm{\mu}}(\bm{\theta})=\frac{1}{\lambda}\bm{\Sigma}\nabla_{\bm{\mu}}J(\bm{\theta}),$ (33) $\displaystyle\bm{\delta}_{\bm{\Sigma}}(\bm{\theta})=\frac{2\bm{\Sigma}\nabla_{\bm{\Sigma}}J(\bm{\theta})}{\lambda\bm{\Sigma}^{-1}-2\nabla_{\bm{\Sigma}}J(\bm{\theta})}.$ (34) $$

### A.2 Determine the optimal solution

Note that we have

$$ $\displaystyle\mathrm{KL}(p_{\bm{\theta}+\bm{\delta}}\|p_{\bm{\theta}})$ $\displaystyle=\frac{1}{2}\left[\text{tr}(\bm{\Sigma}^{-1}(\bm{\Sigma}+\bm{\delta}_{\bm{\Sigma}}))+\bm{\delta}_{\bm{\mu}}^{\top}\bm{\Sigma}^{-1}\bm{\delta}_{\bm{\mu}}+\log\frac{|\bm{\Sigma}|}{|(\bm{\Sigma}+\bm{\delta}_{\bm{\Sigma}})|}-d\right]$ (35) $\displaystyle=\frac{1}{2}\left[\text{tr}(I+\frac{2}{\lambda}\nabla_{\bm{\Sigma}}J(\bm{\theta})(\bm{\Sigma}+\bm{\delta}_{\bm{\Sigma}}))+\bm{\delta}_{\bm{\mu}}^{\top}\bm{\Sigma}^{-1}\bm{\delta}_{\bm{\mu}}+\log\frac{|\bm{\Sigma}|}{|(\bm{\Sigma}+\bm{\delta}_{\bm{\Sigma}})|}-d\right]$ (36) $\displaystyle=\frac{1}{2}\left[\text{tr}(\frac{2}{\lambda}\nabla_{\bm{\Sigma}}J(\bm{\theta})\bm{\delta}_{\bm{\Sigma}})+\bm{\delta}_{\bm{\mu}}^{\top}\bm{\Sigma}^{-1}\bm{\delta}_{\bm{\mu}}+Q\right],$ (37) $$

where

$$ $\displaystyle Q=|\frac{2}{\lambda}\bm{\Sigma}\nabla_{\bm{\Sigma}}J(\bm{\theta})|+\log(I-\frac{2}{\lambda}\bm{\Sigma}\nabla_{\bm{\Sigma}}J(\bm{\theta})).$ (38) $$

Since $\bm{\Sigma}$ and $\nabla_{\bm{\Sigma}}J(\bm{\theta})$ are both diagonal matrix, we denote $\mathrm{diag}(\bm{\Sigma}\nabla_{\bm{\Sigma}}J(\bm{\theta}))=(v^{1},\dots,v^{d})$, then we have

$$ $Q=\log(\prod_{i=1}^{d}(1-\frac{2}{\lambda}v^{i}))+\sum_{i=1}^{d}\frac{2}{\lambda}v^{i}=\sum_{i=1}^{d}\Big{(}\log(1-\frac{2}{\lambda}v^{i})+\frac{2}{\lambda}v^{i}\Big{)}=\sum_{i=1}^{d}-\frac{2}{\lambda^{2}}(v^{i})^{2}-\mathcal{O}(\frac{1}{\lambda^{3}}(v^{i})^{3}).$ (39) $$

We denote $\mathrm{diag}(\nabla_{\bm{\Sigma}}J(\bm{\theta}))=(\hat{v}^{1},\dots,\hat{v}^{d})$, then we have $\hat{v}^{i}\sigma^{i}=v^{i}$ and

$$ $\text{tr}(\frac{2}{\lambda}\nabla_{\bm{\Sigma}}J(\bm{\theta})\bm{\delta}_{\bm{\Sigma}})=\sum_{i=1}^{d}\frac{2}{\lambda}\hat{v}^{i}(\frac{1}{(\sigma^{i})^{-1}-\frac{2}{\lambda}\hat{v}^{i}}-\sigma^{i})=\sum_{i=1}^{d}\frac{2}{\lambda}\hat{v}^{i}\sigma^{i}(\frac{2}{\lambda}\sigma^{i}\hat{v}^{i})+\mathcal{O}(\frac{4}{\lambda^{2}}(\sigma^{i}\hat{v}^{i})^{2}).$ (40) $$

Then substituting $\bm{\delta}_{\bm{\mu}}(\bm{\theta})$ and $\bm{\delta}_{\bm{\Sigma}}(\bm{\theta})$ into the inequality $\mathrm{KL}(p_{\bm{\theta}+\bm{\delta}}\|p_{\bm{\theta}_{t}})\leq\rho^{2}$, we can obtain that

$$ $\displaystyle\mathrm{KL}(p_{\bm{\theta}+\bm{\delta}}\|p_{\bm{\theta}})$ $\displaystyle=\frac{1}{2}\left[\frac{2}{\lambda^{2}}\|\bm{\Sigma}\nabla_{\bm{\Sigma}}J(\bm{\theta})\|_{\mathrm{F}}^{2}+\frac{1}{\lambda^{2}}\|\bm{\Sigma}^{\frac{1}{2}}\nabla_{\bm{\mu}}J(\bm{\theta})\|^{2}_{2}\right]+\epsilon\leq\rho^{2},$ (41) $$

where $\epsilon=\mathcal{O}(\frac{4(\sigma^{i}\hat{v}^{i})^{2}}{\lambda^{2}})=\mathcal{O}(\frac{1}{\lambda^{2}})$. Let the equality holds and solve Eq. ([41](#A1.E41)), we have

$$ $\displaystyle\lambda\approx\frac{1}{\rho}\sqrt{\|\bm{\Sigma}\nabla_{\bm{\Sigma}}J(\bm{\theta})\|_{\mathrm{F}}^{2}+0.5\|\bm{\Sigma}^{\frac{1}{2}}\nabla_{\bm{\mu}}J(\bm{\theta})\|^{2}_{2}}.$ (42) $$

### A.3 Update rule

Note that we have

$$ $\begin{split}&\left\langle\bm{\theta}-\bm{\theta}_{t},\nabla_{\bm{\theta}}J(\bm{\theta}_{t}+\bm{\delta}_{t})\right\rangle+\frac{1}{\beta_{t}}\mathrm{KL}(p_{\bm{\theta}}\|p_{\bm{\theta}_{t}})\\ &=(\bm{\mu}-\bm{\mu}_{t})^{\top}\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})+\text{tr}((\bm{\Sigma}-\bm{\Sigma}_{t})\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t}+\bm{\delta}_{t}))\\ &\ ~{}~{}~{}~{}~{}+\frac{1}{2\beta_{t}}\left[\text{tr}(\bm{\Sigma}_{t}^{-1}\bm{\Sigma})+(\bm{\mu}-\bm{\mu}_{t})^{\top}\bm{\Sigma}_{t}^{-1}(\bm{\mu}-\bm{\mu}_{t})+\log\frac{|\bm{\Sigma}_{t}|}{|\bm{\Sigma}|}-d\right],\end{split}$ $$

where $\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$ and $\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$ denotes the derivative w.r.t. $\bm{\mu}$ and $\bm{\Sigma}$ taking at $\bm{\mu}=\bm{\mu}_{t}+\bm{\delta}_{\bm{\mu}_{t}}$ and $\bm{\Sigma}=\bm{\Sigma}_{t}+\bm{\delta}_{\bm{\Sigma}_{t}}$, respectively. We can see the above problem is convex with respect to $\bm{\mu}$ and $\bm{\Sigma}$. Taking the derivative w.r.t. $\bm{\mu}$ and $\bm{\Sigma}$ and setting them to zero, we can obtain that

$$ $\displaystyle\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})+\frac{1}{\beta_{t}}\bm{\Sigma}_{t}^{-1}(\bm{\mu}-\bm{\mu}_{t})=0,$ (43) $\displaystyle\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t}+\bm{\delta}_{t})+\frac{1}{2\beta_{t}}[\bm{\Sigma}_{t}^{-1}-\bm{\Sigma}^{-1}]=0.$ (44) $$

Therefore, we obtain the update rule for $\bm{\theta}_{t}$ as

$$ $\displaystyle\bm{\mu}_{t+1}$ $\displaystyle=\bm{\mu}_{t}-\beta_{t}\bm{\Sigma}_{t}\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t}),$ (45) $\displaystyle\bm{\Sigma}_{t+1}^{-1}$ $\displaystyle=\bm{\Sigma}_{t}^{-1}+2\beta_{t}\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t}+\bm{\delta}_{t}).$ (46) $$

### A.4 The gradients of the reparameterized objective

To obtain the gradients of the reparameterized objective $J(\bm{\theta})$ w.r.t. $\bm{\mu}$ and $\bm{\Sigma}$, we use the following theorem to show that only function queries are needed.

###### Theorem A.1 .

[]
The gradient of the expectation of an integrable function $F(\bm{x})$ under a Gaussian distribution $p_{\bm{\theta}}:=\mathcal{N}(\bm{\mu},\bm{\Sigma})$ with respect to the mean $\bm{\mu}$ and the covariance $\bm{\Sigma}$ can be expressed as

$$ $\displaystyle\nabla_{\bm{\mu}}J(\bm{\theta})$ $\displaystyle=\mathbb{E}_{p_{\bm{\theta}}}[\bm{\Sigma}^{-1}(\bm{x}-\bm{\mu})F(\bm{x})],$ (47) $\displaystyle\begin{split}\nabla_{\bm{\Sigma}}J(\bm{\theta})&=\frac{1}{2}\mathbb{E}_{p_{\bm{\theta}}}[\big{(}\bm{\Sigma}^{-1}(\bm{x}-\bm{\mu})(\bm{x}-\bm{\mu})^{\top}\bm{\Sigma}^{-1}-\bm{\Sigma}^{-1}\big{)}F(\bm{x})].\end{split}$ (48) $$

###### Theorem A.1 .

## Appendix B Mini-batch SABO

For the proposed SABO with a mini-batch function query, the stochastic approximation of the gradients $\nabla_{\bm{\mu}}J(\bm{\theta}_{t})$ and $\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})$ using Monte Carlo sampling are given as

$$ $\displaystyle\bm{g}_{t}^{\prime}$ $\displaystyle=\frac{1}{NM}\sum_{j=1}^{N}\sum_{i=1}^{M}\bm{\Sigma}_{t}^{-1}(\bm{x}_{j}^{\prime}-\bm{\mu}_{t})\big{(}F(\bm{x}_{j}^{\prime};(X_{i}^{\prime},y_{i}^{\prime}))-F(\bm{\mu}_{t};(X_{i}^{\prime},y_{i}^{\prime}))\big{)},$ (49) $\displaystyle\bm{G}_{t}^{\prime}$ $\displaystyle=\frac{1}{2NM}\sum_{j=1}^{N}\sum_{i=1}^{M}\mathrm{diag}\Big{[}\bm{\Sigma}_{t}^{-1}\big{[}\mathrm{diag}\big{(}(\bm{x}_{j}^{\prime}-\bm{\mu}_{t})(\bm{x}_{j}^{\prime}-\bm{\mu}_{t})^{\top}\bm{\Sigma}_{t}^{-1}-\bm{I}\big{)}$ $\displaystyle\ ~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}\times(F(\bm{x}_{j}^{\prime};(X_{i}^{\prime},y_{i}^{\prime}))-F(\bm{\mu}_{t};(X_{i}^{\prime},y_{i}^{\prime})))\big{]}\Big{]},$ (50) $$

where $\bm{x}_{j}^{\prime}$ denotes the $j$-th sample sampled from the distribution $\mathcal{N}(\bm{x}\mid\bm{\mu}_{t},\bm{\Sigma}_{t})$, and $(X_{i}^{\prime},y_{i}^{\prime})$ denotes the $i$-th data in the mini-batch dataset $\mathcal{B}^{\prime}$. Then by setting $\bm{\delta}_{t}=\{\frac{1}{\lambda}\bm{\Sigma}_{t}\bm{g}_{t}^{\prime},\frac{2\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}}{\lambda\bm{\Sigma}_{t}^{-1}-2\bm{G}_{t}^{\prime}}\},$ the gradients $\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$ and $\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$ can also be approximated as

$$ $\displaystyle\bm{g}_{t}$ $\displaystyle=\frac{1}{NM}\sum_{j=1}^{N}\sum_{i=1}^{M}\widehat{\bm{\Sigma}}_{t}^{-1}(\bm{x}_{j}-\widehat{\bm{\mu}}_{t})\big{(}F(\bm{x}_{j};(X_{i},y_{i}))-F(\widehat{\bm{\mu}}_{t};(X_{i},y_{i}))\big{)},$ (51) $\displaystyle\bm{G}_{t}$ $\displaystyle=\frac{1}{2NM}\sum_{j=1}^{N}\sum_{i=1}^{M}\mathrm{diag}\Big{[}\widehat{\bm{\Sigma}}_{t}^{-1}\big{[}\mathrm{diag}\big{(}(\bm{x}_{j}-\widehat{\bm{\mu}}_{t})(\bm{x}_{j}-\widehat{\bm{\mu}}_{t})^{\top}\widehat{\bm{\Sigma}}_{t}^{-1}-\bm{I}\big{)}$ $\displaystyle\ ~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}\times(F(\bm{x}_{j};(X_{i},y_{i}))-F(\widehat{\bm{\mu}}_{t};(X_{i},y_{i})))\big{]}\Big{]},$ (52) $$

where $\widehat{\bm{\Sigma}}_{t}=\bm{\Sigma}_{t}+\bm{\delta}_{\bm{\Sigma}_{t}}$, $\widehat{\bm{\mu}}_{t}=\bm{\mu}_{t}+\bm{\delta}_{\bm{\mu}_{t}}$, $\bm{x}_{j}$ denotes the $j$-th sample sampled from the distribution $\mathcal{N}(\bm{x}\mid\widehat{\bm{\mu}}_{t},\widehat{\bm{\Sigma}}_{t})$, and $(X_{i},y_{i})$ denotes the $i$-th data in the mini-batch dataset $\mathcal{B}$. Then we can update $\bm{\theta}_{t}$ by Eq. ([22](#S3.E22)). The entire algorithm of SABO with a mini-batch function query is shown in Algorithm [2](#alg2).

Figure: Algorithm 2 Mini-batch SABO

## Appendix C Technical Lemmas

In this section, we introduce the following technical lemmas for analysis. The proof of all technical lemmas is put in Appendix [E](#A5).

###### Lemma C.1 .

Suppose $\bm{\Sigma}$ and $\hat{\bm{\Sigma}}$ are two $d$-dimensional diagonal matrix and $\bm{z}$ is a $d$-dimensional vector, then we have $\|\bm{\Sigma}\bm{z}\|\leq\|\bm{\Sigma}\|_{\mathrm{F}}\|\bm{z}\|$ and $\|\bm{\Sigma}\hat{\bm{\Sigma}}\|_{\mathrm{F}}\leq\|\bm{\Sigma}\|_{\mathrm{F}}\|\hat{\bm{\Sigma}}\|_{\mathrm{F}}$.

###### Lemma C.2 .

Given a convex function $F(\bm{x})$, for Gaussian distribution with parameters ${\bm{\theta}}:=\{{\bm{\mu}},\bm{\Sigma}^{\frac{1}{2}}\}$, let $J({\bm{\theta}}):=\mathbb{E}_{p({\bm{x}};{\bm{\theta}})}[F({\bm{x}})]$. Then $J({\bm{\theta}})$ is a convex function with respect to ${\bm{\theta}}$.

###### Lemma C.3 .

Suppose that the gradient $\bm{G}_{t}$ are positive semi-definite matrix and satisfies $\xi\bm{I}\preceq\bm{G}_{t}\preceq b\bm{I}$. Then for algorithm 1 and 2, we have the following results.

- (a)
The (diagonal) covariance matrix $\bm{\Sigma}_{T}$ satisfies $\frac{1}{2b\sum_{t=1}^{T}\beta_{t}\bm{I}+\bm{\Sigma}_{0}^{-1}}\preceq\bm{\Sigma}_{T}\preceq\frac{1}{2\xi\sum_{t=1}^{T}\beta_{t}\bm{I}+\bm{\Sigma}_{0}^{-1}}$.
- (b)
The Frobenius norm for the covariance matrix $\bm{\Sigma}_{t}$ satisfies $\|\bm{\Sigma}_{t}\|_{\mathrm{F}}\leq\frac{\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}}$.

###### Lemma C.4 .

For given $\bm{\theta}_{t}$, denote the approximated gradients of $\nabla_{\bm{\mu}}J(\bm{\theta}_{t})$ and $\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})$ by $\bm{g}_{t}^{\prime}$ and $\bm{G}_{t}^{\prime}$, respectively. Then for Algorithm [1](#alg1) and Algorithm [2](#alg2), if $\rho\leq\frac{\sqrt{d}}{2}$, then the perturbation satisfies, $\|\bm{\delta}_{\bm{\mu}_{t}}\|\leq\rho\sqrt{2}\|\bm{\Sigma}^{\frac{1}{2}}_{t}\|_{\mathrm{F}}$ and $\|\bm{\delta}_{\bm{\Sigma}_{t}}\|_{\mathrm{F}}\leq 2\rho r\|\bm{\Sigma}_{t}\|_{\mathrm{F}}$, where $r$ is a positive constant.

###### Lemma C.5 .

In Algorithm [1](#alg1), suppose the gradient estimator $\bm{g}_{t}^{\prime}$ in $t$-th iteration as

$$ $\displaystyle\bm{g}_{t}^{\prime}=\bm{\Sigma}_{t}^{-\frac{1}{2}}\bm{z}\big{(}F(\bm{\mu}_{t}+\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{z})-F(\bm{\mu}_{t})\big{)},$ (53) $$

where $\bm{z}\sim\mathcal{N}(0,I)$. Then $\bm{g}_{t}^{\prime}$ is an unbiased estimator of the gradient $\nabla_{\bm{\mu}}\mathbb{E}_{p_{\bm{\theta}_{t}}}[F(\bm{x})]$.

###### Lemma C.6 .

In Algorithm [1](#alg1), suppose the gradient estimator $\bm{g}_{t}$ in $t$-th iteration as

$$ $\displaystyle\bm{g}_{t}=\widehat{\bm{\Sigma}}_{t}^{-\frac{1}{2}}\bm{z}\big{(}F(\widehat{\bm{\mu}}_{t}+\widehat{\bm{\Sigma}}_{t}^{\frac{1}{2}}\bm{z})-F(\widehat{\bm{\mu}}_{t})\big{)},$ (54) $$

where $\bm{z}\sim\mathcal{N}(0,I)$. Suppose that Assumption [4.1](#S4.Thmtheorem1) holds, $\rho<\frac{\sqrt{d}}{2}$, the gradient $\bm{G}_{t}$ is positive semi-definite matrix and satisfies $\xi\bm{I}\preceq\bm{G}_{t}\preceq b\bm{I}$ and $\bm{\Sigma}_{0}\preceq R\bm{I}$, where $\xi,b,R\geq 0$. Then we have

$$ $\displaystyle\mathbb{E}[\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}]\leq\frac{H^{2}(1+2\rho r^{\prime})(d+4)^{2}}{2\xi(\sum_{k=1}^{t}\beta_{k})},$ (55) $$

where $r^{\prime}$ is a positive constant.

###### Lemma C.7 .

In Algorithm [2](#alg2), suppose assumption [4.4](#S4.Thmtheorem4) holds and the gradient estimator $\bm{g}_{t}^{\prime}$ in $t$-th iteration as

$$ $\displaystyle\bm{g}_{t}^{\prime}=\bm{\Sigma}_{t}^{-\frac{1}{2}}\bm{z}\big{(}F(\bm{\mu}_{t}+\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{z};\mathcal{B})-F(\bm{\mu}_{t};\mathcal{B})\big{)},$ (56) $$

where $\bm{z}\sim\mathcal{N}(0,I)$. Then $\bm{g}_{t}^{\prime}$ is an unbiased estimator of the gradient $\nabla_{\bm{\mu}}\mathbb{E}_{p_{\bm{\theta}_{t}}}[F(\bm{x})]$.

###### Lemma C.8 .

In Algorithm [2](#alg2), suppose the gradient estimator $\bm{g}_{t}$ in $t$-th iteration as

$$ $\displaystyle\bm{g}_{t}=\widehat{\bm{\Sigma}}_{t}^{-\frac{1}{2}}\bm{z}\big{(}F(\widehat{\bm{\mu}}_{t}+\widehat{\bm{\Sigma}}_{t}^{\frac{1}{2}}\bm{z};\mathcal{B})-F(\widehat{\bm{\mu}}_{t};\mathcal{B})\big{)},$ (57) $$

where $\bm{z}\sim\mathcal{N}(0,I)$. Suppose that Assumption [4.1](#S4.Thmtheorem1) and Assumption [4.4](#S4.Thmtheorem4) hold, $\rho<\frac{\sqrt{d}}{2}$, the gradient $\bm{G}_{t}$ is positive semi-definite matrix and satisfies $\xi\bm{I}\preceq\bm{G}_{t}\preceq b\bm{I}$ and $\bm{\Sigma}_{0}\preceq R\bm{I}$, where $\xi,b,R\geq 0$. Then we have

$$ $\displaystyle\mathbb{E}[\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}]\leq\frac{H^{2}(1+2\rho r^{\prime})(d+4)^{2}}{6\xi(\sum_{k=1}^{t}\beta_{k})}+\frac{2(d+4)\varepsilon^{2}}{3},$ (58) $$

where $r^{\prime}$ is a positive constant, $\varepsilon^{2}=2\varepsilon^{2}_{\mathrm{B}}+2\varepsilon^{2}_{\mathrm{D}}$.

###### Lemma C.1 .

###### Lemma C.2 .

###### Lemma C.3 .

###### Lemma C.4 .

###### Lemma C.5 .

###### Lemma C.6 .

###### Lemma C.7 .

###### Lemma C.8 .

## Appendix D Proof of the Result in Section 4

In this section, we provide the proof of the result in Section [4](#S4).

### D.1 Proof of the Proposition 4.2

Note that $F(\bm{x})$ is a convex function, we have

$$ $\displaystyle F(\bm{\mu})=F(\mathbb{E}_{\bm{x}\sim\mathcal{N}(\bm{\mu},\bm{\Sigma})}[\bm{x}])\leq\mathbb{E}_{\bm{x}\sim\mathcal{N}(\bm{\mu},\bm{\Sigma})}[F(\bm{x})]=J(\bm{\mu},\bm{\Sigma}).$ (59) $$

Since $F(\bm{\mu}^{*})=J(\bm{\mu}^{*},\bm{0})$, it follows that

$$ $\displaystyle F(\bm{\mu})-F(\bm{\mu}^{*})\leq J(\bm{\mu},\bm{\Sigma})-J(\bm{\mu}^{*},\bm{0}),$ (60) $$

where we reach the conclusion.

### D.2 Proof of Theorem 4.3

The update rule of $\bm{\mu}$ can be represented as $\bm{\mu}_{t+1}=\bm{\mu}_{t}-\beta_{t}\bm{\Sigma}_{t}\bm{g}_{t}$. Since $F(\bm{x})$ is convex function, we have $J(\bm{\theta})$ is convex w.r.t. $\bm{\theta}=\{\bm{\mu},\bm{\Sigma}^{\frac{1}{2}}\}$ by Lemma [C.2](#A3.Thmtheorem2), together with $J(\bm{\theta})$ is $c$-strongly convex w.r.t. $\bm{\mu}$ we obtain

$$ $\displaystyle J(\bm{\theta}_{t})\leq J(\bm{\mu}^{*},0)+\nabla_{\bm{\mu}}J(\bm{\theta}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})+\nabla_{\bm{\Sigma}^{\frac{1}{2}}}J(\bm{\theta}_{t})^{\top}\bm{\Sigma}_{t}^{\frac{1}{2}}-\frac{c}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}.$ (61) $$

Note that $\nabla_{\bm{\Sigma}^{\frac{1}{2}}}J(\bm{\theta}_{t})=\bm{\Sigma}_{t}^{\frac{1}{2}}\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})+\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})\bm{\Sigma}_{t}^{\frac{1}{2}}$, we have

$$ $\displaystyle J(\bm{\theta}_{t})\leq J(\bm{\mu}^{*},0)+\nabla_{\bm{\mu}}J(\bm{\theta}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})+2\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})\bm{\Sigma}_{t}-\frac{c}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}.$ (62) $$

Let $A_{t}=J(\bm{\mu}_{t},\bm{\Sigma}_{t})-J(\bm{\mu}^{*},0)$, we have

$$ $\displaystyle\beta_{t}\mathbb{E}[A_{t}]\leq\beta_{t}\mathbb{E}_{\bm{z}}[\nabla_{\bm{\mu}}J(\bm{\theta}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})]+2\beta_{t}\mathbb{E}_{\bm{z}}[\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})^{\top}\bm{\Sigma}_{t}]-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}.$ (63) $$

Note that

$$ $\displaystyle\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}$ $\displaystyle=\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t}-\beta_{t}\bm{\Sigma}_{t}\bm{g}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}$ (64) $\displaystyle=\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\big{(}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-2\beta_{t}\left<\bm{\mu}_{t}-\bm{\mu}^{*},\bm{g}_{t}\right>+\beta_{t}^{2}\left<\bm{\Sigma}_{t}\bm{g}_{t},\bm{g}_{t}\right>\big{)}$ (65) $\displaystyle=2\beta_{t}\bm{g}_{t}^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})-\beta_{t}^{2}(\bm{\Sigma}_{t}\bm{g}_{t})^{\top}\bm{g}_{t}.$ (66) $$

Therefore we have

$$ $\displaystyle\bm{g}_{t}^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})$ $\displaystyle=\frac{1}{2\beta_{t}}\big{(}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}\big{)}+\frac{\beta_{t}}{2}(\bm{\Sigma}_{t}\bm{g}_{t})^{\top}\bm{g}_{t}.$ (67) $$

According to Lemma [C.5](#A3.Thmtheorem5), we have $\mathbb{E}\bm{g}_{t}=\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$. Then we have

$$ $\displaystyle\mathbb{E}\big{[}(\nabla_{\bm{\mu}}J(\bm{\theta}_{t})-\bm{g}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}$ $\displaystyle=\mathbb{E}\big{[}(\nabla_{\bm{\mu}}J(\bm{\theta}_{t})-\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t}))^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}+\mathbb{E}\big{[}(\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})-\bm{g}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}$ (68) $\displaystyle\leq\|\bm{\mu}_{t}-\bm{\mu}^{*}\|\|\nabla_{\bm{\mu}}J(\bm{\theta}_{t})-\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})\|$ (69) $\displaystyle\leq DL\|\bm{\delta}_{t}\|$ (70) $\displaystyle\leq DL\rho_{t}U_{t},$ (71) $$

where the first inequality is due to the Cauchy-Schwarz inequality, the second inequality is due to $\|\bm{\mu}_{t}-\bm{\mu}^{*}\|\leq D$ and Assumption [4.1](#S4.Thmtheorem1), the last inequality is due to Lemma [C.4](#A3.Thmtheorem4), and $U_{t}=\max(2\sqrt{2}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}},4r\|\bm{\Sigma}_{t}\|_{\mathrm{F}})$. Then we have

$$ $\displaystyle\mathbb{E}_{\bm{z}}\big{[}\nabla_{\bm{\mu}}J(\bm{\theta}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}$ (72) $\displaystyle=\mathbb{E}_{\bm{z}}\big{[}\bm{g}_{t}^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})+(\nabla_{\bm{\mu}}J(\bm{\theta}_{t})-\bm{g}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}$ (73) $\displaystyle\leq\frac{1}{2\beta_{t}}\mathbb{E}_{\bm{z}}[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}]+\beta_{t}(\bm{\Sigma}_{t}\bm{g}_{t})^{\top}\bm{g}_{t}]+DL\rho_{t}U_{t},$ (74) $$

where the inequality is due to Eq. ([67](#A4.E67)), Eq. ([71](#A4.E71)), and $\beta_{t}\geq 0$. Note that

$$ $\displaystyle\mathbb{E}\big{[}\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})^{\top}\bm{\Sigma}_{t}\big{]}\leq\mathbb{E}\big{[}\|\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})\|\big{]}\|\bm{\Sigma}_{t}\|_{\mathrm{F}}\leq H\|\bm{\Sigma}_{t}\|_{\mathrm{F}},$ (75) $$

where the first inequality is due to Lemma [C.1](#A3.Thmtheorem1) and the second inequality is due to the Lipschitz continuous assumption of the function $J(\bm{\theta})$.
Then substituting Eq. ([74](#A4.E74)) and Eq. ([75](#A4.E75)) into Eq. ([63](#A4.E63)) and multiplying $\beta_{t}$ on both sides of the inequality, we have

$$ $\begin{split}\beta_{t}\mathbb{E}[A_{t}]&\leq\frac{1}{2}\mathbb{E}[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}]-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}+2\beta_{t}^{2}H\|\bm{\Sigma}_{t}\|_{\mathrm{F}}\\ &\ ~{}~{}~{}~{}~{}~{}~{}~{}~{}+DL\rho_{t}U_{t}+\beta_{t}^{2}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2},\end{split}$ (76) $$

We further obtain that

$$ $\displaystyle\sum_{t=0}^{T-1}\left[\frac{1}{2}\mathbb{E}[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}]-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}\right]$ (77) $\displaystyle\leq\frac{1}{2}\sum_{t=0}^{T-1}\left[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t-1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t-1}^{-1}}^{2}-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}\right]$ (78) $\displaystyle\ ~{}~{}~{}~{}~{}~{}~{}~{}+\frac{1}{2}\left[\|\bm{\mu}_{0}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{0}^{-1}}^{2}-\|\bm{\mu}_{T}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{T-1}^{-1}}^{2}\right]$ (79) $\displaystyle\leq\frac{1}{2}\sum_{t=0}^{T-1}\left[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{2\beta_{t}\bm{G}_{t}}^{2}-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}\right]+\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}$ (80) $\displaystyle\leq\frac{1}{2}\sum_{t=0}^{T-1}\left[\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}\right]+\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}$ (81) $\displaystyle=\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2},$ (82) $$

where the second inequality is due to the update rule of $\bm{\Sigma}_{t}$ and $\|\bm{\mu}_{t}-\bm{\mu}^{*}\|\leq D$, and the third inequality is due to Cauchy-Schwarz inequality and $\bm{G}_{t}\preceq\frac{c}{4}\bm{I}$. Then we have

$$ $\displaystyle\sum_{t=0}^{T-1}\beta_{t}\mathbb{E}[A_{t}]$ $\displaystyle\leq\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}+\sum_{t=0}^{T-1}\Bigg{(}2H\beta_{t}^{2}\|\bm{\Sigma}_{t}\|_{\mathrm{F}}+DL\rho_{t}U_{t}+\beta_{t}^{2}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}\Bigg{)}$ (83) $\displaystyle\leq\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}+\sum_{t=0}^{T-1}\Bigg{(}\frac{2H\beta_{t}^{2}\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}}+DL\rho_{t}U_{t}+\frac{H^{2}R^{\frac{1}{2}}(1+2\rho r^{\prime})(d+4)^{2}\beta_{t}^{2}}{2\xi(\sum_{k=1}^{t}\beta_{k})}\Bigg{)},$ (84) $$

where the first inequality is due to Eq. ([82](#A4.E82)) and Eq. ([76](#A4.E76)), and the second inequality is due to Lemma [C.6](#A3.Thmtheorem6) and $\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}\leq R^{\frac{1}{2}}$.
According to [C.3](#A3.Thmtheorem3) (b), we have $U_{t}\leq\max(\frac{2\sqrt{2}d^{\frac{1}{4}}}{(2\xi\sum_{k=1}^{t}\beta_{k})^{\frac{1}{2}}},\frac{4r\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}})$. Therefore, there exists a constant $t^{*}$, when $t>t^{*}-1$, $U_{t}\leq\frac{2\sqrt{2}d^{\frac{1}{4}}}{(2\xi\sum_{k=1}^{t}\beta_{k})^{\frac{1}{2}}}$. Denote $\sum_{t=0}^{t^{*}-1}\frac{4r\rho_{t}\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}}$ by a constant $\Gamma$. Then if $T>t^{*}$, we have

$$ $\begin{split}\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}[A_{t}]\leq\frac{\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}}{T\beta_{t}}&+\frac{1}{T}\sum_{t=0}^{T-1}\Bigg{(}\frac{2H\beta_{t}\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}}+\frac{H^{2}R^{\frac{1}{2}}(1+2\rho r^{\prime})(d+4)^{2}\beta_{t}}{2\xi\sum_{k=1}^{t}\beta_{k}}\Bigg{)}\\ &+\frac{1}{T}\sum_{t=t*}^{T-1}\Bigg{(}\frac{2\sqrt{2}DL\rho_{t}d^{\frac{1}{4}}}{\beta_{t}(2\xi\sum_{k=1}^{t}\beta_{k})^{\frac{1}{2}}}\Bigg{)}+\frac{DL\Gamma}{T\beta_{t}},\end{split}$ (85) $$

Let $\beta_{t}=\beta$ and $\rho_{t}=\frac{\rho_{0}}{\sqrt{t}}$, we can obtain that

$$ $\displaystyle\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}[A_{t}]\leq\frac{\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}}{T\beta}+\frac{1}{T}\sum_{t=0}^{T-1}\Bigg{(}\frac{C}{2\xi t}\Bigg{)}+\frac{1}{T}\sum_{t=t^{*}}^{T-1}\Bigg{(}\frac{2\sqrt{2}DL\rho_{0}d^{\frac{1}{4}}}{\sqrt{2\xi}\beta^{\frac{3}{2}}t}\Bigg{)}+\frac{DL\Gamma}{T\beta},$ (86) $$

where $C=2H\sqrt{d}+H^{2}R^{\frac{1}{2}}(1+2\rho r^{\prime})(d+4)^{2}$. Since we have $\sum_{t=1}^{T}\frac{1}{t}\leq 1+\log(T)$, we obtain

$$ $\displaystyle\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}\left[J(\bm{\mu}_{t+1},\bm{\Sigma}_{t})-J(\bm{\mu}^{*},0)\right]=\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}[A_{t}]=\mathcal{O}\left(\frac{1}{T}+\frac{\log T}{T}\right),$ (87) $$

where we reach the conclusion.

###### Remark D.1 .

In Theorem [4.3](#S4.Thmtheorem3), if we set $\beta=\mathcal{O}(1)$ and $\rho=\mathcal{O}(\frac{1}{\sqrt{T}})$, then the third term in right-hand side of Eq. ([85](#A4.E85)) has a convergence rate of $\mathcal{O}(\frac{\log T}{T})$. If we set $\beta=\mathcal{O}(1)$ and $\rho=\mathcal{O}(1)$, then the third term in right-hand side of Eq. ([85](#A4.E85)) has a convergence rate of $\mathcal{O}(\frac{1}{\sqrt{T}})$, since we have $\sum_{t=1}^{T}\frac{1}{\sqrt{t}}\leq 2\sqrt{T}$. Then we obtain

$$ $\displaystyle\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}\left[J(\bm{\mu}_{t+1},\bm{\Sigma}_{t})-J(\bm{\mu}^{*},0)\right]=\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}[A_{t}]=\mathcal{O}\left(\frac{1}{T}+\frac{\log T}{T}+\frac{1}{\sqrt{T}}\right).$ (88) $$

This implies $\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}\left[J(\bm{\mu}_{t+1},\bm{\Sigma}_{t})-J(\bm{\mu}^{*},0)\right]=\mathcal{O}\left(\frac{1}{\sqrt{T}}\right)$.

###### Remark D.1 .

### D.3 Proof of the Proposition 4.5

Since the datasets $\mathcal{D}$ and $\mathcal{B}$ are i.i.d. sampled from one data distribution $P(X,y)$. We have $\mathbb{E}_{(X,y)\sim P}[F(\bm{x};(X,y))]=\mathbb{E}[F(\bm{x};\mathcal{B})]=\mathbb{E}[F(\bm{x};\mathcal{D})]$. We obtain that

$$ $\|F(\bm{x};\mathcal{B})-F(\bm{x};\mathcal{D})\|_{2}^{2}\leq 2\|F(\bm{x};\mathcal{B})-\mathbb{E}F(\bm{x};\mathcal{B})\|_{2}^{2}+2\|F(\bm{x};\mathcal{D})-\mathbb{E}F(\bm{x};\mathcal{D})\|_{2}^{2}=2(\varepsilon_{\mathrm{B}}^{2}+\varepsilon_{\mathrm{D}}^{2})$ (89) $$

### D.4 Proof of Theorem 4.6

Note that for SABO with a mini-batch function query, the update rule of $\bm{\mu}$ can be represented as $\bm{\mu}_{t+1}=\bm{\mu}_{t}-\beta_{t}\bm{\Sigma}_{t}\bm{g}_{t}$.
Let $B_{t}=J(\bm{\mu}_{t},\bm{\Sigma}_{t})-J(\bm{\mu}^{*},0)$, then by using Eq. ([62](#A4.E62)), we have

$$ $\displaystyle\beta_{t}\mathbb{E}[B_{t}]$ $\displaystyle\leq\beta_{t}\mathbb{E}[\nabla_{\bm{\mu}}J(\bm{\theta}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})]+2\beta_{t}\mathbb{E}[\nabla_{\bm{\Sigma}}J(\bm{\theta}_{t})^{\top}\bm{\Sigma}_{t}]-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}$ (90) $\displaystyle\leq\beta_{t}\mathbb{E}[\nabla_{\bm{\mu}}J(\bm{\theta}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})]+2\beta_{t}H\|\bm{\Sigma}_{t}\|_{\mathrm{F}}-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2},$ (91) $$

where the second inequality is due to Lemma [C.1](#A3.Thmtheorem1) and Lipschitz continuous assumption of the function $J(\bm{\theta})$. Note that

$$ $\displaystyle\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}$ $\displaystyle=\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t}-\beta_{t}\bm{\Sigma}_{t}\bm{g}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}$ (92) $\displaystyle=\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\big{(}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-2\beta_{t}\left<\bm{\mu}_{t}-\bm{\mu}^{*},\bm{g}_{t}\right>+\beta_{t}^{2}\left<\bm{\Sigma}_{t}\bm{g}_{t},\bm{g}_{t}\right>\big{)}$ (93) $\displaystyle=2\beta_{t}\bm{g}_{t}^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})-\beta_{t}^{2}(\bm{\Sigma}_{t}\bm{g}_{t})^{\top}\bm{g}_{t}.$ (94) $$

It follows that

$$ $\displaystyle\bm{g}_{t}^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})$ $\displaystyle=\frac{1}{2\beta_{t}}\big{(}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}\big{)}+\frac{\beta_{t}}{2}(\bm{\Sigma}_{t}\bm{g}_{t})^{\top}\bm{g}_{t}$ (95) $\displaystyle\leq\frac{1}{2\beta_{t}}\big{(}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}\big{)}+\beta_{t}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2},$ (96) $$

where the inequality is due to $\beta_{t}\geq 0$ and Lemma [C.1](#A3.Thmtheorem1). According to Lemma [C.7](#A3.Thmtheorem7), we have $\mathbb{E}\bm{g}_{t}=\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})$. Therefore

$$ $\displaystyle\mathbb{E}\big{[}(\nabla_{\bm{\mu}}J(\bm{\theta}_{t})-\bm{g}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}$ $\displaystyle=\mathbb{E}_{\bm{z}}\big{[}(\nabla_{\bm{\mu}}J(\bm{\theta}_{t})-\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t}))^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}+\mathbb{E}_{\bm{z}}\big{[}(\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})-\bm{g}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}$ (97) $\displaystyle\leq\|\bm{\mu}_{t}-\bm{\mu}^{*}\|\|\nabla_{\bm{\mu}}J(\bm{\theta}_{t})-\nabla_{\bm{\mu}}J(\bm{\theta}_{t}+\bm{\delta}_{t})\|$ (98) $\displaystyle\leq DL\|\bm{\delta}_{t}\|$ (99) $\displaystyle\leq DL\rho_{t}U_{t},$ (100) $$

where the first inequality is due to the Cauchy-Schwarz inequality, the second inequality is due to $\|\bm{\mu}_{t}-\bm{\mu}^{*}\|\leq D$ and smoothness assumption of the function $J(\bm{\theta})$, the last inequality is due to Lemma [C.4](#A3.Thmtheorem4), and $U_{t}=\max(2\sqrt{2}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}},4r\|\bm{\Sigma}_{t}\|_{\mathrm{F}})$. Then we have

$$ $\displaystyle\mathbb{E}_{\bm{z}}\big{[}\nabla_{\bm{\mu}}J(\bm{\theta}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}$ (101) $\displaystyle=\mathbb{E}_{\bm{z}}\big{[}\bm{g}_{t}^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})+(\nabla_{\bm{\mu}}J(\bm{\theta}_{t})-\bm{g}_{t})^{\top}(\bm{\mu}_{t}-\bm{\mu}^{*})\big{]}$ (102) $\displaystyle\leq\frac{1}{2\beta_{t}}\mathbb{E}_{\bm{z}}[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}]+\beta_{t}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}+DL\rho_{t}U_{t},$ (103) $$

where the inequality is due to Eq. ([96](#A4.E96)) and Eq. ([100](#A4.E100)).
Then substituting Eq. ([103](#A4.E103)) into Eq. ([91](#A4.E91)) and multiplying $\beta_{t}$ on both sides of the inequality, we have

$$ $\begin{split}\beta_{t}\mathbb{E}[B_{t}]&\leq\frac{1}{2}\mathbb{E}[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}]-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}+2\beta_{t}^{2}H\|\bm{\Sigma}_{t}\|_{\mathrm{F}}\\ &\ ~{}~{}~{}~{}~{}~{}~{}~{}~{}+DL\rho_{t}U_{t}+\beta_{t}^{2}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2},\end{split}$ (104) $$

We further obtain that

$$ $\displaystyle\sum_{t=0}^{T-1}\left[\frac{1}{2}\mathbb{E}[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t+1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}]-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}\right]$ (105) $\displaystyle\leq\frac{1}{2}\sum_{t=0}^{T-1}\left[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t}^{-1}}^{2}-\|\bm{\mu}_{t-1}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{t-1}^{-1}}^{2}-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}\right]$ (106) $\displaystyle\ ~{}~{}~{}~{}~{}~{}~{}~{}+\frac{1}{2}\left[\|\bm{\mu}_{0}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{0}^{-1}}^{2}-\|\bm{\mu}_{T}-\bm{\mu}^{*}\|_{\bm{\Sigma}_{T-1}^{-1}}^{2}\right]$ (107) $\displaystyle\leq\frac{1}{2}\sum_{t=0}^{T-1}\left[\|\bm{\mu}_{t}-\bm{\mu}^{*}\|_{2\beta_{t}\bm{G}_{t}}^{2}-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}\right]+\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}$ (108) $\displaystyle\leq\frac{1}{2}\sum_{t=0}^{T-1}\left[\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}-\frac{c\beta_{t}}{2}\|\bm{\mu}_{t}-\bm{\mu}^{*}\|^{2}\right]+\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}$ (109) $\displaystyle=\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2},$ (110) $$

where the second inequality is due to the update rule of $\bm{\Sigma}_{t}$ and $\|\bm{\mu}_{t}-\bm{\mu}^{*}\|\leq D$, and the third inequality is due to Cauchy-Schwarz inequality and $\bm{G}_{t}\preceq\frac{c}{4}\bm{I}$. Then we have

$$ $\displaystyle\sum_{t=0}^{T-1}\beta_{t}\mathbb{E}[B_{t}]$ $\displaystyle\leq\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}+\sum_{t=0}^{T-1}\Bigg{(}2H\beta_{t}^{2}\|\bm{\Sigma}_{t}\|_{\mathrm{F}}+DL\rho_{t}U_{t}+\beta_{t}^{2}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}\Bigg{)}$ (111) $\displaystyle\leq\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}+\sum_{t=0}^{T-1}\Bigg{(}\frac{2H\beta_{t}^{2}\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}}+DL\rho_{t}U_{t}+\frac{2\beta_{t}^{2}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}(d+4)\varepsilon^{2}}{3}$ $\displaystyle\ ~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}+\frac{H^{2}(1+2\rho r^{\prime})(d+4)^{2}\beta_{t}^{2}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}}{6\xi(\sum_{k=1}^{t}\beta_{k})}\Bigg{)}$ (112) $\displaystyle\leq\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}+\sum_{t=0}^{T-1}\Bigg{(}\frac{2H\beta_{t}^{2}\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}}+DL\rho_{t}U_{t}+\frac{2\beta_{t}^{2}(d+4)d^{\frac{1}{4}}\varepsilon^{2}}{3(2\xi\sum_{k=1}^{t}\beta_{k})^{\frac{1}{2}}}$ $\displaystyle\ ~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}+\frac{H^{2}(1+2\rho r^{\prime})(d+4)^{2}\beta_{t}^{2}R^{\frac{1}{2}}}{6\xi(\sum_{k=1}^{t}\beta_{k})}\Bigg{)},$ (113) $$

where the first inequality is due to Eq. ([104](#A4.E104)) and Eq. ([110](#A4.E110)), the second inequality is due to Lemma [C.3](#A3.Thmtheorem3) (b) and Lemma [C.8](#A3.Thmtheorem8), and the third inequality is due to $\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}\leq R^{\frac{1}{2}}$ and Lemma [C.3](#A3.Thmtheorem3) (b). According to [C.3](#A3.Thmtheorem3) (b), we have $U_{t}\leq\max(\frac{2\sqrt{2}d^{\frac{1}{4}}}{(2\xi\sum_{k=1}^{t}\beta_{k})^{\frac{1}{2}}},\frac{4r\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}})$. Therefore, there exists a constant $t^{*}$, when $t>t^{*}-1$, $U_{t}\leq\frac{2\sqrt{2}d^{\frac{1}{4}}}{(2\xi\sum_{k=1}^{t}\beta_{k})^{\frac{1}{2}}}$. Denote $\sum_{t=0}^{t^{*}-1}\frac{4r\rho_{t}\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}}$ by a constant $\Gamma$. Then if $T>t^{*}$, we have

$$ $\begin{split}\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}[B_{t}]&\leq\frac{\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}}{T\beta_{t}}+\frac{1}{T}\sum_{t=0}^{T-1}\Bigg{(}\frac{2H\beta_{t}\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}}+\frac{H^{2}R^{\frac{1}{2}}(1+2\rho r^{\prime})(d+4)^{2}\beta_{t}}{6\xi\sum_{k=1}^{t}\beta_{k}}\Bigg{)}\\ &\ ~{}~{}~{}~{}+\frac{1}{T}\sum_{t=t*}^{T-1}\Bigg{(}\frac{2\sqrt{2}DL\rho_{t}d^{\frac{1}{4}}}{\beta_{t}(2\xi\sum_{k=1}^{t}\beta_{k})^{\frac{1}{2}}}\Bigg{)}+\frac{1}{T}\sum_{t=0}^{T-1}\Bigg{(}\frac{2\beta_{t}(d+4)d^{\frac{1}{4}}\varepsilon^{2}}{3(2\xi\sum_{k=1}^{t}\beta_{k})^{\frac{1}{2}}}\Bigg{)}+\frac{DL\Gamma}{T\beta_{t}},\end{split}$ (114) $$

Let $\beta_{t}=\beta$ and $\rho_{t}=\rho$, we have

$$ $\begin{split}\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}[B_{t}]&\leq\frac{\|\bm{\Sigma}_{0}^{-1}\|_{\mathrm{F}}D^{2}}{T\beta}+\frac{1}{T}\sum_{t=0}^{T-1}\Bigg{(}\frac{2H\sqrt{d}}{2\xi t}+\frac{H^{2}R^{\frac{1}{2}}(1+2\rho r^{\prime})(d+4)^{2}}{6\xi t}\Bigg{)}\\ &\ ~{}~{}~{}~{}~{}~{}~{}~{}+\frac{1}{T}\sum_{t=t*}^{T-1}\Bigg{(}\frac{2\sqrt{2}DL\rho d^{\frac{1}{4}}}{\sqrt{2\xi}\beta^{\frac{3}{2}}\sqrt{t}}\Bigg{)}+\frac{1}{T}\sum_{t=0}^{T-1}\Bigg{(}\frac{2(d+4)d^{\frac{1}{4}}\varepsilon^{2}}{3\sqrt{2\xi}\sqrt{t}}\Bigg{)}+\frac{DL\Gamma}{T\beta}.\end{split}$ (115) $$

Since we have $\sum_{t=1}^{T}\frac{1}{t}\leq 1+\log(T)$ and $\sum_{t=1}^{T}\frac{1}{\sqrt{t}}\leq 2\sqrt{T}$, we obtain

$$ $\displaystyle\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}\left[J(\bm{\mu}_{t+1},\bm{\Sigma}_{t})-J(\bm{\mu}^{*},0)\right]=\frac{1}{T}\sum_{t=0}^{T-1}\mathbb{E}[B_{t}]=\mathcal{O}\left(\frac{1}{T}+\frac{\log T}{T}+\frac{1}{\sqrt{T}}\right),$ (116) $$

where we reach the conclusion.

### D.5 Proof of Theorem 4.8

Using PAC-Bayesian bound [] and following [], for any prior distribution $p$ and any posterior distribution $q$ that may be dependent on finite dataset $\mathcal{S}$, where data set $S$ with $M$ i.i.d. samples drawn from data distribution $P(X,y)$, we have

$$ $\displaystyle\forall q,\ ~{}~{}~{}\mathbb{E}_{q(\bm{x})}\big{[}\mathbb{E}_{P(X,y)}[F(\bm{x};(X,y))]\big{]}\leq\mathbb{E}_{q(\bm{x})}\big{[}F(\bm{x};\mathcal{S})\big{]}+\sqrt{\frac{\mathrm{KL}(q||p)+\log(\frac{M}{\kappa})}{2(M-1)}}.$ (117) $$

Set the posterior distribution as $q=p_{\bm{\theta}}:=\mathcal{N}(\bm{\mu},\bm{\Sigma})$. Denote the set $\mathcal{M}(\bm{\theta})=\{\bm{\delta}\mid\mathrm{KL}(p_{\bm{\theta}+\bm{\delta}}\|p_{\bm{\theta}})+\mathrm{KL}(p_{\bm{\delta}}\|p_{\bm{\theta}+\bm{\delta}})\leq\rho^{2}\}$. We can choose $p\in\mathcal{M}(\bm{\theta})$. Then, we know that for the prior distribution $p$, the following inequality holds with a probability at least $1-\kappa$,

$$ $\displaystyle\mathbb{E}_{p_{\bm{\theta}}(\bm{x})}\big{[}\mathbb{E}_{P(X,y)}[F(\bm{x};(X,y))]\big{]}$ $\displaystyle\leq\mathbb{E}_{p_{\bm{\theta}}(\bm{x})}\big{[}F(\bm{x};\mathcal{S})\big{]}+\sqrt{\frac{\mathrm{KL}(p_{\bm{\theta}}||p)+\log(\frac{M}{\kappa})}{2(M-1)}}.$ (118) $$

Note that for any density $p,q$, we have $\mathrm{KL}(p||q)\geq 0$. Thus, we know $\mathcal{M}(\bm{\theta})\subset\mathcal{C}(\bm{\theta})$.
It follows that

$$ $\displaystyle\mathbb{E}_{p_{\bm{\theta}}(\bm{x})}\big{[}\mathbb{E}_{P(X,y)}[F(\bm{x};(X,y))]\big{]}$ $\displaystyle\leq\mathbb{E}_{p_{\bm{\theta}}(\bm{x})\in\mathcal{C}(\bm{\theta})}\big{[}F(\bm{x};\mathcal{S})\big{]}+\sqrt{\frac{\mathrm{KL}(p_{\bm{\theta}}||p)+\log(\frac{M}{\kappa})}{2(M-1)}}$ (119) $\displaystyle\leq\max_{\bm{\delta}\in\mathcal{C}(\bm{\theta})}\mathbb{E}_{p_{\bm{\theta}+\bm{\delta}}}\big{[}F(\bm{x};\mathcal{S})\big{]}+\max_{\bm{\delta}\in\mathcal{M}(\bm{\theta})}\sqrt{\frac{\mathrm{KL}(p_{\bm{\theta}}||p_{\theta+\delta})+\log(\frac{M}{\kappa})}{2(M-1)}}$ (120) $\displaystyle\leq\max_{\bm{\delta}\in\mathcal{C}(\bm{\theta})}\mathbb{E}_{p_{\bm{\theta}+\bm{\delta}}}\big{[}F(\bm{x};\mathcal{S})\big{]}+\sqrt{\frac{\rho^{2}+\log(\frac{M}{\kappa})}{2(M-1)}}.$ (121) $$

Note that $F(\bm{x};(X,y))$ is convex function w.r.t. $\bm{x}$, we know that $\mathbb{E}_{P(X,y)}[F(\bm{x};(X,y))]$ is a convex function w.r.t. $\bm{x}$. It follows that

$$ $\displaystyle\mathbb{E}_{P(X,y)}[F(\bm{\mu};(X,y))]$ $\displaystyle=\mathbb{E}_{P(X,y)}[F(\mathbb{E}_{p_{\bm{\theta}}(\bm{x})}[\bm{x}];(X,y))]$ (122) $\displaystyle\leq\mathbb{E}_{p_{\bm{\theta}}(\bm{x})}\big{[}\mathbb{E}_{P(X,y)}[F(\bm{x};(X,y))]\big{]}.$ (123) $$

Finally, we know that with a probability of at least $1-\kappa$, the following inequality holds.

$$ $\displaystyle\mathbb{E}_{P(X,y)}[F(\bm{\mu};\bm{z},y)]\leq\max_{\bm{\delta}\in\mathcal{C}(\bm{\theta})}\mathbb{E}_{p_{\bm{\theta}+\bm{\delta}}}\big{[}F(\bm{x};\mathcal{S})\big{]}+\sqrt{\frac{\rho^{2}+\log(\frac{M}{\kappa})}{2(M-1)}}.$ (124) $$

## Appendix E Proof of Technical Lemmas

In this section, we provide the proof of lemmas in Appendix [C](#A3). Note that Lemma [C.1](#A3.Thmtheorem1) and Lemma [C.2](#A3.Thmtheorem2) can be directly obtained by Lemma B.1. and Lemma B.2. in [], respectively.

### E.1 Proof of Lemma C.3

(a): Since we have $\bm{\Sigma}_{t+1}^{-1}=\bm{\Sigma}_{t}^{-1}+2\beta_{t}\bm{G}_{t}$. We can obtain that

$$ $\displaystyle\bm{\Sigma}_{t}^{-1}+2b\beta_{t}\bm{I}\succeq\bm{\Sigma}_{t+1}^{-1}\succeq\bm{\Sigma}_{t}^{-1}+2\xi\beta_{t}\bm{I}.$ (125) $$

Summing up it over $t=0,\dots,T-1$, we have

$$ $\displaystyle\bm{\Sigma}_{0}^{-1}+2b\sum_{t=1}^{T}\beta_{t}\bm{I}\succeq\bm{\Sigma}_{T}^{-1}\succeq\bm{\Sigma}_{0}^{-1}+2\xi\sum_{t=1}^{T}\beta_{t}\bm{I}.$ (126) $$

Therefore, we have

$$ $\displaystyle\frac{1}{2b\sum_{t=1}^{T}\beta_{t}\bm{I}+\bm{\Sigma}_{0}^{-1}}\preceq\bm{\Sigma}_{T}\preceq\frac{1}{2\xi\sum_{t=1}^{T}\beta_{t}\bm{I}+\bm{\Sigma}_{0}^{-1}}.$ (127) $$

(b): We have

$$ $\displaystyle\|\bm{\Sigma}_{t}\|_{\mathrm{F}}\leq\left\|\frac{1}{2\xi\sum_{k=1}^{t}\beta_{k}\bm{I}+\bm{\Sigma}_{0}^{-1}}\right\|_{\mathrm{F}}\leq\left\|\frac{1}{2\xi\sum_{k=1}^{t}\beta_{k}\bm{I}}\right\|_{\mathrm{F}}=\frac{\sqrt{d}}{2\xi\sum_{k=1}^{t}\beta_{k}}.$ (128) $$

### E.2 Proof of Lemma C.4

In Algorithm [1](#alg1) and Algorithm [2](#alg2), for given $\bm{\theta}_{t}$, the perturbation $\bm{\delta}_{t}$ satisfies

$$ $\bm{\delta}_{\bm{\mu}_{t}}=\frac{1}{\lambda}\bm{\Sigma}_{t}\bm{g}_{t}^{\prime}=\frac{\rho\bm{\Sigma}_{t}\bm{g}_{t}^{\prime}}{\sqrt{\|\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}\|_{\mathrm{F}}^{2}+0.5\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}^{\prime}\|^{2}_{2}}}\leq\frac{\rho\sqrt{2}\bm{\Sigma}^{\frac{1}{2}}_{t}\bm{\Sigma}^{\frac{1}{2}}_{t}\bm{g}_{t}^{\prime}}{\sqrt{\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}^{\prime}\|^{2}_{2}}}.$ (129) $$

Therefore, we have $\|\bm{\delta}_{\bm{\mu}_{t}}\|\leq\rho\sqrt{2}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\|_{\mathrm{F}}$. For $\bm{\delta}_{\bm{\Sigma}_{t}}$, we have

$$ $\bm{\delta}_{\bm{\Sigma}_{t}}=\frac{2\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}}{\lambda\bm{\Sigma}_{t}^{-1}-2\bm{G}_{t}^{\prime}}=\bm{\Sigma}_{t}\frac{\frac{2}{\lambda}\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}}{\bm{I}-\frac{2}{\lambda}\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}}.$ (130) $$

Note that $\frac{2}{\lambda}\leq\frac{2\rho}{\|\bm{\Sigma}\bm{G}_{t}^{\prime}\|_{\mathrm{F}}}$. Therefore, we can obtain that

$$ $\displaystyle\|\bm{\delta}_{\bm{\Sigma}_{t}}\|_{\mathrm{F}}\leq\|\bm{\Sigma}_{t}\|_{\mathrm{F}}\frac{2\rho}{\|\bm{I}-\frac{2}{\lambda}\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}\|_{\mathrm{F}}}.$ (131) $$

If $\rho<\frac{\sqrt{d}}{2}$, there exist a constant $r>0$, $\|\bm{I}-\frac{2}{\lambda}\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}\|_{\mathrm{F}}>\frac{1}{r}$ holds. Therefore we have $\|\bm{\delta}_{\bm{\Sigma}_{t}}\|_{\mathrm{F}}\leq 2\rho r\|\bm{\Sigma}_{t}\|_{\mathrm{F}}$.

### E.3 Proof of Lemma C.5

We have

$$ $\displaystyle\mathbb{E}_{\bm{z}}[\bm{g}_{t}^{\prime}]$ $\displaystyle=\mathbb{E}_{\bm{z}}[\bm{\Sigma}_{t}^{-\frac{1}{2}}{\bm{z}}F(\bm{\mu}_{t}+\bm{\Sigma}_{t}^{\frac{1}{2}}{\bm{z}})]-\mathbb{E}_{\bm{z}}[\bm{\Sigma}_{t}^{-\frac{1}{2}}{\bm{z}}F(\bm{\mu}_{t})]$ (132) $\displaystyle=\mathbb{E}_{\bm{z}}[\bm{\Sigma}_{t}^{-\frac{1}{2}}{\bm{z}}F(\bm{\mu}_{t}+\bm{\Sigma}_{t}^{\frac{1}{2}}{\bm{z}})]$ (133) $\displaystyle=\mathbb{E}_{\bm{x}\sim\mathcal{N}(\bm{\mu}_{t},\bm{\Sigma}_{t})}[\bm{\Sigma}_{t}^{-1}(\bm{x}-\bm{\mu}_{t})F(\bm{x})]$ (134) $\displaystyle=\nabla_{\bm{\mu}}\mathbb{E}_{p_{\bm{\theta}_{t}}}[F(\bm{x})],$ (135) $$

where we reach the conclusion.

### E.4 Proof of Lemma C.6

Denote the diagonal elements of $\bm{\Sigma}$ and $\widehat{\bm{\Sigma}}$ by $\bm{\sigma}$ and $\widehat{\bm{\sigma}}$, respectively. Then we have

$$ $\displaystyle\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t,j}\|^{2}=\|\bm{\sigma}_{t}^{\frac{1}{2}}\odot\widehat{\bm{\sigma}}_{t}^{-\frac{1}{2}}\odot\bm{z}_{j}(F(\widehat{\bm{\mu}}_{t}+\widehat{\bm{\sigma}}_{t}^{\frac{1}{2}}\odot\bm{z}_{j})-F(\widehat{\bm{\mu}}_{t}))\|^{2}.$ (136) $$

Note that $\widehat{\bm{\Sigma}}_{t+1}^{-1}=\bm{\Sigma}_{t}^{-1}-2\lambda\bm{G}_{t}^{\prime}$, we have $\bm{\sigma}_{t}\leq\widehat{\bm{\sigma}}_{t}$.
Since we have

$$ $\displaystyle\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}_{2}=\mathbb{E}\left\|\frac{1}{N}\sum_{j=1}^{N}\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t,j}\right\|^{2}\leq\frac{1}{N}\sum_{j=1}^{N}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t,j}\|^{2}.$ (137) $$

It follows that

$$ $\displaystyle\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}_{2}$ $\displaystyle\leq\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\|\widehat{\bm{\sigma}}_{t}^{\frac{1}{2}}\odot\widehat{\bm{\sigma}}_{t}^{-\frac{1}{2}}\odot\bm{z}_{j}(F(\widehat{\bm{\mu}}_{t}+\widehat{\bm{\sigma}}_{t}^{\frac{1}{2}}\odot\bm{z}_{j})-F(\widehat{\bm{\mu}}_{t}))\|^{2}$ (138) $\displaystyle\leq\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[\|\bm{z}_{j}\|^{2}(F(\widehat{\bm{\mu}}_{t}+\widehat{\bm{\sigma}}_{t}^{\frac{1}{2}}\odot\bm{z}_{j})-F(\widehat{\bm{\mu}}_{t}))^{2}\right]$ (139) $\displaystyle\leq\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[H^{2}\|\bm{z}_{j}\|^{4}\|\widehat{\bm{\sigma}}_{t}\|_{\infty}\right].$ (140) $$

Note that

$$ $\displaystyle\widehat{\bm{\Sigma}}_{t}=\bm{\Sigma}_{t}+\bm{\delta}_{\bm{\Sigma}_{t}}=\bm{\Sigma}_{t}+\bm{\Sigma}_{t}\frac{\frac{2}{\lambda}\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}}{\bm{I}-\frac{2}{\lambda}\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}}.$ (141) $$

Note that $\frac{2}{\lambda}\leq\frac{2\rho}{\|\bm{\Sigma}\bm{G}_{t}^{\prime}\|_{\mathrm{F}}}$. Then if $\rho<\frac{\sqrt{d}}{2}$, there exist a constant $r^{\prime}>0$, the inequality $\widehat{\bm{\sigma}}_{t}\leq\bm{\sigma}_{t}+2\rho r^{\prime}\bm{\sigma}_{t}$ holds. Therefore, we have

$$ $\displaystyle\|\widehat{\bm{\sigma}}_{t}\|_{\infty}\leq(1+2\rho r^{\prime})\|\bm{\sigma}_{t}\|_{\infty}\leq\frac{1+2\rho r^{\prime}}{\|\bm{\sigma}_{0}^{-1}\|_{min}+2(\sum_{k=1}^{t}\beta_{k})\xi},$ (142) $$

where $\|\cdot\|_{min}$ denotes the minimum element in the input. Noticed that $\mathbb{E}_{\bm{z}}[\|{\bm{z}}\|^{2}]\leq d+4$ as shown in [], we obtain

$$ $\displaystyle\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}_{2}$ $\displaystyle\leq\frac{H^{2}(1+2\rho r^{\prime})(d+4)^{2}}{2\xi(\sum_{k=1}^{t}\beta_{k})},$ (143) $$

where we reach the conclusion.

### E.5 Proof of Lemma C.7

We have

$$ $\displaystyle\mathbb{E}[\bm{g}_{t}^{\prime}]$ $\displaystyle=\mathbb{E}[\bm{\Sigma}_{t}^{-\frac{1}{2}}{\bm{z}}F(\bm{\mu}_{t}+\bm{\Sigma}_{t}^{\frac{1}{2}}{\bm{z}};\mathcal{B})]-\mathbb{E}[\bm{\Sigma}_{t}^{-\frac{1}{2}}{\bm{z}}F(\bm{\mu}_{t};\mathcal{B})]$ (144) $\displaystyle=\mathbb{E}_{\bm{z}}[\bm{\Sigma}_{t}^{-\frac{1}{2}}{\bm{z}}F(\bm{\mu}_{t}+\bm{\Sigma}_{t}^{\frac{1}{2}}{\bm{z}};\mathcal{D})]-\mathbb{E}_{\bm{z}}[\bm{\Sigma}_{t}^{-\frac{1}{2}}{\bm{z}}F(\bm{\mu}_{t};\mathcal{D})]$ (145) $\displaystyle=\mathbb{E}_{\bm{z}}[\bm{\Sigma}_{t}^{-\frac{1}{2}}{\bm{z}}F(\bm{\mu}_{t}+\bm{\Sigma}_{t}^{\frac{1}{2}}{\bm{z}})]$ (146) $\displaystyle=\mathbb{E}_{\bm{x}\sim\mathcal{N}(\bm{\mu}_{t},\bm{\Sigma}_{t})}[\bm{\Sigma}_{t}^{-1}(\bm{x}-\bm{\mu}_{t})F(\bm{x})]$ (147) $\displaystyle=\nabla_{\bm{\mu}}\mathbb{E}_{p_{\bm{\theta}_{t}}}[F(\bm{x})],$ (148) $$

where the second equality is due to Proposition [4.5](#S4.Thmtheorem5).

### E.6 Proof of Lemma C.8

Denote the diagonal elements of $\bm{\Sigma}$ and $\widehat{\bm{\Sigma}}$ by $\bm{\sigma}$ and $\widehat{\bm{\sigma}}$, respectively. Then we have

$$ $\displaystyle\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t,j}\|^{2}=\|\bm{\sigma}_{t}^{\frac{1}{2}}\odot\widehat{\bm{\sigma}}_{t}^{-\frac{1}{2}}\odot\bm{z}_{j}(F(\widehat{\bm{\mu}}_{t}+\widehat{\bm{\sigma}}_{t}^{\frac{1}{2}}\odot\bm{z}_{j};\mathcal{B})-F(\widehat{\bm{\mu}}_{t});\mathcal{B})\|^{2}$ (149) $$

Note that $\widehat{\bm{\Sigma}}_{t+1}^{-1}=\bm{\Sigma}_{t}^{-1}-2\lambda\bm{G}_{t}^{\prime}$, we have $\bm{\sigma}_{t}\leq\widehat{\bm{\sigma}}_{t}$.
Since we have

$$ $\displaystyle\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}_{2}=\mathbb{E}\left\|\frac{1}{N}\sum_{j=1}^{N}\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t,j}\right\|^{2}\leq\frac{1}{N}\sum_{j=1}^{N}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t,j}\|^{2}.$ (150) $$

It follows that

$$ $\displaystyle\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}_{2}$ $\displaystyle\leq\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\|\widehat{\bm{\sigma}}_{t}^{\frac{1}{2}}\odot\widehat{\bm{\sigma}}_{t}^{-\frac{1}{2}}\odot\bm{z}_{j}(F(\widehat{\bm{\mu}}_{t}+\widehat{\bm{\sigma}}_{t}^{\frac{1}{2}}\odot\bm{z}_{j};\mathcal{B})-F(\widehat{\bm{\mu}}_{t};\mathcal{B}))\|^{2}$ (151) $\displaystyle\leq\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[\|\bm{z}_{j}\|^{2}(F(\widehat{\bm{\mu}}_{t}+\widehat{\bm{\sigma}}_{t}^{\frac{1}{2}}\odot\bm{z}_{j};\mathcal{B})-F(\widehat{\bm{\mu}}_{t};\mathcal{B}))^{2}\right]$ (152) $\displaystyle\leq\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[\|\bm{z}_{j}\|^{2}(\frac{1}{3}H^{2}\|\widehat{\bm{\sigma}}_{t}^{\frac{1}{2}}\odot\bm{z}_{j}\|^{2}+\frac{2}{3}\varepsilon^{2})\right],$ (153) $$

where the third inequality is due to Proposition [4.5](#S4.Thmtheorem5).
Since $\widehat{\bm{\Sigma}}_{t}=\bm{\Sigma}_{t}+\bm{\Sigma}_{t}\frac{\frac{2}{\lambda}\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}}{\bm{I}-\frac{2}{\lambda}\bm{\Sigma}_{t}\bm{G}_{t}^{\prime}}$, and $\frac{2}{\lambda}\leq\frac{2\rho}{\|\bm{\Sigma}\bm{G}_{t}^{\prime}\|_{\mathrm{F}}}$. Then if $\rho<\frac{\sqrt{d}}{2}$, there exist a constant $r^{\prime}>0$, the inequality $\widehat{\bm{\sigma}}_{t}\leq\bm{\sigma}_{t}+2\rho r^{\prime}\bm{\sigma}_{t}$ holds. Therefore, we have

$$ $\displaystyle\|\widehat{\bm{\sigma}}_{t}\|_{\infty}\leq(1+2\rho r^{\prime})\|\bm{\sigma}_{t}\|_{\infty}\leq\frac{1+2\rho r^{\prime}}{\|\bm{\sigma}_{0}^{-1}\|_{min}+2(\sum_{k=1}^{t}\beta_{k})\xi},$ (154) $$

where $\|\cdot\|_{min}$ denotes the minimum element in the input. Noticed that $\mathbb{E}_{\bm{z}}[\|{\bm{z}}\|^{2}]\leq d+4$ as shown in [], we obtain

$$ $\displaystyle\mathbb{E}\|\bm{\Sigma}_{t}^{\frac{1}{2}}\bm{g}_{t}\|^{2}_{2}$ $\displaystyle\leq\frac{1}{N}\sum_{j=1}^{N}\mathbb{E}\left[\frac{(1+2\rho r^{\prime})H^{2}}{3}\big{\|}\bm{\sigma}_{t}\big{\|}_{\infty}\times\big{\|}\bm{z}_{j}\big{\|}^{4}+\frac{2\varepsilon^{2}}{3}\|\bm{z}_{j}\|^{2}\right]$ (155) $\displaystyle\leq\frac{H^{2}(1+2\rho r^{\prime})(d+4)^{2}}{6\xi(\sum_{k=1}^{t}\beta_{k})}+\frac{2(d+4)\varepsilon^{2}}{3},$ (156) $$

where we reach the conclusion.

## Appendix F Additional Materials for Section 6

In this section, we provide additional experiments and more implementation details of Section [6](#S6).

### F.1 Synthetic Problems

The four numerical benchmark test functions employed in Section [6.1](#S6.SS1) are listed as follows:

$$ $\displaystyle F(\bm{x})$ $\displaystyle=\sum\nolimits_{i=1}^{d}10^{\frac{2(i-1)}{d-1}}\bm{x}_{i}^{2},$ (157) $\displaystyle F(\bm{x})$ $\displaystyle=\sum\nolimits_{i=1}^{d}10^{\frac{2(i-1)}{d-1}}|\bm{x}_{i}|^{\frac{1}{2}},$ (158) $\displaystyle F(\bm{x})$ $\displaystyle=\sqrt{\sum\nolimits_{i=1}^{d}|x_{i}|^{2+4\frac{i-1}{d-1}}},$ (159) $\displaystyle\begin{split}F(\bm{x})&=\sin^{2}(\pi\omega_{1})+\sum_{i=1}^{d-1}(\omega_{i}-1)^{2}(1+10\sin^{2}(\omega_{i}\pi+1))+(\omega_{d}-1)^{2}(1+\sin^{2}(2\pi\omega_{d})),\\ &\ ~{}~{}~{}~{}~{}~{}\text{where}\ ~{}\omega_{i}=1+\frac{x_{i}-1}{4},i\in\{1,\dots,d\}.\end{split}$ (160) $$

Test functions ([157](#A6.E157))-([160](#A6.E160)) are called the ellipsoid function, $l_{\frac{1}{2}}$-ellipsoid function, different powers function and Levy function, respectively.

#### Implementation Details.

For all the methods, we initialize $\bm{\mu}_{0}$ from the uniform distribution $\text{Uni}[0,1]$, and set $\bm{\Sigma}_{0}=\bm{I}$. For the INGO, BES, and SABO methods, we use a fixed step size of $\beta=0.1$. According to our assumption in Theorem [4.3](#S4.Thmtheorem3), we set $\rho=100/\sqrt{T+1}$ for the proposed SABO method. We set the spacing $c=1$ for the BES method and employ the default hyperparameter setting from [] for the MMES method. We then assess these methods using varying dimensions, i.e., $d\in\{200,500,1000\}$. For $d=200$, we assess these methods using varying sample sizes, i.e., $N\in\{10,50,100\}$. The mean value of $\mathcal{E}$ over 3 independent runs is reported.

#### Results.

Figure [2](#A6.F2) and [3](#A6.F3) show the results on four test functions with problem dimensions $d=200$ and $d=1000$, respectively. Figure [4](#A6.F4) and [5](#A6.F5) show the results on $200$-dimensional test functions with sample size $N=10$ and $N=100$, respectively. Combining these results with the result from Figure [1](#S6.F1), we observe consistent performance for the proposed ASBO method. It achieves a similar convergence result to the INGO method, as they have the same theoretical convergent rate. In some cases, i.e., Figure [3](#A6.F3) (d) and Figure [4](#A6.F4) (d), it converges slightly faster than INGO. The CMA-ES method and MMES method both work for ellipsoid and different powers functions, but fail in $l_{\frac{1}{2}}$-ellipsoid and Levy functions. In most cases, they converge slower than INGO and SABO. With a large sample size, i.e., $N=100$, the CMA-ES method can maintain a fast converge rate according to Figure [4](#A6.F4) (a). The BES method fails to achieve high precision in all test functions. It diverges in ellipsoid, $l_{\frac{1}{2}}$-ellipsoid, Levy functions, and only achieves a low precision in the different power functions. These results demonstrate the superiority of the SABO method in optimizing high-dimensional problems, and verify our theoretical convergence results in Section [4](#S4).

Figure: (a) Ellipsoid.
Refer to caption: /html/2410.12457/assets/x5.png

Figure: (a) Ellipsoid.
Refer to caption: /html/2410.12457/assets/x9.png

Figure: (a) Ellipsoid.
Refer to caption: /html/2410.12457/assets/x13.png

Figure: (a) Ellipsoid.
Refer to caption: /html/2410.12457/assets/x17.png

**Table 2: Performance ($\%$) on SST-2, AG’s News, MRPC, RTE, SNLI, and Yelp P. datasets. We report the mean and standard deviation over $3$ random seeds. The best result across all groups is highlighted in bold and the best result in each group is marked with underlined.**
| Methods | SST-2 | AG’s News | MRPC | RTE | SNLI | Yelp P. |
| --- | --- | --- | --- | --- | --- | --- |
| Zero-shot | $79.82$ | $76.96$ | $67.40$ | $51.62$ | $38.82$ | $89.64$ |
| Dimension $d=200$ |  |  |  |  |  |  |
| GIBO | $83.53_{\pm 0.15}$ | $75.79_{\pm 0.08}$ | $79.21_{\pm 0.09}$ | $53.07_{\pm 0.29}$ | $38.73_{\pm 0.09}$ | $89.63_{\pm 0.03}$ |
| TurBo | $83.30_{\pm 0.19}$ | $79.01_{\pm 1.68}$ | $69.59_{\pm 8.51}$ | $46.57_{\pm 2.30}$ | $40.27_{\pm 0.69}$ | $90.16_{\pm 0.19}$ |
| SABO | $\underline{87.88}_{\pm 0.53}$ | $\underline{82.22}_{\pm 0.41}$ | $\underline{79.35}_{\pm 0.12}$ | $\underline{\textbf{53.67}}_{\pm 0.17}$ | $\underline{40.72}_{\pm 0.15}$ | $\underline{91.50}_{\pm 0.13}$ |
| Dimension $d=500$ |  |  |  |  |  |  |
| GIBO | $83.49_{\pm 0.09}$ | $75.70_{\pm 0.05}$ | $79.03_{\pm 0.08}$ | $52.95_{\pm 0.17}$ | $38.71_{\pm 0.16}$ | $89.65_{\pm 0.02}$ |
| TurBo | $84.52_{\pm 0.65}$ | $80.03_{\pm 1.97}$ | $75.30_{\pm 2.34}$ | $48.01_{\pm 0.59}$ | $38.82_{\pm 0.34}$ | $90.20_{\pm 0.45}$ |
| SABO | $\underline{87.31}_{\pm 0.38}$ | $\underline{82.65}_{\pm 0.59}$ | $\underline{79.62}_{\pm 0.07}$ | $\underline{53.55}_{\pm 0.17}$ | $\underline{\textbf{42.29}}_{\pm 2.48}$ | $\underline{91.83}_{\pm 0.16}$ |
| Dimension $d=1000$ |  |  |  |  |  |  |
| GIBO | $83.45_{\pm 0.11}$ | $75.67_{\pm 0.10}$ | $79.15_{\pm 0.0}$ | $52.95_{\pm 0.17}$ | $38.87_{\pm 0.18}$ | $89.65_{\pm 0.04}$ |
| TurBo | $85.90_{\pm 0.95}$ | $82.36_{\pm 0.21}$ | $77.30_{\pm 0.86}$ | $50.30_{\pm 1.11}$ | $39.87_{\pm 1.07}$ | $90.14_{\pm 0.20}$ |
| SABO | $\underline{\textbf{87.96}}_{\pm 0.83}$ | $\underline{\textbf{82.77}}_{\pm 0.41}$ | $\underline{\textbf{79.68}}_{\pm 0.23}$ | $\underline{53.31}_{\pm 0.17}$ | $\underline{40.32}_{\pm 0.27}$ | $\underline{\textbf{91.96}}_{\pm 0.41}$ |

### F.2 Additional Results on Black-box Prompt Fine-tuning

We conduct an additional experiment on the black-box prompt fine-tuning task to compare the proposed SABO method with high-dimensional BO methods, i.e., TuRBO [] and GIBO [] methods discussed in Section [5](#S5). We employ the default setting of TuRBO-1 from [] for the TuRBO method, and the default setting from [] for the GIBO method. The results on six benchmark datasets are reported in Table [2](#A6.T2). According to the results, the SABO method consistently outperforms these two baselines, highlighting its effectiveness.

### F.3 Synthetic Image Classification Problem

We conduct experiments on a synthetic image classification problem to empirically evaluate the performance of the proposed mini-batch SABO method.
Specifically, we apply black-box optimization methods to train a model to classify the images accurately. Moreover, following [], we conduct additional experiments on the noisy label setting. Particularly, we train the model on a corrupted version of the Fashion-MNIST dataset, where some of its training labels are randomly flipped according to different noise rates, while the testing set is clean.
To construct problems with different dimensions, we first adopt UMAP [] to reduce the dimension of Fashion-MNIST to $8$ while preserving the class discriminability [], and then employ a fixed randomly initialized matrix $\bm{P}\in\mathbb{R}^{8\times\widetilde{d}}$ to project the extracted features to a $\widetilde{d}$-dimensional space. After preprocessing, a linear layer is used as a classifier. Therefore, the total number of the trainable parameters is $d=10\times\widetilde{d}$.

#### Datasets.

Fashion-MNIST [] is a image classifications dataset. It contains $60,000$ training samples and $10,000$ test samples, each representing a $28\times 28$-pixel grayscale image of fashion items from $10$ different classes.

#### Implementation Details.

For a fair comparison, we set the same population size, number of batch samples, and the initialization for CMA-ES [], INGO [], BES [], and SABO. The population size and the number of batch samples are set to $N=100$ and $M=2048$, respectively. The Gaussian distributions are initialized with $\bm{\mu}_{0}=\bm{0}$ and $\bm{\Sigma}_{0}=0.5\bm{I}$.
For BES, INGO and SABO methods, we search the learning rate $\beta$ over $\{0.1,0.5,1,5\}$. Moreover, we perform grid-search on $\rho$ over $\{100,500,1000,5000\}$ for SABO.
We employ the default hyperparameter setting from [] for the MMES method.
For the BES method, we perform grid-search on the spacing $c$ over $\{0.1,1,10\}$. The cross-entropy loss is used as the training objective.
All experiments are repeatedly run with three independent seeds and the mean and standard deviation are reported.

#### Results.

Table [3](#A6.T3) shows experimental results on Fashion-MNIST dataset with different $\widetilde{d}$ and different noise rates.
We can see that the SABO method consistently outperforms all baselines across different noise rates and dimensions, highlighting its effectiveness. These results show that the proposed SABO method can achieve good robustness performance and demonstrate the effectiveness of the proposed SABO method in improving model generalization performance.

**Table 3: Performance ($\%$) on Fashion-MNIST dataset with noise labels. The best result across all groups is highlighted in bold and the best result in each group is marked with underlined.**
| Methods | noise=$0\%$ | noise=$20\%$ | noise=$40\%$ | noise=$60\%$ | noise=$80\%$ |
| --- | --- | --- | --- | --- | --- |
| Dimension $d=100$ |  |  |  |  |  |
| CMA-ES | $76.68\pm_{0.44}$ | $77.25\pm_{1.26}$ | $75.71\pm_{1.15}$ | $70.81\pm_{1.01}$ | $57.86\pm_{3.65}$ |
| MMES | $76.20\pm_{1.45}$ | $74.68\pm_{0.33}$ | $74.46\pm_{0.32}$ | $71.12\pm_{2.24}$ | $65.92\pm_{1.30}$ |
| BES | $76.26\pm_{0.21}$ | $73.60\pm_{0.50}$ | $67.60\pm_{0.72}$ | $55.14\pm_{6.26}$ | $45.65\pm_{4.35}$ |
| INGO | $76.68\pm_{0.90}$ | $75.16\pm_{0.55}$ | $71.26\pm_{0.74}$ | $61.70\pm_{3.89}$ | $46.75\pm_{1.35}$ |
| SABO | $\underline{\textbf{78.10}}\pm_{0.90}$ | $\underline{77.41}\pm_{0.68}$ | $\underline{\textbf{77.10}}\pm_{1.69}$ | $\underline{73.41}\pm_{0.33}$ | $\underline{\textbf{66.45}}\pm_{0.43}$ |
| Dimension $d=1000$ |  |  |  |  |  |
| CMA-ES | $76.75\pm_{0.54}$ | $74.85\pm_{0.06}$ | $73.21\pm_{0.52}$ | $67.32\pm_{0.84}$ | $50.00\pm_{0.98}$ |
| MMES | $75.85\pm_{0.43}$ | $72.39\pm_{3.59}$ | $73.73\pm_{0.50}$ | $72.51\pm_{1.47}$ | $56.60\pm_{2.72}$ |
| BES | $76.18\pm_{0.63}$ | $72.45\pm_{0.44}$ | $66.64\pm_{2.09}$ | $58.74\pm_{1.50}$ | $46.99\pm_{4.30}$ |
| INGO | $74.73\pm_{0.62}$ | $74.87\pm_{1.72}$ | $69.75\pm_{0.73}$ | $63.18\pm_{4.04}$ | $41.80\pm_{1.36}$ |
| SABO | $\underline{76.87}\pm_{1.41}$ | $\underline{\textbf{77.73}}\pm_{0.85}$ | $\underline{74.97}\pm_{1.12}$ | $\underline{\textbf{73.81}}\pm_{1.26}$ | $\underline{61.96}\pm_{0.08}$ |

## Appendix G Discussion with Offline Black-box Optimization

The typical black-box optimization problem we studied in this work can also be called the online black-box optimization problem. Since we have access to the objective function $F$, the problem can be solved in an online iterative manner, where in each iteration the solver proposes new $\bm{x}$ and queries the objective function for feedback in order to inform better solution proposals at the next iteration.

The offline black-box optimization [, ] is different from the research line of standard online black-box optimization. In offline black-box optimization, access to the true objective $F$ is not available. Instead, the offline black-box optimization algorithm is provided access to a static dataset $\mathcal{D}=\{\bm{x_{i}},F(\bm{x}_{i})\}$ of the variable $\bm{x_{i}}$ and its corresponding objective value $F(\bm{x}_{i})$. Therefore, the basic settings of offline black-box optimization and online black-box optimization are different.

Moreover, the challenges of offline black-box optimization and online black-box optimization are also different. Offline black-box optimization focuses on producing query candidates by a surrogate model trained with a prior static dataset [, ]. The goal is to produce a good query set based on the fixed model at one time without considering query feedback update, exploration and exploitation balance in the long term. Along this line, [] proposed a model-based offline optimization by training a conditional generative model that conditions the objective value.
In [], the authors formulated the problem as a non-zero-sum game and proposed an alternating ascent-descent algorithm for model-based offline optimization.
[] proposed the conservative objective models, which presents a technique similar to adversarial training that avoids overestimation of out-of-distribution inputs. In contrast to offline black-box optimization, standard online black-box optimization needs to balance exploration and exploitation, which focuses on long-term convergence performance.
As a result, offline black-box optimization is not suitable for our online prompt fine-tuning tasks.

## Peer Reviews

confidence: 3
contribution: 2
initial_rating: 5
presentation: 3
questions: See Weakness
soundness: 2
strengths: 1. The paper is clearly written, well organized and easy to follow.
2. I think the whole topic is intriguing and worthy to probe, and the authors give some interesting insights.
3. It is valuable to find that the proposed method could be effective in prompt tuning.
summary: The authors strive to introduce SAM to the black-box optimization algorithm, and propose a new algorithm called Sharpness-Aware Black-box Optimization (SABO). SABO reparameterizes the objective function using its expectation over a Gaussian distribution and iteratively updates the parameterized distribution based on approximated stochastic gradients of the maximum objective value within a small neighborhood. The authors provide theoretical guarantees for SABO's convergence rate and generalization bounds. Also, they give empirical experiments to validate its effectiveness in improving model generalization performance on black-box prompt fine-tuning tasks and synthetic numerical problems.
weaknesses: 1. My primary concern is the motivation regarding the introduction of SAM in black-box optimization. In other words, can we safely use the proposed method in black-box optimization? The implicit premise of SAM relies on the smoothness of the mapping function. However, black-box optimization may deal with highly complex, non-smooth, or potentially noisy objective landscapes, scenarios much more complicated than optimizing neural networks. In my opinion, it is crucial for the authors to clarify the scope of which minimizing the maximization of the neighborhood region is appropriate for black-box optimization.

2. I noticed that in the algorithm, the authors use two different samples to separately calculate the ascent and descent gradients. Could the authors give some explanations here given that we always keep the same batch in SAM?

3. The essence of SAM involves two main steps: first, a maximization step using gradient ascent based on the current parameters with a step size of rho/||g||; and second, a minimization step using gradient descent based on the maximized parameters with a step size equal to the learning rate. What are the results of applying this procedure in the current black-box optimization?

4. The convergence analysis mostly relies on the assumption of convexity of F(x), which I think could be less useful.

5. Some minor mistakes: (1) Line 121, SAM is proposed by [Foretet al., 2021] not [Kwon et al., 2021].
(2) Eq 6. nabla_F(x + epsilon_t) should be rewrite as nabla_F(x)|x = epsilon_t

---

confidence: 2
contribution: 3
initial_rating: 10
presentation: 4
questions: - in the setup on line 177, you claim that for computational consideration you can assume the search distribution to be Gaussian. Is that a fair simplification? How sensitive is SABO to this choice of distribution and would it make more sense to sample the local parameter space in some way that captures the nature of the posterior rather than with a naive Gaussian?
 - related to above (line 180), is it fair to assume that the covariance matrix is diagonal. This seems optimistic and is essentially assuming no correlation between parameters, correct? And then your perturbations $\delta$ also have diagonal covariance. I understand that, as you point out in line 195, removing this assumption would make the problem computationally substantially more difficult.  Theoretically, how much of a sacrifice does this assumption make? and are there any other means you can use to not lose this amount of information?
soundness: 3
strengths: This is an interesting paper that combines ideas from sharpness aware minimization (SAM) with black box optimization. SAM has already proven to be a very useful technique to improve generalization in training ML models and the authors present an intuitive motivation for applying these ideas to BBO. SABO is a creative combination of SAM techniques and black box optimization methods. 

This paper has a strong theoretical analysis of SABO, including convergence rate and generalization error bounds. These theoretical results are then empirically validated on synthetic and real-worl black-box prompt fine-tuning models. 

The paper is very well-written and logically structured. The authors do a good job motivating their approach and everything is presented in a clear understandable manner.
summary: This paper combines ideas from sharpness aware minimization (SAM) in machine learning with black-box optimization problems to propose a new algorithm called Sharpness-Aware Black-box Optimization (SABO). The goal of SABO is to improve the generalization performance of models by finding flat minima in the loss landscape, which are believed to lead to better generalization. The authors provide theoretical analysis of SABO's convergence rate and generalization bound, and they demonstrate its effectiveness through experiments on synthetic problems and black-box prompt fine-tuning tasks.
weaknesses: Overall, I see no major weaknesses for this paper and some of these points raised here are covered in the questions below which might clear up these apparent weaknesses. For example, in eqn (8) there is a dependence on the parameter $\rho$. The performance of SABO would be sensitive to the neighborhood size $\rho$ in (8) and tuning this parameter would be critical for successful application of this technique. It would be nice to see methods which can automatically adapt this parameter in a useful way to reduce dependence on manual tuning. 

The authors do not mention in detail the increased computational cost of the sharpness-aware updates in SABO. One way to decrease computational cost, the authors assume the covariance matrix $\Sigma$ is diagonal. This seems like a fairly naive assumption and one I would assume is too simplistic in general. To address this problem in SAM for traditional ML settings, variants of SAM have been introduced. Could similar techniques be applied here?

The authors focus on Gaussian distributions for the search distribution. Is that a fair simplification? How effective or sensitive is SABO with other types of distributions that are likely to occur for parameters; e.g. heavy-tail distributions, or other posteriors that are more informative.

---

confidence: 4
contribution: 3
initial_rating: 6
presentation: 3
questions: Please see weakness
soundness: 3
strengths: - This paper is novel in that it extends sharpness-aware minimization (SAM) principles, typically used in gradient-based settings, to black-box optimization where gradients are unavailable.
- This paper is supported by rigorous theoretical analysis, including convergence rates and generalization bounds in both full-batch and mini-batch settings. 
- The empirical results highlight SABO’s practical value and potential to improve generalization in real-world applications where direct gradient access is unavailable.
summary: The paper introduces the Sharpness-Aware Black-Box Optimization (SABO) algorithm, enhancing model generalization in black-box optimization tasks by incorporating sharpness-aware minimization principles. Unlike traditional methods that may lead to poor generalization by focusing solely on training loss, SABO optimizes a Gaussian parameter distribution, finding "flat" minima that improve robustness. Through reparameterization and stochastic gradient approximation, SABO adapts Sharpness-Aware Minimization (SAM) for use without gradient information, making it suitable for black-box settings. Theoretical results demonstrate SABO's convergence and generalization bounds, and empirical evaluations on synthetic problems and prompt fine-tuning tasks show that SABO outperforms existing methods in high-dimensional settings, positioning it as a powerful tool for black-box optimization in machine learning.
weaknesses: - Due to its reliance on iterative updates to a Gaussian distribution, along with stochastic gradient approximations and KL divergence constraints, this method may incur higher computational costs per iteration compared to some existing black-box optimization methods, especially in high-dimensional settings. A discussion of the time complexity for various black-box optimization algorithms would strengthen the paper.

- This paper introduces additional hyperparameters (e.g., neighborhood size for KL divergence) that require careful tuning to achieve optimal performance. Discussing the sensitivity of these hyperparameters would improve the analysis.

- Although the proposed method is evaluated on prompt fine-tuning tasks, the backbone model used in the experiments may not be sufficiently representative. Evaluating it on a larger model, such as LLaMA-7B, would provide a stronger assessment of its effectiveness.

---

confidence: 4
contribution: 2
initial_rating: 3
presentation: 3
questions: - Could you provide a table showing how selecting the neighborhood size $\rho$ affects both convergence and performance?
- What specific $\rho$ would you recommend to ensure stability?
- Could the authors include a comparison of computational usage against other methods in Table 1?
- Could the authors provide a comparison of your theoretical and empirical results with the paper Lin et al. [2022]?
- See the Weakness part.

Tianyi Lin, Zeyu Zheng, and Michael Jordan. Gradient-free methods for deterministic and stochastic nonsmooth nonconvex optimization. Advances in Neural Information Processing Systems, 35:26160–26175, 2022.
soundness: 2
strengths: - The paper is well-written and easy to follow.
- The proposal of integrating SAM techniques into black-box optimization is interesting. They show clearly how to achieve the update algorithm. The paper provides proofs for the convergence rate and generalization bounds.
- Comprehensive experiments are conducted on synthetic benchmarks and real-world tasks (black-box prompt fine-tuning for large language models). Results show that SABO outperforms other methods.
summary: This paper studies Sharpness-Aware Minimization in the context of black-box optimization. It introduces Sharpness-Aware Black-Box Optimization (SABO), which is the first to design a sharpness-aware minimization strategy specifically for black-box optimization. The key contributions include the formulation of the SABO algorithm, convergence results, generalization bounds, and extensive empirical validation through experiments on both synthetic benchmark functions and black-box prompt fine-tuning tasks.
weaknesses: - The paper lacks evidence that SABO can find flat minima, even though this is their motivation for improving model generalization.
- The paper lacks a strong motivation for why we need Sharpness-Aware Minimization for black-box optimization, aside from the fact that it may help improve generalization. Could the authors identify the problems that existing black-box optimization methods encounter and demonstrate how SABO addresses these issues?
- Assumption 4.1 is narrow because it includes both L-smoothness and H-Lipschitz conditions; typically, we only assume one of them. For example, the simple function $f(x) = x^2$ satisfies L-smoothness but does not satisfy H-Lipschitz.
- In Theorem 4.3, the assumption is too strong for real applications in neural networks. The result in Lin et al. [2022] is more general.
- Theorem 4.6 is not strong because the assumptions are overly restrictive, which limits the practical applicability of this theorem.
- The empirical convergence results of SABO in synthetic problems do not show much improvement over INGO.
- The empirical results in Table 1 are not very significant, showing only around a 0-1\% improvement, even for datasets with low accuracy such as RTE and SNLI.

Tianyi Lin, Zeyu Zheng, and Michael Jordan. Gradient-free methods for deterministic and stochastic nonsmooth nonconvex optimization. Advances in Neural Information Processing Systems, 35:26160–26175, 2022.

Meta-review: This paper extends the idea of Sharpness-Aware Minimization (SAM) to black-box optimization settings, where an objective function is optimized using only function queries. Extending SAM to this setting is non-trivial because the core idea of SAM involves computing the worst-case loss perturbation via a gradient step, but in black-box optimization, only function values can be queried, not gradients. The key idea behind the proposed approach is to consider a probability distribution over the queries, such as a normal distribution centered at each query point, and then work with the expected function value instead of the actual queried value. By differentiating the expected function, one can obtain an approximation of the objective function's gradient using only the available query points. The paper also presents a convergence analysis for the convex case, as well as empirical results on black-box prompt fine-tuning tasks. Given the prevalence of black-box optimization in practical settings like reinforcement learning and prompt fine-tuning, the proposed idea is practically relevant and potentially impactful.

Out of four reviews, three reviewers rated the paper on the accept side (including one strong accept), and one reviewer rated the paper as marginally below the acceptance threshold. The reviewers raised some concerns about the limitations of the theoretical results (specially the convergence proof relying on convexity assumption), as well as some clarifying questions about the connection between the proposed method and SAM and the computational cost of the proposed method. The authors responded to these issues, and as a result, some reviewers increased their initial scores.

In concordance with the reviewers, I believe the paper provides solid and potentially impactful contributions to the field, and thus I recommend acceptance.

Decision: Accept (Poster)