OptEx: Expediting First-Order Optimization with Approximately Parallelized Iterations

Authors: yao  shu # †
jiongfeng  fang
ying tiffany he
fei  richard
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 05:06:57

Abstract: First-order optimization (FOO) algorithms are pivotal in numerous computational domains, such as reinforcement learning and deep learning. However, their application to complex tasks often entails significant optimization inefficiency due to their need of many sequential iterations for convergence. In response, we introduce first-order optimization expedited with approximately parallelized iterations (OptEx), the first general framework that enhances the time efficiency of FOO by leveraging parallel computing to directly mitigate its requirement of many sequential iterations for convergence. To achieve this, OptEx utilizes a kernelized gradient estimation that is based on the history of evaluated gradients to predict the gradients required by the next few sequential iterations in FOO, which helps to break the inherent iterative dependency and hence enables the approximate parallelization of iterations in FOO. We further establish theoretical guarantees for the estimation error of our kernelized gradient estimation and the iteration complexity of SGD-based OptEx, confirming that the estimation error diminishes to zero as the history of gradients accumulates and that our SGD-based OptEx enjoys an effective acceleration rate of Θ(√N ) over standard SGD given parallelism of N, in terms of the sequential iterations required for convergence. Finally, we provide extensive empirical studies, including synthetic functions, reinforcement learning tasks, and neural network training on various datasets, to underscore the substantial efficiency improvements achieved by our OptEx in practice.

Main paper: 
 introduction: 
 First-order optimization (FOO) algorithms, such as stochastic gradient descent (SGD) [1] , Nesterov Accelerated Gradient (NGA) [2] , AdaGrad [3] , Adam [4] etc., have already been the cornerstone of many computational disciplines, driving advancements in areas ranging from reinforcement learning [5] to machine learning [6] . These algorithms, which are widely known for their straightforward form of iterative gradient-based updates, are fundamental in solving both simple and intricate optimization problems. However, their applications usually encounter substantial optimization inefficiency, especially when addressing complex functions that not only are expensive in evaluating their function values and gradients but also necessitate a large number of sequential iterations to converge in practice, e.g., deep reinforcement learning [7] and neural network training [8] .To this end, parallel computing has been widely used in the literature to considerably enhance the optimization (e.g., time) efficiency of FOO by reducing the evaluation cost of function and gradient per iteration in FOO [9] . For instance, in the field of neural network training, techniques that are based # Correspondence to: Yao Shu <shuyao@gml.ac.cn> on parallel computing, e.g., data parallelism [8, 10, 11, 12] , model parallelism [13] , and pipeline parallelism [14, 15] , have been employed to reduce the evaluation time of loss function and parameter gradient by processing multiple input samples and network components concurrently. However, to the best of our knowledge, few efforts have been devoted to leveraging parallel computing to reduce the number of sequential iterations required for convergence to mitigate the optimization inefficiency in FOO. Different from the methods of reducing the evaluation time per iteration during optimization, which requires specialized human efforts in a specific domain (e.g., neural network training), the reduction of sequential iterations is likely to be more general since no such specialized domain efforts are required and thus shall enjoy a wider application in practice. This underscores the need to explore the potential of parallelizing sequential iterations in standard FOO.However, the inherent iterative dependency in FOO where the output of each iteration servers as the input of the next iteration, poses a significant barrier to independent and concurrent iteration execution, thereby making it nearly impossible to realize iteration parallelism within standard FOO. To this end, we develop a novel framework called first-order optimization expedited with approximately parallelized iterations (OptEx) that is capable of bypassing the challenge of inherent iterative dependency in standard FOO and therefore make parallelized iterations in FOO possible. Specifically, our framework begins with a novel kernelized gradient estimation strategy, which uses the history of gradients during optimization to predict the gradients for any input within the domain such that these estimated gradients can be used in standard FOO algorithms to determine the inputs for the next few iterations. We further introduce the techniques of separable kernel function and local history of gradients to enhance the computational efficiency of this gradient estimation (Sec. 4.1). We then apply standard FOO algorithms with this kernelized gradient estimation to determine the inputs for the next N sequential iterations efficiently (namely proxy updates), aiming to approximate the ground-truth sequential updates and bypass the iteration dependency in standard FOO (Sec. 4.2). Lastly, we complete our approximately parallelized iterations for standard FOO by leveraging parallel computing with parallelism of N to concurrently execute standard FOO algorithms over these N inputs obtained from our proxy updates using the ground-truth gradients (Sec. 4.3) .Apart from proposing our innovative OptEx framework, we further establish rigorous theoretical guarantees and extensive empirical studies underpinning its efficacy. Specifically, we give a theoretical bound for the estimation error of our kernelized gradient estimation. Remarkably, this error approaches zero asymptotically as the number of historical gradients increases, ranging across a broad spectrum of kernel functions. This suggests that our kernelized gradient estimation can facilitate effective proxy updates to help parallelize sequential iterations in FOO (Sec.5.1). Building on this, we delineate both upper and lower bounds for the sequential iteration complexity of our SGD-based OptEx, showing that our SGD-based OptEx is able to reduce the sequential iteration complexity of standard FOO algorithms at a rate of Θ( √ N ) with parallelism of N (Sec.5.2). Finally, through extensive empirical studies, including the optimization of synthetic functions, reinforcement learning tasks, and neural network training on both image and text datasets, we demonstrate the consistent advantages of our OptEx in expediting existing FOO algorithms (Sec. 6).To summarize, our contribution to this work includes:• To the best of our knowledge, we are the first to develop a general framework (i.e., OptEx) that can leverage parallel computing to approximately parallelize the sequential iterations in FOO, thereby considerably reducing the sequential iteration complexity of FOO algorithms. • We provide the first upper and lower iteration complexity bound for SGD-based OptEx, which gives an effective acceleration rate of Θ( √ N ) with parallelism of N . • We conduct extensive empirical studies, including the optimization of synthetic function, reinforcement learning tasks, and neural network training on both image and text datasets, to support the efficacy of our OptEx framework. 
related work: 
 Reduction of Iteration Complexity. In the literature, various techniques have been developed to enhance the optimization efficiency of FOO by improving their sequential iteration complexity. For example, variance reduction strategies [16, 17, 18] have been proposed to accelerate stochastic optimization by effectively reducing the gradient variance and therefore aligning the iteration complexity of SGD with that of gradient descent (GD) in expectation. These strategies usually yield significant improvements in high-variance problems whereas their compelling performance is hard to extend to low-variance scenarios and deterministic contexts. Meanwhile, adaptive gradient methods, e.g., AdaGrad [3] , Adam [4] , and AdaBelief [19] , have been introduced to employ an adaptive learning rate for a better-performing optimization where fewer iterations are required for convergence. Furthermore, acceleration techniques like the Nesterov method [2] and momentum-based updates [20] have also been proven to be capable of reducing the sequential iteration complexity for GD and SGD efficiently. Orthogonal to these established methodologies, our paper introduces parallel computing as a distinct and innovative strategy to further decrease the sequential iteration complexity of FOO.Of note, such an approach not only stands independently but also offers potential for synergistic integration with existing methods, promising enhanced optimization outcomes.Reduction of Time Complexity Per Iteration using Parallel Computing. In the realm of enhancing the computational efficiency of FOO, parallel computing has emerged as a rescue by reducing the time complexity per iteration in FOO. Particularly in the field of neural network training, data parallelism [8, 10, 11, 12] has been introduced to evaluate the gradients of model parameters w.r.t mini-batch input samples simultaneously. In addition to data parallelism, model parallelism [13] has been developed to process various neural network components concurrently. Furthermore, pipeline parallelism [14, 15] divides the neural network into stages and assigns each stage to a different device, allowing different stages of the computation to be executed in parallel across the pipeline.However, the tailored nature of these methods constrains their application to wider contexts. Contradictory to these case-specified solutions, this paper proposes a general framework that can leverage parallel computing to enhance the optimization efficiency of FOO in wide practical applications. 
problem setup: 
 In this paper, we aim to enhance the optimization efficiency of the following stochastic minimization problem by leveraging parallel computing with parallelism of N : minEQUATIONHere, ∇f (θ) is assumed to follow a specific Gaussian distribution, i.e., ∇f (θ) ∼ N (∇F (θ), σ 2 I) for any θ ∈ R, which has already been widely used in the literature [21, 22, 23] . Besides, we adopt a common assumption that ∇F is sampled from a Gaussian process, i.e., ∇F ∼ GP(0, K(•, •)) [24, 25, 26] . Of note, (1) has found extensive applications in practice, e.g., neural network training [27] and reinforcement learning [28] . Importantly, although our primary focus is on this stochastic optimization, our method can also be applied to deterministic optimization (evidenced in Sec. 6.1).Standard FOO algorithms commonly optimize (1) in an iterative and sequential manner:EQUATION) where t is the iteration number. Ideally, if parallel computing can be used to parallelize the sequential iterations in FOO (i.e., to execute several sequential iterations simultaneously), it will be able to lead to a noticeable improvement in its optimization efficiency since fewer sequential iterations will be required for convergence. Unfortunately, there is an inherent iterative dependency in standard FOO, that is, the output of each iteration t (e.g., θ t ) is the input of the next iteration t + 1. Such an iterative and sequential process makes it nearly impossible to attain θ t and θ t+1 concurrently, and therefore parallelize the iterations for established FOO algorithms. 
the optex framework: 
 To this end, we introduce the first general framework in Algo. 1 with a detailed illustration in Fig. 1 , namely first-order optimization expedited with approximately parallelized iterations (OptEx), to overcome the aforementioned inherent iterative dependency in FOO and facilitate the realization of parallelized iterations therein. To achieve this, we first propose a kernelized gradient estimation with the technique of separable kernel function and local history of gradient to efficiently and effectively estimate the gradient at any input in the domain (Sec. 4.1). We then follow standard FOO algorithms with this kernelized gradient estimation to approximate the inputs for the next N sequential iterations to be parallelized (Sec.4.2), aiming to overcome the inherent iterative dependency in FOO. Lastly, we finish our approximately parallelized iterations by leveraging parallel computing to run standard FOO algorithms on these N inputs concurrently using the ground-truth gradients (Sec. 4.3). Input: FO-OPT, k(•, •), θ 0 , T , N , G = ∅ for sequential iteration t ∈ [T ] do Initialization: θ t,0 ← θ t-1 Update µ t (θ) using (4.1) with G for proxy step s ∈ [N -1] do θ t,s ← FO-OPT(θ t,s-1 , µ t (θ t,s-1 )) for process i ∈ [N ] in parallel do Sample f to evaluate ∇f (θ t,i-1 ) θ (i) t ← FO-OPT(θ t,i-1 , ∇f (θ t,i-1 )) G ← G ∪ {(θ t,i-1 , ∇f (θ t,i-1 ))} θ t ← θ (N ) tAs mentioned in our Sec. 3, ∇F is assumed to be sampled from a Gaussian process, i.e., ∇F ∼ GP(0, K(•, •)) with kernel function K. Then, for every sequential iteration t of Algo. 1, conditioned on the history of gradients during optimizationG ≜ {(θ τ , ∇f (θ τ )} N (t-1) τ =11 , ∇F then follows the posterior Gaussian process: ∇F ∼ GP µ t (•), Σ 2 t (•, •) with the mean function µ t (•) and the covariance function Σ 2 t (•, •) defined as below [24] :EQUATIONwhere vec(•) vectorizes a matrix into a column vector,G t ≜ [∇f (θ τ )] N (t-1) τ =1 is a d × N (t -1)- dimensional matrix, V ⊤ t (θ) ≜ [K(θ, θ τ )] N (t-1) τ =1is a d × N (t -1)d-dimensional matrices, andU t ≜ K(θ τ , θ τ ′ ) N (t-1) τ,τ ′ =1 is a N (t -1)d × N (t -1)d-dimensional matrices.We therefore propose to use µ t (•) to estimate the gradient at any input θ ∈ R d , that is,EQUATION)and covariance Σ 2 (θ) ≜ Σ 2 (θ, θ) to measure the quality of this gradient estimation in a principled way, which will be further theoretically supported in our Sec. 5.1.However, for every sequential iteration t of Algo. 1 with (3), it will incur a computational complexity of O(N 3 (t -1) 3 d 3 ), along with a space complexity of O(N (t -1)d). Practically, this presents a significant challenge in scenarios with a large input dimension d or requiring a substantial number T of sequential iterations for convergence, such as in neural network training [8] . To mitigate these complexity issues, we introduce two techniques: the separable kernel function and the local history of gradients, to reduce both the computational and space complexities associated with our kernelized gradient estimation, thereby enhancing its efficiency and practical applicability.Separable Kernel Function. Let K(•, •) = k(•, •) I where k(•,•) produces a scalar value and I is a d × d identity matrix, and define theN (t -1)-dimensional vector k ⊤ t (θ) ≜ [k(θ, θ τ )] N (t-1) τ =1, andN (t -1) × N (t -1)-dimensional matrix K t ≜ [k(θ τ , θ τ ′ )] N (t-1) τ =τ ′ =1, we can prove that the Gaussian process in (3) can be simplified as the Gaussian process in Prop. µ t (θ) = k ⊤ t (θ) K t + σ 2 I -1 G t ⊤ , Σ 2 t (θ, θ ′ ) = k(θ, θ ′ ) -k ⊤ t (θ) K t + σ 2 I -1k t (θ ′ ) I . 1 We slightly abuse the notation f to denote the different functions that are randomly sampled per iteration and (θ τ , ∇f (θ τ ) to denote a historical evaluation till sequential iteration t -1 with parallelism of N .Its proof is in Appx. A.1. Prop. 4.1 shows that with a separable kernel function K(•, •) = k(•, •) I, the multi-output Gaussian process in a d-dimensional space can be effectively decoupled into d independent single-output Gaussian processes. Each of these processes results from the same scalar kernel function k, leading to a uniform posterior form shared by all these processes, i.e., the expressionk ⊤ t (θ) K t + σ 2 I -1 in µ t (θ) and k(θ, θ ′ )-k ⊤ t (θ) K t + σ 2 I -1 k t (θ ′ ) in Σ 2 t (θ, θ ′ ).This thus considerably diminishes the computational complexity, now quantified as O(N 3 (t-1) 3 +N (t-1)d), resulting in a more computationally efficient gradient estimation in practice.Local History of Gradients. Conventional FOO algorithms predominantly operate by optimizing within a localized region neighboring the initial input θ 0 [29] . This therefore indicates that our Algo. 1 only requires precise gradient estimation within a local region. In this context, the use of a local gradient history is posited as sufficiently informative for effective kernelized gradient estimation, which can be supported by the theoretical results in [30] and the empirical evidence in our Sec. 6. As a result, rather than relying on a complete gradient history, we propose to use a localized gradient history of size T 0 that neighbors θ to estimate the gradient at θ. This strategic modification results in a substantial reduction of computational complexity to O(T 3 0 + T 0 d) as well as a corresponding decrease in space complexity to O(T 0 d), which is especially beneficial in the situations where T 0 is considerably smaller than N (t -1) for t ∈ [T ].The ability of our kernelized gradient estimation to provide gradient estimation at any input θ then enables the application of a multi-step gradient estimation. This helps to approximate the inputs for the next N sequential iterations {θ τ +i } N -1 i=0 to be parallelized in standard FOO, given θτ . Specifically, in the context of our Algo. 1, for every sequential iteration t ∈ [T ], by employing a first-order optimizer (FO-OPT), we can approximate the inputs required by our parallelized iteration in Sec. 4.3 sequentially as below through our multi-step proxy updates (line 4-5 of Algo. 1).θ t,s = FO-OPT(θ t,s-1 , µ t (θ t,s-1 )), ∀s ∈ [N -1] .(5) Intuitively, these proxy updates imitate the sequential iterations in standard FOO by using only the estimated gradients in our Sec. 4.1. We will show that these proxy updates can provide a reasonably good approximation of the ground-truth updates in Sec. 5.1. Meanwhile, despite the iterative and sequential nature of (5), our proxy updates based on operations on relatively small-sized matrices (refer to the Prop. 4.1) will still be able to provide significantly enhanced computational efficiency compared to the ground-truth updates based on expensive evaluation of function values and gradients in complex tasks, like neural network training. This effectiveness and efficiency of (5) thus render it an essential foundation for achieving parallelized iterations and improved the optimization efficiency in FOO.Upon obtaining the inputs {θ t,s-1 } N s=1 for the next N sequential iterations to be parallelized, we then finish our approximately parallelized iteration by executing standard FOO algorithms over each of {θ t,s-1 } N s=1 based on the ground-truth gradients {∇f (θ t,s-1 )} N s=1 in parallel (line 6-9 of Algo. 1, see also the processes in Fig. 1 ). That is,EQUATION)After that, the final input θ t = θ (N ) t will be used in the next sequential iteration (line 10 of Algo. 1). Of note, central to the approximately parallelized iterations in our OptEx framework is the necessity of evaluating the gradients ∇f (θ t,s-1 ) N s=1 in our Algo. 1. These evaluations in fact play pivotal roles in reducing the estimation error of our kernelized gradient estimation and consequently improving the performance of our OptEx by augmenting the gradient history near the input θ t with N more evaluations, which will be supported by the theoretical results in our Sec. 5 and the empirical evidence in our Appx. B.3. 
theoretical results: 
 To begin with, we formally present the assumptions mentioned in our Sec. 3 as below.Assumption 1. ∇f (θ) -∇F (θ) follows N 0, σ 2 I for any θ ∈ R d . Assumption 2. ∇F is sampled from a Gaussian process GP (0, K(•, •)) where K(•, •) = k(•, •) I and |k(θ, θ)| ≤ κ for any θ ∈ R d .Note that the Assump. 1 has already been widely employed in the literature [21, 22, 23] . Meanwhile, it is also common to assume that F is sampled from a Gaussian process [24, 31] , implying that ∇F follows a Gaussian process as well [24, 25, 26] (Assump. 2), i.e., ∇F can be any function in this prior. The inclusion of a separable kernel function in Assump. 2 aims to enhance the efficiency of our kernelized gradient estimation in Sec. 4.1 and simplify our theoretical analyses below, whereas our conclusions apply to non-separable kernel functions as well by following our proof techniques.Following the principled idea in kernelized bandit [32, 33] and Bayesian Optimization [34, 31] , we define the maximal information gain as belowEQUATIONwhereI(vec(G n ); vec(∇ n )) is the mutual information between G n ≜ [∇f (θ i )] n i=1 and ∇ n ≜ [∇F (θ i )] n i=1 .In essence, γ n encapsulates the maximum amount of information about ∇F that can be gleaned from observing any set of n evaluated gradients, represented as G n , which is known to be problem dependent measure that is highly related to the kernel function k(•, •) [32] . Built on this notation, we then provide the following theoretical result for our gradient estimation.Theorem 1 (Gradient Estimation Error). Let δ ∈ (0, 1) and α ≜ d + ( √ d + 1) ln(1/δ). Given Assump. 1 and 2, let |G| = T 0 -1 for any sequential iteration t in Algo. 1, then for any θ ∈ R d , t > 0, with a probability of at least 1 -δ, ∥∇F (θ) -µ t (θ)∥ ≤ α Σ 2 (θ) where κ (κ + 1/σ 2 ) T 0 -1 ≤ Σ 2 (θ) ≤ 4 max{κ, σ 2 }γ T 0 T 0 d .The proof is in Appx. A.2. It is important to note that since FOO pertains to local optimization, the global fulfillment of Assump. 2 is not a prerequisite. That is, the assumption that ∇F is sampled from GP(0, K) within a local region will already be sufficient for our kernelized gradient estimation in Sec. 4.1 to achieve accurate gradient estimation in practice. Our Sec. 6 will later evidence this empirically. Thm. 1 with upper bound on Σ 2 (θ) illustrates that the efficacy of our kernelized gradient estimation in the worst case will enjoy a polynomial error rate of O γ T 0 /T 0 . This means that if γ T 0 /T 0 will asymptotically approach zero w.r.t. T 0 , the error of our kernelized gradient estimation method will become significantly small given a large number of evaluated gradients T 0 . This consequently facilitates the effectiveness of our proxy updates in (5) built on our kernelized gradient estimation to approximate the ground-truth updates when |G| is sufficiently large. Meanwhile, Thm.1 with lower bound on Σ 2 (θ) illustrates that our kernelized gradient estimation in the best case may achieve an exponential error rate of O κ/(κ + 1/σ 2 ) T 0 -1 , which thus further elaborates the efficacy of kernelized gradient estimation in Sec. 4.1 and proxy updates in Sec. 4.2.It is important to note that the ratio γ T 0 /T 0 has been demonstrated to asymptotically approach zero for a range of kernel functions, as evidenced in existing literature [35] . This therefore underpins the establishment of concrete error bounds for our kernelized gradient estimation where notation O is applied to hide the logarithmic factors, delineated as follows:Corollary 1 (Concrete Error Bounds). Let k(•, •) be the radial basis function (RBF) kernel, then∥∇F (θ) -µ t (θ)∥ = O T -1/2 0 . Let k(•,•) be the Matérn kernel where ν is the smoothness parameter, then∥∇F (θ) -µ t (θ)∥ = O T -ν/(2ν+d(d+1)) 0 .Cor. 1 elucidates that with kernel functions such as RBF and Matérn kernel, the error in our kernelized gradient estimation indeed will diminish asymptotically w.r.t. T 0 . That is, as T 0 increases, the estimation error ∥∇F (θ) -µ t (θ)∥ decreases and consequently the proxy updates in (5) become closer to the ground-truth updates. It is important to note that this reduction typically follows a nonlinear trajectory, suggesting that the effect of an increasing T 0 on our kernelized gradient estimation diminishes when T 0 is reasonably large. This consequently affirms the reasonability of our utility of local history for gradient estimation in Sec. 4.1, which leads to not only accurate but also efficient gradient estimations.We first introduce Assump. 3, which has been widely applied in stochastic optimization [16, 36] , to underpin the analysis of sequential iteration complexity of our OptEx framework.Assumption 3. F is L-Lipschitz smooth: ∇F (θ) -∇F (θ ′ ) ≤ L θ -θ ′ for any θ, θ ′ ∈ R d .To simplify the analysis, we primarily prove the sequential iteration complexity of our SGD-based OptEx where we use min τ ∈[N T ) ∥∇F (θ τ )∥ 2 to denote the minimal gradient norm we can achieve within the whole optimization process when applying our OptEx with T sequential iterations and parallelism of N for a clear and fair comparison with standard SGD. Notably, our analysis can also be extended to other FOO-based OptEx by following similar proof idea.Theorem 2 (Upper Bound). Let δ ∈ (0, 1), ∆ ≜ F (θ) -inf θ F (θ), β ≜ max{κ, σ 2 } and ρ ≜ (1 -1 N ) 4βγ T 0 σ 2 T 0 + 1 N . Under Assump. 1-3, by choosing T ≥ 2∆L N σ 2 ρ , η = 2∆ N LT σ 2 ρand |G| = T 0 -1 for our SGD-based Algo. 1, with a probability of at least 1 -δ,EQUATION)The proof of Thm. 2 is in Appx. A.3. Of note, our Thm. 2 with N = 1 aligns with the established upper bound for standard SGD, as discussed in [36] . Importantly, our Thm. 2 elucidates that with parallelism N > 1, our SGD-based OptEx algorithm can expedite the standard SGD by a factor of at least N/ρ, where 1/ρ quantifies the impact of the error introduced by our kernelized gradient estimation. This efficiency gain can be further amplified as the accuracy of our kernelized gradient estimation increases (i.e., a decrease in ρ), which can be achieved by augmenting the number T 0 as discussed in our Sec. 5.1. In addition, Thm. 2 also demonstrates that for a fixed learning rate η, there exists a constant N max , e.g., N max = 2∆/(η 2 LT σ 2 ρ) in Thm. 2, the parallelism N should roughly remain below to ensure the fastest convergence of function F to a stationary point. In contrast, if N exceeds N max , our SGD-based OptEx will underperform due to the increased gradient estimation error. This observation is further supported by the results presented in Appx. B.3. However, when the learning rate η is relatively small (e.g., during fine-tuning in practice), the parallelism N can be significantly larger to achieve a further improved speedup.Theorem 3 (Lower Bound). Let δ ∈ (0, 1), ∆ ≜ F (θ) -inf θ F (θ), β ≜ max{κ, σ 2 }, and β ≜ min{κ/(κ + 1/σ 2 ) T 0 -1 , σ 2 }. Then, for any L > 0, ∆ > 0, N ≥ 1, T ≥ 1 and η ∈ [0, 1/L), there exists a F on R d (∀d > d 0 = O β/(∆L 2 ) ln N T /δ ) satisfying Assump.EQUATION)The proof of Thm. 3 is in Appx. A.4. Note that when N = 1, Thm. 3 aligns with the recognized lower bound for SGD, as elucidated in [37] . Thm. 3 illustrates that, with parallelism of N , our SGD-based OptEx can potentially accelerate standard SGD by up to √ N /(κ/(σ 2 (1 + 1/σ 2 ) T 0 -1 )), under the condition that κ/(1 + 1/σ 2 ) T 0 -1 ≤ min{∆L, 1, σ 2 }. This upper limit in fact corresponds with the lower bound of the variance in our kernelized gradient estimation, as established in Thm. 1. Essentially, the agreement between Thm. 2 and Thm. 3, in the aspect of parallelism N , demonstrates the tightness of our sequential complexity analysis for SGD-based Algo. 1. Finally, the combination of Thm. 2 and Thm. 3 enables us to specify the effective acceleration that can be achieved by our SGD-based OptEx tightly, as shown in our Cor. 2 below.  
experiments: 
 In this section, we use extensive experiments to show that our OptEx framework can considerably enhance the efficiency of FOO with parallel computing, including synthetic experiments (Sec. 6.1), reinforcement learning (Sec. 6.2) and neural network training on various datasets (Sec. 6.3).Here, we utilize synthetic functions to demonstrate the enhanced performance of our OptEx framework compared to existing baselines, including the standard FOO algorithm, namely Vanilla, and FOO with ideally parallelized iterations, namely Target, which ideally but impractically utilizes the ground-truth gradient to obtain the inputs for the iterations to be parallelized. More specifically, the Vanilla baseline is equivalent to Algo. 1 with parallelism of N = 1, and the Target baseline is equivalent to Algo. 1 with µ t (θ t,s-1 ) being replaced with ∇f (θ t,s-1 ), indicating the desired parallelized iteration we aim to approximate. We have also provided a comprehensive illustration of these baselines in Appx. B.1 and detailed experimental setup applied here in Appx. B.2.1.The results in Fig. 2 with σ 2 = 0 and N = 5 have demonstrated the efficacy of our OptEx framework for deterministic optimization (i.e., σ 2 = 0). Specifically, Fig. 2 shows that OptEx consistently achieves a notable speedup in optimization efficiency measured by the number of sequential iterations, which is at least 2× more efficient than the Vanilla baseline, when optimizing with parallelism of N = 5 to reach an equivalent level of optimality gap. This is roughly in line with the result of our Cor. 2, implying the validity of our Cor. 2. Meanwhile, although our OptEx framework slightly underperforms the Target baseline, such a phenomenon is in fact quite reasonable since the Target baseline can leverage the ground-truth gradient whereas OptEx relies on the kernelized gradient estimation with estimation error bounded in Thm. 1 to parallelize sequential iterations. This also aligns with the insight from our iteration complexity analysis in Thm. 2. Overall, the results in Fig. 2 have provided strong empirical support for the efficacy of our OptEx in expediting FOO, as theoretically justified in our Sec. 5.2. We also present a number of ablation studies as well as analyses in Appx. B.3 to examine the effects of different components in our proposed OptEx framework on its effectiveness.We proceed to compare our OptEx framework with previously established baselines under various reinforcement learning tasks with different parameter dimension d from the OpenAI Gym suite [38] , with the deployment of DQN agents [39] . Here, the parallelism parameter is set to be N = 4 and a detailed experimental setup is provided in Appx. B.2.2. The results are presented in Fig. 3 . As illustrated in Fig. 3 , the integration of parallel computing techniques, including Target and OptEx, considerably outperforms the traditional Vanilla baseline in terms of the optimization efficiency quantified by the number of sequential iterations. More importantly, amongst these methodologies, OptEx consistently demonstrates a more stable and superior improvement on the optimization effi- is applied. This is likely because the gradient variance (i.e., ∥Σ 2 (θ)∥) in our OptEx framework can asymptotically approach zero by using a large number of history of gradient (refer to our Sec. 4.1), whereas the gradient variance in the Target baseline remains the same.At last, we examine the efficacy of our OptEx in expediting the optimization (i.  
conclusion: 
 In conclusion, our OptEx framework represents a significant advancement in FOO. By leveraging kernelized gradient estimation to enable approximately parallelized iterations, OptEx effectively re-duces the number of sequential iterations required for convergence and thus addresses the traditional inefficiencies of FOO. Theoretical analyses and extensive empirical studies validate the reliability and efficacy of OptEx, confirming its potential to expedite optimization processes across various applications. Of note, a limitation of OptEx is the additional storage and computational cost introduced by the kernelized gradient estimation, which we aim to mitigate further in the future work. 
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Answer: [Yes]Justification: See our Sec. B for the details.Guidelines:• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper). 
code of ethics: 
 Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines?Answer: [Yes] Justification: The research conducted in the paper indeed conforms with the NeurIPS Code of Ethics.Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction). 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?Answer: [NA] Justification: We do not see any societal impact of the work performed.Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.• Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)? Answer: [NA] Justification: This paper poses no such risks. Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort. 12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected? Answer: [Yes] Justification: All the data and codes used in the paper are open-sourced. Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided. • If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. 
 
 Appendix: A.1 Proof of Proposition 4.1Recall that we have definedk ⊤ t (θ) ≜ [k(θ, θ τ )] N (t-1) τ =1, andK t ≜ [k(θ τ , θ τ ′ )] N (t-1) τ =τ ′ =1. Let ⊗ denote the Kronecker product, by introducing the fact that K(•, •) = k(•, •) I into V ⊤ t (θ) and U t from the Gaussian process posterior (3), we have thatEQUATION)Similarly,EQUATION)By introducing the results above into the posterior mean and variance in (3), we haveEQUATION)where (c) come from the bi-linearity of the Kronecker product, i.e., (A + B)⊗ C = A ⊗ C + B ⊗ C while (d) is from the inverse of the Kronecker product, i.e., (A ⊗ B) -1 = A -1 ⊗ B -1 .In addition, (e) is due to the mixed-product property of the Kronecker product, i.e., (A ⊗ B)(C ⊗ D) = (AC) ⊗ (BD), and (f ) results from the mixed Kronecker matrix-vector product of the Kronecker product, i.e., (A ⊗ B)vec(C) = vec(BCA ⊤ ).Similarly,EQUATIONwhere (b) comes from the result in ( 12), (c) results from the mixed-product property of the Kronecker product and the fact thatk ⊤ t (θ) K t + σ 2 I -1 k t (θ ′) is a scalar. This finally concludes our proof.To begin with, we introduce the following lemmas: 43] ). Let ζ ∼ N (0, I d ) and δ ∈ (0, 1) thenLemma A.1 ([EQUATION)Lemma A.2 (Lemma 2 in Appx. B of [34] ). For any σ ∈ R and any matrix A, the following holdEQUATION)Lemma A.3 (Sherman-Morrison formula). For any invertible square matrix A and column vectors u, v, suppose A + uv ⊤ is invertible, then the following holdsEQUATION)Lemma A.4 (Non-Increasing Variance Norm). Define variance Σ 2 n (θ) ≜ Σ 2 n (θ, θ) with n being the number of gradients employed to evaluate the mean and covariance in Prop. 4.1. Then for any θ ∈ R d and n ≥ 1,EQUATION)Proof. We follow the idea in [34] and [25] to prove it. Specifically, we firstly define k(θ, θ ′ ) = ϕ(θ) ⊤ ϕ(θ ′ ) and ϕ n ≜ [ϕ(θ i )] n i=1 . Then the matrix K n in Prop. 4.1 can be reformulated asEQUATION)and based on the definition ofEQUATIONwhere (c) comes from Lemma A.2 by replacing the matrix A in Lemma A.2 with the matrix ϕ ⊤ n .As a result,Σ 2 n (θ) (a) = σ 2 ϕ(θ) ⊤ Φ -1 n ϕ(θ) I (b) = σ 2 ϕ(θ) ⊤ ϕ n-1 ϕ ⊤ n-1 + σ 2 I + ϕ(θ n )ϕ(θ n ) ⊤ -1 ϕ(θ) I (c) = σ 2 ϕ(θ) ⊤ Φ n-1 + ϕ(θ n )ϕ(θ n ) ⊤ -1 ϕ(θ) I (d) = σ 2 ϕ(θ) ⊤ Φ -1 n-1 ϕ(θ) -σ 2 1 + ϕ(θ n ) ⊤ Φ -1 n-1 ϕ(θ n ) -1 ϕ(θ) ⊤ Φ -1 n-1 ϕ(θ n )ϕ(θ n ) ⊤ Φ -1 n-1 ϕ(θ) I (e) = Σ 2 n-1 (θ) -σ 2 1 + ϕ(θ n ) ⊤ Φ -1 n-1 ϕ(θ n ) -1 ϕ(θ) ⊤ Φ -1 n-1 ϕ(θ n )ϕ(θ n ) ⊤ Φ -1 n-1 ϕ(θ) I (f ) ≼ Σ 2 n-1 (θ) (20) where (b) is due to the fact that Φ n Φ ⊤ n = Φ n-1 Φ ⊤ n-1 + ϕ(θ n )ϕ(θ n ) ⊤ ,and (d) is from Lemma A.3. Finally, (f ) derives from the positive semi-definite property of Φ -1 n-1 ϕ(θ t )ϕ(θ t ) ⊤ Φ -1 n-1 and Φ -1 n-1 , leading to the conclusion of our proof.EQUATION)Proof. Again, we follow the idea in [34] and [25] to prove it. we first prove the following inequalityEQUATIONwhere (c) comes from the fact thatΦ n-1 = Φ -1 n-2 + ϕ(θ n-1 ϕ(θ n-1 ) ⊤ ≽ Φ n-2 ≽ • • • ≽ σ 2 I and (e) is due to the fact that ϕ(θ n ) ⊤ ϕ(θ n ) = k(θ n , θ n ) ≤ κ.Then, based on the reformulation of Σ 2 n (θ) in ( 20), we have thatEQUATIONwhere (c) comes from (22) . This finally concludes our proof.Lemma A.6 (Information Gain). Define G n ≜ [∇f (θ i )] n i=1 , ∇ n ≜ [∇F (θ i )] n i=1 , and K n ≜ k(θ i , θ j ) n i,j=1. The information gain I(vec(G n ); vec(∇ n )) has the following form with Assump. 1, 2:EQUATION)Proof. Based on our Assump. 1, 2, the following holds respectively:EQUATION)Due to the fact that H(z) = 1 2 ln(det(2πeΣ)) if z ∼ N (µ, Σ), the following holds This then concludes our proof.EQUATIONLemma A.7 (Sum of Variance). Define the maximal information gainEQUATION)the following then holdsEQUATION)Proof. To begin with, we show the following inequalities resulting from the matrix determinant lemma:EQUATIONEQUATION)where (a) follows from the reformulation of Σ 2 i (θ) in ( 19), (b) results from the fact that x/2 ≤ ln(1 + x) for any x ∈ (0, 1), (c) derives from ( 29), (d) is from the telescoping sum, (e) is due to the fact that det(Φ 0 ) = det(σ 2 I), (f ) is from the fact that det(AB -1 ) = det(A)/ det(B), (g) comes from the Sylvester's determinant identity, i.e., det(Φ i ) = det(K i +σ 2 I i ) according to the definition of Φ i , and (h) results from the fact that K n = ϕ ⊤ n ϕ n in (18) , the conclusion in Lemma A.6, and the definition of γ n .Following the same idea, given κ > σ 2 , we haveEQUATION)Combining the results in ( 30) and ( 31), we conclude our proof byEQUATION)Proof of our Thm. 1. Since Σ -1 n (θ) [µ n (θ) -∇F (θ)] ∼ N (0, I d ), according to Lemma A.1, for any δ ∈ (0, 1) and α ≜ d + 2( √ d + 1) ln(1/δ), with a probability of at least 1 -δ,EQUATIONwhere (a) is from CauchySchwarz inequality and (b) is from Lemma A.1. By introducing the results in Lemma A.5 and Lemma A.7 into the result above and letting T 0 = n + 1, we conclude our proof.In general, we follow the idea in [36] to give a high probability convergence for our OptEx algorithm.To begin with, we introduce the following definition and lemma.EQUATION)Proof.EQUATIONwhere (c) comes from the definition of σ-sub-Gaussian random variable.EQUATION)andEQUATION)According to our Lemma A.8 and the fact that each dimension of ε t (θ t,s ) follows an independent Gaussian distribution given Assump. 2, the following holdsEQUATION)Based on Markov inequality, we have thatP exp T t=1 N s=1 X t,s > 1 2δ exp T t=1 N s=1 wη 2 L 2 σ 2 (θ t,s-1 ) ≤ E exp T t=1 N s=1 X t,s exp T t=1 N s=1 wη 2 Lσ 2 (θ t,s-1 )/2 /(2δ) ≤ δ . (39)Therefore, with a probability of at least 1 -δ,EQUATION)which leads to the following inequality with w = wEQUATION)Of note, for every proxy step based on SGD:EQUATIONwhere (a) derives from the Lipschitz smoothness of function F (i.e., Assump. 3), (b) comes from the standard SGD update and the definition of ε t (θ t,s ), and (d) is a rearrangement of the results in (c).By introducing the results above into (42) and choosing w -1 = 2βη with β ≜ max{κ, σ 2 }, we haveEQUATIONwhere (a) comes from η ≤ 1/L, (b) is due to the fact that σ 2 (θ t,s-1 ) ≤ max{κ, σ 2 } = β , and (c) is due to the fact thatEQUATION)By rearranging the result in (43) and defining ρ ≜ (1EQUATIONwhere the last inequality comes from the fact that F (θ T ) ≤ inf θ F (θ) and σ 2 (θ t,s-1 ) ≤ 4βγ T 0 /(T 0 d) in (33) .By choosing T ≥ 2∆L N σ 2 ρ and η = 2∆ N T Lσ 2 ρwhere ∆ ≜ F (θ 0 ) -inf θ F (θ), we conclude our proof.Remark 1. The speedup achieved by OptEx matches that of basic sample averaging (i.e., data parallelism) for stochastic optimization with noisy gradients. However, the speedup from OptEx comes from reduced sequential iterations (first term on the RHS in ( 45)), while sample averaging derives from reduced gradient variance (second term on the RHS in ( 45)). When gradient noise is already small or in deterministic optimization (e.g., the experiments in Sec. 6.1), data parallelism may not provide noticeable speedup, but OptEx can still contribute significantly. Overall, OptEx works in a complementary direction to existing parallelization methods, including sample averaging, to speed up first-order optimization, especially when other methods are not applicable or underperforming as discussed in our Sec. 2).We follow the idea in [37] to prove our Thm. 3. We first introduce the following lemma:Lemma A.9 (Lemma B.12 in [44] ). Let X i ∼ N (0, 1) independently, Z ≜ n i=1 X 2 i , and ϵ ∈ (0, 1) thenEQUATION)EQUATIONwhere θ 0 is initialized with N (0, ∆ L I). We abuse ε(θ τ ) to denote the ε(θ t,s-1 ) defined in our (36) . Based on the update rule of stochastic gradient descent, we then have thatEQUATION)Since ε(θ i ) follows N (0, σ 2 (θ i )) independently where we abuse σ 2 (θ i ) to denote the σ 2 (θ t,s-1 ) defined in (36) and θ 0 is initialized with N (0, ∆ L I), we then have thatEQUATION)Let δ ∈ (0, 1) andβ ≜ min{1/(1 + 1/σ 2 ) T 0 κ, σ 2 }, since ∥∇F (θ τ )∥ 2 = L 2 ∥θ τ ∥ 2 ,the results above into Lemma A.9 with ϵ = 1/2 and d ≥ d 0 ≜ 24 ln(N T /δ), with a probability of at least 1 -δ,EQUATIONWhen η ∈ 0, 1/( √ N T L) , we consider the functionEQUATIONwhere θ 0 is initialized withθ ⊤ 0 = d∆ max 1/L, N T τ =1 η , 0, • • • , 0 .Similarly, we haveθ τ =   1 - 1 2 max 1/L, N T τ =1 η   τ θ 0 + τ -1 i=0 η   1 - 1 2 max 1/L, N T τ =1 η   τ -i-1 ε(θ i ) ,(52) andEQUATION)Therefore, let θ (1) τ denote the first element of θ τ , we have E θ (1) τ(a) =   1 - 1 2 max 1/L, N T τ =1 η   τ d∆ max 1/L, N T τ =1 η (b) ≥ d∆ max 1/L, N T τ =1 η exp   ln 1 2 N T τ =1 η max 1/L, N T τ =1 η   (c) ≥ 1 2 d∆ max 1/L, N T τ =1 η (54)where (b) comes from the fact that 1 -z/2 ≥ exp(ln(1/2)z) for all z ∈ [0, 1].In addition, we have var θ (1) EQUATIONEQUATIONwhere β ≜ max{κ, σ 2 }.Let Ψ denote the CDF of standard normal distribution and follow the idea in [37] , by choosingEQUATION   θ (1) τ -E θ (1) τ var θ (1) τ ≥ -1 4 d∆ max 1/L, N T τ =1 η β/L 2     ≥ P    θ (1) τ -E θ (1) τ var θ (1) τ ≥ -EQUATION)That is, with a probability of at least 1 -δ/(N T ), θ (1) τ ≥EQUATION) Since ∥∇F (θ τ )∥ 2 = 1 2 max{1/L, ∑ N T τ =1 η} 2 θ ⊤ e 1 2, we conclude our proof by applying union bound on (58) as belowmin τ ∈[N T ] ∥∇F (θ τ )∥ 2 = 1 N T N T τ =1 ∥∇F (θ τ )∥ 2 ≥ d 0 ∆ 4 max 1/L, N T τ =1 η ≥ d 0 ∆L 4 √ N T (59)where the last inequality comes from the fact that η ∈ [0, 1/( √ N T L)]. This finally concludes our proof.Target In this section, we provide an illustrated comparison between our OptEx and all the baselines at iteration t in Fig. 5 . Notably, the Target baseline represents an ideal parallelization of the Vanilla baseline. However, this is impractical because the ground-truth gradient (i.e., ∇f (•)) required by a process i ∈ [N ] to produce the update can not be obtained before the start of this process. More specifically, this gradient is the outcome at the end of the corresponding process. In contrast, our Op-tEx framework makes use of the kernelized gradient estimation (i.e., µ t (•)) to achieve approximated parallelized iterations for FOO as illustrated in our Fig. 5 , which is more practical and useful.𝜃 !"# Process 1 𝜃 ! 𝜃 !,# ← 𝜃 !$% 𝜃 !,% ⋯ 𝜃 !,& 𝜃 !,'$& 𝜃 !,'$% 𝜃 ! (%) 𝜃 ! (&) 𝜃 ! (*) 𝜃 ! ('$%) 𝜃 ! (')Let input θ = [θ i ] d i=1 , the Ackley, Sphere, and Rosenbrock functions applied in our synthetic experiments are given below, which have been slightly modified compared with the standard ones.F (θ) = -20 exp   -0.2 1 d d i=1 θ 2 i   -exp( 1 d d i=1cos (2πθ i )) + 20 + exp(1), (Ackley)F (θ) = 1 d d i=1θ 2 i , (Sphere)EQUATIONNote that both Ackley and Sphere function achieve their minimum (i.e., min F (θ) = 0) at θ * = 0, whereas Rosenbrock function achieves its minimum (i.e., min F (θ) = 0) at θ * = 1.In this experiment, the parallelism of N = 5 is applied and all the baselines introduced in Sec. 6.1 as well as our OptEx are based on Adam [4] with a learning rate of 0.1, β 1 = 0.9, and β 2 = 0.999.In addition, we employ a Matérn kernel-based gradient estimation in our OptEx with T 0 = 20.Our experimental framework is built on the Deep Q-Network (DQN) algorithm, as outlined in [39] , and implemented within the OpenAI Gym environment [38] . This study investigates the effectiveness of different optimizer configurations across classical discrete control tasks provided by Gym. Each trial is conducted on a dedicated CPU to maintain consistency in computational conditions. The DQN architecture consists of dual fully connected layers, with 64 or 128 neurons tailored to each task's requirements. Hyperparameters, including a learning rate of 0.001, a reward discount factor of 0.95, and a batch size of 256, are applied for fairness and consistency across experiments.Performance evaluation of the optimizer-enhanced DQN agents is systematically carried out over 100 to 200 episodes per game, employing an ϵ-greedy policy with a minimum epsilon of 0.1 and an exponential epsilon decay with a rate of 2 -1 1500 . A preliminary warm-up phase of either 30 or 50 episodes, depending on the task, is incorporated to stabilize initial learning dynamics. Besides, all baselines introduced in Sec.6.1 and our OptEx are based on Adam [4] with a learning rate of 0.001, β 1 = 0.9, and β 2 = 0.999. For OptEx, we utilize a Mat'ern kernel-based gradient estimation, with T 0 = 150 to accommodate the variance in RL tasks, and parallelism of N = 4 is applied.In this experiment, we compared our OptEx with other baselines using both image classification and text autoregression tasks. Here, we simply make use of the jax.vmap function to simulate parallel computing and measure the wallclock time for each sequential iteration. We believe that the time efficiency of our OptEx can be further improved when it is more properly implemented on a parallel computing platform. Besides, to reduce the computational cost of our kernelized gradient estimation in these high-dimensional optimization problems, we propose to use a randomly sampled subset of dimensions (e.g., d = 10 4 for image classification and d = 10 5 for text autoregression) from the total d dimension to compute the kernel value k(•, •) in each sequential iteration of our OptEx.Image Classification. In this image classification task, we train a 9-layer MLP (including input and output layer) with skip connections on MNIST [45] , Fashion MNIST [46] and a 10-layer MLP (including input and output layer) with skip connections on CIFAR-10 [41] datasets, which have a parameter size of d = 978186 for (fashion-)MNIST and d = 2412298 for CIFAR-10. Both our OptEx and other baselines are based on SGD [1] with a learning rate of 0.001, a batch size of 512, and parallelism of N = 4. For OptEx, we employ a Matérn kernel-based gradient estimation with T 0 = 6.Text Autoregression. In addition, we further train a simple transformer from Haiku library [42] with a parameter size of d = 1626496 on the corpus of "Harry Potter and the Sorcerers Stone" and a subset work from Shakespeare. In both tasks, all the baselines introduced in Sec. 6.1 and our OptEx are based on SGD [1] with a learning rate of 0.01, batch size of 256 and parallelism of N = 4. For OptEx, we employ a Matérn kernel-based gradient estimation, where T 0 = 10.Ablation Studies on Synthetic Function. To better understand our OptEx algorithm, we have conducted a number of ablation studies on the Rosenbrock synthetic function with a dimension of d = 10 5 . The results are in Fig. 6 , in which there are 4 different types of comparisons: (a) We have compared our OptEx with vs. without evaluating the intermediate gradients, i.e., {∇f (θ t,i-1 )} N -1 i=1 at every iteration t, denoting as parallel and sequential respectively in Fig. 6 (a), which aims to • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors). • It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text. 


Decision: The paper presents OptEx, a first-order optimization framework that enhances the efficiency of first-order optimization algorithms (FOO) through approximately parallelized iterations. OptEx mitigates FOO's reliance on numerous sequential iterations by using kernelized gradient estimation to predict required gradients based on historical data, allowing for approximate parallelization.

The paper establishes theoretical guarantees for the estimation error of the kernelized gradient method and the iteration complexity of SGD-based OptEx, showing that the estimation error approaches zero as gradient history accumulates. Numerical examples are provided to showcase the performance of the method in various settings.

Overall, most reviewers find the paper’s contributions to be novel, interesting, and valuable. Some questions, concerns, and suggestions have arisen during the review process, which the authors have adequately addressed. They are strongly encouraged to incorporate these into the final version of the paper.

Reviews: review: paper_summary: paper_summary: This paper introduces a new approach for parallelizing stochastic gradient descent for unconstrained, nonconvex, smooth stochastic optimization. The approach is based on building a Gaussian Process surrogate for the true gradient (based on a history of observed stochastic gradients), and uses this surrogate to suggest potential 'future' iterates, by iterating gradient descent using the surrogate. Iterates are then adjusted based on new stochastic gradient estimates, computed in parallel.

They present worst-case complexity theory that says that being able to compute $N$ stochastic gradients in parallel improves the iteration complexity by a factor of $\sqrt{N}$, with numerical results on a selection of synthetic problems, reinforcement learning problems, and neural network training problems.
strength_weakness: strength_weakness: 1. There are gaps in the literature survey of parallelized SGD methods. For example, there is no mention of the highly cited and relevant papers:
- Recht et al. Hogwild!: A Lock-Free Approach to Parallelizing Stochastic Gradient Descent. NeurIPS 2011
- Mnih et al. Asynchronous Methods for Deep Reinforcement Learning. ICML 2016
- Yu et al. Parallel Restarted SGD with Faster Convergence and Less Communication. AAAI, 2019

2. It is unclear to me whether the convergence theory (Theorem 2) actually supports the assertion of a $\sqrt{N}$ speedup from $N$ processes. The RHS gradient bound indeed is $O(1/\sqrt{N})$, but you are measuring the smallest gradient over $O(N)$ outer iterations. So, doesn't this say that you decrease your gradient by a factor of $O(1/\sqrt{N})$ after doing $N$ times extra work, which is exactly the convergence rate of SGD under your sub-Gaussian stochastic gradient assumption [33]? The work would also benefit from a direct comparison with just sample averaging of stochastic gradients (see question below).

3. It appears to me that Theorem 3 appears to contradict standard theory (e.g. Bottou, Curtis & Nocedal. Optimization Methods for Large-Scale Machine Learning. SIAM Review, 2018). You claim that the complexity bound for (nonconvex) SGD is tight when applied to the spherical function eq (47), which is strongly convex. However (stochastic) gradient descent with properly chosen stepsizes should converge linearly for this problem, not sublinearly as claimed. Can you explain this potential mismatch in results?

4. For very large problems, it is unrealistic that parallelism will be available, since storing the model and training data may already require parallelism (i.e. a single stochastic gradient evaluation is already computed in parallel). It would be helpful to include a comment about the specific regimes (dimension, dataset size, level of parallelism) where this method is intended to be used.
questions: questions: How is the next iterate of OptEx defined? Line 10 of Algorithm 1 says it is the final iterate that is chosen, but Figure 1 suggests that the next iterate is somehow taken from a combination of the $\theta_t^1, \ldots, \theta_t^N$ values. Can you please clarify? The latter interpretation would make more sense, otherwise all but one of the parallel evaluations are discarded immediately.

Is the $O(1/\sqrt{N})$ speedup exactly the same as what you would get by just doing basic sample averaging? That is, at each iteration, just get $N$ different stochastic gradient estimates at the same point and average them (reducing the variance of the estimates). How does your bound compare with this simple use of parallelism?

You say that "$\nabla F$ is assumed to be sampled from the Gaussian Process, $\nabla F \sim GP(0, K(\cdot,\cdot))$". Can you please clarify? I read this to mean that this is the assumed prior of your GP surrogate for grad F, but as written it suggests that your actual objective function has - on average - zero gradient, meaning we would expect the 'average point' to be almost stationary?

You say the "Target baseline is equivalent to Algo. 1 with [GP mean] being replaced with [stochastic gradient], indicating the desired parallelized iteration we aim to approximate". Does the 'target' implementation also have the parallel step? If so, we would expect that to give an unfair advantage to 'target', since it does $2N$ stochastic gradient evaluations per outer iteration rather than $N$ for OptEx.
limitations: limitations: Limitations appear to have been addressed. The theoretical assumptions (e.g. smoothness) are clearly articulated and some practical limitations are addressed in Section 7.

No potential negative societal impact of this work.
paper_summary: paper_summary: This paper presents OptEx, an approach to parallelize optimization methods by using gradients from previous iterations to predict gradient for subsequent iterations which in turn breaks the serial nature of standard stochastic optimization thereby enabling approximately parallel iterations. The gradient prediction is done using Kernel methods. The paper also presents experimental evidence towards the effectiveness of the proposed approach.
strength_weakness: strength_weakness: -- The question of gradient estimation converging to the true gradient seems fairly far-fetched, and there is every reason to believe modeling gradients would essentially carry a constant bias without very strong regularity and realizability conditions. 

-- Furthermore, when applied to problems in RL, this becomes even more confounding because of the necessity to explore and using data collected with exploration to in turn construct gradient estimates. I am not sure if I follow how this can be achieved within this framework.

-- The paper's experimental comparison doesn't offer any sources of comparisons against the many different approaches to parallelization that have already been studied in the literature, e.g. mini-batching, model averaging etc.
questions: questions: Can the authors clarify the speedup offered by utilizing their framework while using acceleration/momentum based methods? Standard mini-batch methods offer different trade-offs particularly wrt batch sizes used while still achieving linear parallelization speedups [1].

[1] Cotter et al: Better Mini-Batch Algorithms via Accelerated Gradient Methods.
limitations: limitations: N/A.
paper_summary: paper_summary: The paper "OptEx: Expediting First-Order Optimization with Approximately Parallelized Iterations" introduces a novel framework, OptEx, designed to improve the efficiency of first-order optimization (FOO) algorithms by approximately parallelizing their iterations.
strength_weakness: strength_weakness: 1. The kernelized gradient estimation and the associated computational techniques are complex, which might pose challenges for practical implementation and understanding by a broader audience.
2. The efficiency and scalability of OptEx in very large-scale optimization problems or in highly dynamic environments might be limited
3. The theoretical guarantees rely on certain assumptions (e.g., the Gaussian distribution of gradients), which may not hold in all practical scenarios, potentially limiting the generality of the results.
questions: questions: see weakness
limitations: limitations: NA
paper_summary: paper_summary: The paper proposes a general framework for accelerating first-order optimization methods. The framework leverages parallel computing by estimating the gradients using kernel methods and breaking iterative dependencies. The paper establishes theoretical guarantees that the method gives an acceleration rate of $\sqrt{N}$ for $N$ parallelism. Then the paper does extensive empirical studies to show the improvements of this framework.
strength_weakness: strength_weakness: The theoretical guarantees rely on assumptions that might not work for $N$ very large, thus limit the extent to which the framework can accelerate optimization.

In the empirical part, $N$ are not very large, thus unable to show that the framework can accelerate things to great extent in general.
questions: questions: Suppose that the step size is very small, can the $N$ potentially be very large?
limitations: limitations: The authors adequately addressed the limitations.

Comments: title: response to author
Very Good Rebuttal
comment: Dear Reviewer ieTR,

Thank you so much for taking the time to review our paper and for your insightful questions. While we have thoroughly addressed the concerns and questions raised by all the other reviewers, we sincerely hope that our clarifications and discussions above have effectively addressed your concerns as well.

If you have any more questions or need more details, we are happy to answer them promptly within the discussion period.

Best, Authors
Dear Reviewer 1m6n,

Thank you so much for your positive response after our rebuttal! We are so happy to hear that your concerns have been well addressed. We will incorporate all the discussions above into our revised version.

Sincerely, Authors
Thanks for the detailed response from the author, my concern has already been addressed.
Dear Reviewer RYzW,

Thank you so much for your prompt and thoughtful response. We are pleased to hear that our clarifications have improved your opinion of our work, and we will incorporate all the clarifications and discussions mentioned above into our revised manuscript, as per your suggestion.

We also appreciate your valuable insight regarding the parallel computing required by combining sample averaging and OptEx, in relation to your question #2. We agree with you that in this scenario, we will need $N=N_1 \times N_2$ parallel computations, leading to a speedup of $O(1/\sqrt{N})$. The benefits of our OptEx over sample averaging may be more straightforward by looking into our Equation 45, which mirrors the convergence results of sample averaging. Particularly, when the gradient noise $\sigma$ is already decreased to a small value, i.e., the second term on the RHS in Equation 45 is small, the first term on the RHS in Equation 45 that is related to sequential iterations will then dominate the convergence. Consequently, simply applying sample averaging with increased parallel computational resources yields only marginal improvements in convergence. In contrast, OptEx can achieve a noticeable improvement by reducing the first term on the RHS in Equation 45 by $O(1/N)$. We will also include this discussion in our revised paper. Should you have any other questions, we would like to address them.

Sincerely, Authors
Thanks for this information. Regarding question #2, I think this is a very important aspect and warrants a discussion in any new version of the paper. However, I don't think the benefit is as good as suggested: if you want to do $N_1$ parallel samples for each stochastic gradient (to reduce the variance to $\sigma^2/N_1$), then if you run OptEx with $N_2$ parallel runs, then you will need $N_1\times N_2$ stochastic gradients in parallel at each outer loop, not $N_1+N_2$, so your speedup is still $O(1/sqrt{N})$ where $N=N_1\times N_2$ is the number of cores needed.

This has been a helpful discussion, and pending the clarifications in a revised version as agreed by the authors above, I am happy to upgrade my scores of the paper in the initial review above.
Dear Reviewer RYzW,

Thank you so much for your positive feedback and thoughtful suggestions! We sincerely appreciate your recognition and support. We are glad to hear that our clarifications are very helpful. We will add our clarifications above in the writing to make our paper clearer. We would like to address your remaining concerns below.

---

### Questions #2

We totally agree with you that sample averaging applied to standard SGD can achieve similar complexity as our OptEx, without the overhead of managing the GP. As we have justified above, these two speedup approaches stem from different principles:

- OptEx Speedup: This comes from reduced sequential iterations (first term on the RHS in Equation 45).
- Sample Averaging Speedup: This is derived from reduced gradient variance (second term on the RHS in Equation 45).

So, when gradient noise $\sigma$ is already decreased to a small value through sample averaging with parallelism $N'$ (not $N$), increasing $N'$ further (i.e., $N'+N$) provides diminishing returns due to the sublinear rate of $O(1/\sqrt{N'+N})$. For example, the benefit of reducing from $1/\sqrt{10}$ to $1/\sqrt{15}$ is less significant compared to reducing from $1/\sqrt{1}$ to $1/\sqrt{6}$ with an additional parallelism $N=5$.

However, in this scenario of a small gradient noise $\sigma$, our OptEx can still achieve noticeable speedup (e.g., from $1/\sqrt{1}$ to $1/\sqrt{6}$ with the same additional parallelism $N=5$) by leveraging additional parallel computing to reduce the sequential iterations. This reduction in sequential iterations is however unattainable by mini-batch SGD itself.

In light of this, Equation 45, leading to our Theorem 2, in fact, will provide valuable insights into when OptEx is preferable to simple sample averaging. Overall, OptEx complements rather than replaces existing parallelization methods, including sample averaging, to accelerate first-order optimization. This is particularly effective when other methods are **inapplicable** or **underperforming**, i.e., further noticeable improvements cannot be achieved merely by increasing parallel computing resources in those methods, as we have justified above.

### Weaknesses #2
We apologize for the confusion caused by our notation $\tau$. Our intention was to use $\tau$ to denote all gradients evaluated during the optimization process more easily, including those in parallel processes, as mentioned in line 260. In Theorems 2 and 3, the expression $\min_{\tau \leq NT} \left|\nabla F(\theta_{\tau})\right|$ is exactly meant to represent $\min_{t \leq T, s \leq N} \left|\nabla F(\theta_{t,s})\right|$. Thank you for bringing this to our attention. We will correct this notation in our revised paper.

---

We want to thank Reviewer RYzW for your positive feedback and recognition again! We sincerely hope our response has addressed your remaining concerns and can increase your opinions of our work. If you have any other questions, we would like to address them.

Sincerely, Authors
Thank you for these clarifications, they are very helpful.

Pending seeing a new version of the manuscript, I think the responses given for all points above would be very helpful. I have one remaining concern (question Q2) and one point for clarification (weakness Q2) that I still think warrant more consideration.

**Questions #2**

I understand that you would like to potentially use OptEx in the exact gradient case ($\sigma=0$, e.g. Section 6.1), and your results make sense in this context. However, in the SGD context that seems to be your main focus, if I had the ability to compute $N$ stochastic gradients in parallel (each with a given noise level $\sigma$ as you assume), I could either use OptEx or compute $N$ stochastic gradients at the current iterate and average these, giving an effective variance of $\sigma^2/N$. 

In Theorem 1, your complexity bound (when $\sigma>0$) is $O(\sigma/\sqrt{NT})$ for $T$ iterations, which suggests that sample averaging applied to standard SGD (where $\sigma \mapsto \sigma/\sqrt{N}$) would achieve the same complexity as OptEx with essentially no overhead from managing the GP, etc. 

It would be useful if you could either discuss why/when OptEx should be preferred to simple sample averaging, as I'm not sure if your theoretical results show this. See, for example, the discussion "Trade-Offs of (Mini-)Batching" in Section 4 of Bottou, Curtis & Nocedal [your ref 26].

**Weaknesses #2**

I think my confusion arose because your statement of Theorems 2 and 3 are not quite correct. Specifically, the LHS term in eqs (8) and (9) should be $\min_{\tau\leq T, s\leq N} \|\nabla F(\theta_{\tau,s})\|^2$ rather than $\min_{\tau\leq NT} \|\nabla F(\theta_{\tau})\|^2$ as given (i.e. reflecting that you do $T$ outer iterations with $N$ parallel evaluations in each, rather than $NT$ outer iterations with $N$ parallel evaluations in each, as currently stated). This seems to align with what is shown in the proofs.
Dear Reviewer YWBT,

Thank you so much for your prompt and positive response after our rebuttal! We are so happy to hear that your concerns have been well addressed. We will incorporate the discussions above into our revised version.

Sincerely, Authors
Thanks the authors for their well-written rebuttal. The response on my confusion about N is very clear. Thanks a lot.