High-dimensional Problems Research Articles

Accelerated optimization through time-scale analysis of inertial dynamics with asymptotic vanishing and Hessian-driven dampings

ABSTRACT Gradient optimization algorithms can be effectively studied from the perspective of ordinary differential equations (ODEs), where these algorithms are obtained through the temporal discretization of the continuous dynamics. This approach provides a powerful tool to comprehend acceleration phenomena in optimization, thanks to time-scale techniques and Lyapunov analysis. In the context of a Hilbert setting, our study focuses on the rapid optimization properties of inertial dynamics, which combine asymptotically vanishing viscous damping with Hessian-driven damping. These dynamics are derived as high-resolution ODEs from both the Nesterov accelerated gradient method and the Ravine method. For a general differentiable convex function f, choosing the viscous damping coefficient in the form of α / t with 3 $ ]]> α > 3 guarantees an inverse quadratic convergence rate of the values, i.e. o ( 1 / t 2 ) , the weak convergence of trajectories towards optimal solutions, and the fast convergence of gradients towards zero. By carefully scaling the dynamics temporally and based on the tuning of the damping coefficient in front of the Hessian term, we identify the limit dynamic as α becomes large. In particular, we examine the case where the limit dynamic corresponds to the Levenberg-Marquardt regularization of Newton's continuous method. This explanation accounts for the aforementioned fast convergence properties and provides fresh insights into the complexity of these methods, which play a central role in optimization for high-dimensional problems. The interplay between numerical algorithms and continuous-time ODEs plays a fundamental role in our analysis for the design and understanding of accelerated optimization algorithms. Numerical experiments are presented to illustrate and support the theoretical results.

Alternating minimization for generalized rank-1 matrix sensing: sharp predictions from a random initialization

Abstract We consider the problem of estimating the factors of a rank-$1$ matrix with i.i.d. Gaussian, rank-$1$ measurements that are nonlinearly transformed and corrupted by noise. Considering two prototypical choices for the nonlinearity, we study the convergence properties of a natural alternating update rule for this non-convex optimization problem starting from a random initialization. We show sharp convergence guarantees for a sample-split version of the algorithm by deriving a deterministic one-step recursion that is accurate even in high-dimensional problems. Notably, while the infinite-sample population update is uninformative and suggests exact recovery in a single step, the algorithm—and our deterministic one-step prediction—converges geometrically fast from a random initialization. Our sharp, non-asymptotic analysis also exposes several other fine-grained properties of this problem, including how the nonlinearity and noise level affect convergence behaviour. On a technical level, our results are enabled by showing that the empirical error recursion can be predicted by our deterministic one-step updates within fluctuations of the order $n^{-1/2}$ when each iteration is run with $n$ observations. Our technique leverages leave-one-out tools originating in the literature on high-dimensional $M$-estimation and provides an avenue for sharply analyzing complex iterative algorithms from a random initialization in other high-dimensional optimization problems with random data.

Active learning for adaptive surrogate model improvement in high-dimensional problems

This paper investigates a novel approach to efficiently construct and improve surrogate models in problems with high-dimensional input and output. In this approach, the principal components and corresponding features of the high-dimensional output are first identified. For each feature, the active subspace technique is used to identify a corresponding low-dimensional subspace of the input domain; then a surrogate model is built for each feature in its corresponding active subspace. A low-dimensional adaptive learning strategy is proposed to identify training samples to improve the surrogate model. In contrast to existing adaptive learning methods that focus on a scalar output or a small number of outputs, this paper addresses adaptive learning with high-dimensional input and output, with a novel learning function that balances exploration and exploitation, i.e., considering unexplored regions and high-error regions, respectively. The adaptive learning is in terms of the active variables in the low-dimensional space, and the newly added training samples can be easily mapped back to the original space for running the expensive physics model. The proposed method is demonstrated for the numerical simulation of an additive manufacturing part, with a high-dimensional field output quantity of interest (residual stress) in the component that has spatial variability due to the stochastic nature of multiple input variables (including process variables and material properties). Various factors in the adaptive learning process are investigated, including the number of training samples, range and distribution of the adaptive training samples, contributions of various errors, and the importance of exploration versus exploitation in the learning function.

A scalable adaptive sampling approach for surrogate modeling of rigid pavements using machine learning

Rigid pavement design is a high-dimensional optimization problem, involving several variables and design considerations. The existing machine learning (ML) design models are either low-dimensional or less accurate than computational simulations. Surrogate modeling is a powerful tool for approximating the results of high-fidelity computational simulations, which is commonly used to approximate the results of these simulations by sampling their solutions to train models. However, conventional and adaptive sampling methods face the challenge of “curse of dimensionality”, due to the exponential increase in required sample size and computations with the number of variables. To overcome this, we propose a scalable adaptive sampling (SAS) method that uses random samples for testing the model's performance but generates new training samples as a subset of a full factorial design of experiments (DoE). The factorial level increases with each iteration, allowing the algorithm to sample training data at a progressively finer scale, and updating the ML models each time adaptively. The proposed method was tested on several 2D benchmark examples, as well as practical pavement design problems. Our results demonstrate that the proposed technique can develop accurate surrogate models for both low- and high-dimensional inference spaces. For a 4D inference space, the surrogate models derived from the proposed method had an order of magnitude lower error than those derived from “one-shot” conventional sampling with the same sample size. For a 6D inference space, the sample size required by the proposed method was only 5 % of that required by conventional sampling for comparable performance.

Embedded exponential Runge–Kutta–Nyström methods for highly oscillatory Hamiltonian systems

Exponential/extended Runge–Kutta–Nyström (ERKN) methods have been validated by numerous theoretical and numerical results to be more efficient and suitable than classical Runge–Kutta–Nyström (RKN) methods in dealing with highly oscillatory Hamiltonian systems. Once numerical solutions are required to achieve a prescribed local error tolerance for practical problems, the embedded ERKN pairs with adaptive stepsize are needed. Given the higher efficiency of high-order methods than low-order methods, this paper focuses on high-order embedded ERKN pairs. Using the particular mapping from RKN methods into ERKN methods, we establish an approach to constructing high-order embedded ERKN pairs. By diagonalizing the frequency matrix, an implementation algorithm with less calculation amount than the direct calculation procedure is proposed for high-dimensional problems. In addition, the dispersion and dissipation of the ERKN pairs ERKN4(3) and ERKN8(6) are analyzed. Furthermore, the application of ERKN pairs to an important highly oscillatory problem, i.e., the wave equation prescribed with different boundary conditions is discussed. For the periodic boundary condition, a special implementation algorithm incorporated with FFT techniques is presented, which decreases the calculation amount in one time step from O(N2) to O(Nlog⁡N). Finally, numerical results with the FPU problem, the Klein–Gordon equation in the nonrelativistic limit regime, and a two-dimensional wave equation demonstrate the superiority of ERKN pairs over classical RKN pairs.

StocIPNet: A novel probabilistic interpretable network with affine-embedded reparameterization layer for high-dimensional stochastic inverse problems

The stochastic inverse problem (StocIP), which aims to align push-forward and observed output distributions by estimating probability distributions of unknown system inputs, often faces optimization challenges and the curse of dimensionality. A novel deep network called StocIPNet which comprises an affine-embedded reparameterization subnetwork (ReparNet) and a complex system metamodeling subnetwork (MetaNet) is proposed to alleviate these issues. The ReparNet subnetwork embeds the affine transformation to convert the statistical parameters of the physical random vector into learnable weights and biases, effectively implementing the reparameterization trick by separating random and deterministic elements in the stochastic sampling operation to preserve differentiability. In parallel, the MetaNet subnetwork offers a computationally efficient alternative to time-consuming forward solvers, facilitating the generation of push-forward distributions. The entire StocIPNet utilizes the kernel maximum mean discrepancy (MMD) as a distribution-free loss function, quantifying the discrepancy between push-forward and observed output distributions. By leveraging the ReparNet’s advantage of reformulating the sampling process as a differentiable transformation and combining two subnetworks seamlessly, the StocIP is reconfigured into a pure network training paradigm preserving differentiability perfectly, which allows for direct modeling and efficient inference of uncertainty within the network using automatic differentiation, backpropagation and gradient-based optimization methods, enabling ease of scaling to high-dimensional problems. The proposed framework has been theoretically demonstrated to be equivalent to the maximum likelihood method, ensuring its solid probabilistic interpretable foundation. The proposed framework is applied to perform stochastic model updating on a numerical and an experimental structure, which effectively demonstrates the framework’s remarkable effectiveness and high efficiency in treating high-dimensional StocIPs.

A surrogate-assisted evolutionary algorithm with dual restricted Boltzmann machines and reinforcement learning-based adaptive strategy selection

To improve the effectiveness of surrogate-assisted evolutionary algorithms (SAEAs) in solving high-dimensional expensive optimization problems with multi-polar and multi-variable coupling properties, a new approach called DRBM-ASRL is proposed. This approach leverages restricted Boltzmann machines (RBMs) for feature learning and reinforcement learning for adaptive strategy selection. DRBM-ASRL integrates four search strategies based on three heterogeneous surrogate modeling approaches, each catering to different preferences. Two of these strategies focus on generative sampling in the subspaces with varying dimensions, while the other two aim to explore the local and global landscapes in the high-dimensional source space. This allows for more effective tradeoffs between exploration and exploitation in the solution space. Reinforcement learning is employed to adaptively prioritize the search strategies during optimization , based on the online feedback information from the optimal solution. In addition, to enhance the representation of potentially optimal samples in the solution space, two task-driven RBMs are separately trained to construct a feature subspace and reconstruct the features of the source space. DRBM-ASRL has been evaluated on various high-dimensional benchmarks ranging from 50 to 200 dimensions, as well as 14 CEC 2013 complex benchmark problems with 100 dimensions and a power system problem with 118 dimensions. Experimental results demonstrate its superior convergence performance and optimization efficiency compared to eight state-of-the-art SAEAs.

A sequence decomposition genetic algorithm for petrophysical inversion with spatial correlation

Predicting petrophysical properties based on seismic attributes is essential for accurate subsurface characterization. Typically, a forward model is first established to relate the petrophysical properties to geophysical measurements, and then the petrophysical properties are inverted from the actual measurements. Imposing prior constraints of the spatial relationship is also desirable to regularize the ill-posed inverse problem. However, petrophysical inversion involving spatial correlation presents a high-dimensional optimization problem with multiple local extrema. This problem can be challenging to solve, particularly when discrete variables (such as lithofacies) and continuous variables (such as porosity and fluid saturation) are considered together. Therefore, we develop a global optimization algorithm that uses the Markov chain decomposition and a genetic algorithm for a joint discrete-continuous petrophysical inversion with spatial correlation. The innovation of this algorithm is to spatially decouple the sequence optimization problem into multiple low-dimensional subproblems and overcome the initial value dependence to converge to the global optimal solution. Furthermore, an approximation operation is introduced to expedite the optimization algorithm. Numerical experiments show that after algorithmic acceleration, the spatially correlated inversion requires the same computational time as a simple spatially uncorrelated inversion and achieves significantly better accuracy than the latter. Finally, we apply the method to a case study in East China, wherein the lithofacies, clay content, porosity, and water saturation are inferred from the data of the wave velocities and density. The inversion results around the wells generally match the well-logging data, verifying the feasibility of the new method.

Using group level factor models to resolve high dimensionality in model-based sampling.

Joint modeling of decisions and neural activation poses the potential to provide significant advances in linking brain and behavior. However, methods of joint modeling have been limited by difficulties in estimation, often due to high dimensionality and simultaneous estimation challenges. In the current article, we propose a method of model estimation that draws on state-of-the-art Bayesian hierarchical modeling techniques and uses factor analysis as a means of dimensionality reduction and inference at the group level. This hierarchical factor approach can adopt any model for the individual and distill the relationships of its parameters across individuals through a factor structure. We demonstrate the significant dimensionality reduction gained by factor analysis and good parameter recovery, and illustrate a variety of factor loading constraints that can be used for different purposes and research questions, as well as three applications of the method to previously analyzed data. We conclude that this method provides a flexible and usable approach with interpretable outcomes that are primarily data-driven, in contrast to the largely hypothesis-driven methods often used in joint modeling. Although we focus on joint modeling methods, this model-based estimation approach could be used for any high dimensional modeling problem. We provide open-source code and accompanying tutorial documentation to make the method accessible to any researchers. (PsycInfo Database Record (c) 2024 APA, all rights reserved).

A novel hybrid Artificial Gorilla Troops Optimizer with Honey Badger Algorithm for solving cloud scheduling problem

Cloud computing has revolutionized the way a variety of ubiquitous computing resources are provided to users with ease and on a pay-per-usage basis. Task scheduling problem is an important challenge, which involves assigning resources to users’ Bag-of-Tasks applications in a way that maximizes either system provider or user performance or both. With the increase in system size and the number of applications, the Bag-of-Tasks scheduling (BoTS) problem becomes more complex due to the expansion of search space. Such a problem falls in the category of NP-hard optimization challenges, which are often effectively tackled by metaheuristics. However, standalone metaheuristics generally suffer from certain deficiencies which affect their searching efficiency resulting in deteriorated final performance. This paper aims to introduce an optimal hybrid metaheuristic algorithm by leveraging the strengths of both the Artificial Gorilla Troops Optimizer (GTO) and the Honey Badger Algorithm (HBA) to find an approximate scheduling solution for the BoTS problem. While the original GTO has demonstrated effectiveness since its inception, it possesses limitations, particularly in addressing composite and high-dimensional problems. To address these limitations, this paper proposes a novel approach by introducing a new updating equation inspired by the HBA, specifically designed to enhance the exploitation phase of the algorithm. Through this integration, the goal is to overcome the drawbacks of the GTO and improve its performance in solving complex optimization problems. The initial performance of the GTOHBA algorithm tested on standard CEC2017 and CEC2022 benchmarks shows significant performance improvement over the baseline metaheuristics. Later on, we applied the proposed GTOHBA on the BoTS problem using standard parallel workloads (CEA-Curie and HPC2N) to optimize makespan and energy objectives. The obtained outcomes of the proposed GTOHBA are compared to the scheduling techniques based on well-known metaheuristics under the same experimental conditions using standard statistical measures and box plots. In the case of CEA-Curie workloads, the GTOHBA produced makespan and energy consumption reduction in the range of 8.12–22.76% and 6.2–18.00%, respectively over the compared metaheuristics. Whereas for the HPC2N workloads, GTOHBA achieved 8.46–30.97% makespan reduction and 8.51–33.41% energy consumption reduction against the tested metaheuristics. In conclusion, the proposed hybrid metaheuristic algorithm provides a promising solution to the BoTS problem, that can enhance the performance and efficiency of cloud computing systems.

Random projection enhancement: A Novel method for improving performance of surrogate models

Surrogate models serve as powerful tools for engineering application but constructing highly accurate models with limited samples remains a challenge. Inspired by this, this paper proposes a general random projection enhancement method (RPE) aimed at improving surrogate model performance. Taking least squares support vector regression (LSSVR) as an example, we conducted numerical experiments to evaluate the feasibility of the RPE method. The RPE-enhanced LSSVR demonstrates superior predictive accuracy and robustness compared to the unenhanced model and other benchmark models. The results show that the RPE method can not only enhance the predictive accuracy of the model but also improve the robustness of the model. Meanwhile, analysis of optimization experiments indicated that the enhanced LSSVR model obtained better solutions compared to the unenhanced LSSVR model. Even when extended to high-dimensional problems, the RPE method remains effective. Furthermore, the applicability of RPE method extends from the LSSVR model to other models, as demonstrated by the enhancement of the extreme learning machine, regression kriging, and radial basis function neural network. More importantly, the feasibility of the RPE method is also verified in engineering problems, obtaining the same results. This method offers a reliable alternative for real-world engineering optimization problems.

Taking another step: A simple approach to high-dimensional Bayesian optimization

Scaling Bayesian optimization to high-dimensional problems is a meaningful but challenging task. Most of current approaches assume only a few variables are effective to the objective function or the objective function is additively separable, therefore are not suitable when the problem violates these assumptions. In this work, we propose a simple and efficient approach to extend Bayesian optimization to high dimensions. The proposed approach does not make these two assumptions. In each iteration, after locating a candidate point by the global Gaussian process model, we train a local Gaussian process model around the candidate point and locate a new point using the local model. Instead of evaluating the solution located by the global model, we evaluate the solution located by the local model. This simple taking-another-step approach is shown to be able to improve Bayesian optimization's performance significantly on high-dimensional optimization problems. The proposed algorithm also shows competitive performance when compared with four high-dimensional Bayesian optimization algorithms and four surrogate-assisted evolutionary algorithms.

But how can I optimise my high-dimensional problem with only very little data? – A composite manufacturing application

In this research, a Gaussian process (GP) surrogate modelling framework for the forming process of dry carbon-fibre textile was investigated. A particular focus of the work is the development of dimension reduction algorithms, allowing to solve high-dimensional sparse optimisation problems. The concept of active subspace is adopted to find the principal space of the problem. Then, a low-dimensional (i.e., active) subspace can be obtained by selecting the directions with highest explained variance. A kernel-combined GP format is developed. This takes advantage of the active subspace to build a robust, high-dimensional emulator that can be regarded as a special case of multi-fidelity GP. A two-step adaptive sequential design approach is adopted, which further improves the efficiency of data design. Different sequential design strategies are compared. A case study with eight input parameters demonstrates the capability of the proposed approach, where an accurate and robust optimum condition is obtained from only tens of simulations.

A score-based filter for nonlinear data assimilation

We propose a score-based generative sampling method for solving the nonlinear filtering problem with superior accuracy. A major drawback of existing nonlinear filtering methods, e.g., particle filters, is the low accuracy in handling high-dimensional nonlinear problems. To overcome this issue, we incorporate the score-based diffusion model into the recursive Bayesian filter framework to develop a novel score-based filter (SF). The key idea of SF is to store the information of the recursively updated filtering density function in the score function, instead of storing the information in a set of finite Monte Carlo samples (used in particle filters and ensemble Kalman filters). By leveraging the reverse-time diffusion process, SF can generate unlimited samples to characterize the filtering density. An essential aspect of SF is its analytical update step, gradually incorporating data information into the score function. This step is crucial in mitigating the degeneracy issue faced when dealing with very high-dimensional nonlinear filtering problems. Three benchmark problems are used to demonstrate the performance of our method. In particular, SF provides surprisingly impressive performance in reliably capturing/tracking the 100-dimensional stochastic Lorenz system that is a well-known challenging problem for existing filtering methods.

A naïve Bayes regularized logistic regression estimator for low-dimensional classification

To reduce the estimator's variance and prevent overfitting, regularization techniques have attracted great interest from the statistics and machine learning communities. Most existing regularized methods rely on the sparsity assumption that a model with fewer parameters predicts better than one with many parameters. This assumption works particularly well in high-dimensional problems. However, the sparsity assumption may not be necessary when the number of predictors is relatively small compared to the number of training instances. This paper argues that shrinking the coefficients towards a low-variance data-driven estimate could be a better regularization strategy for such situations. For low-dimensional classification problems, we propose a naïve Bayes regularized logistic regression (NBRLR) that shrinks the logistic regression coefficients toward the naïve Bayes estimate to provide a reduction in variance. Our approach is primarily motivated by the fact that naïve Bayes is functionally equivalent to logistic regression if naïve Bayes' conditional independence assumption holds. Under standard conditions, we prove the consistency of the NBRLR estimator. Extensive simulation and empirical experimental results show that NBRLR is a competitive alternative to various state-of-the-art classifiers, especially on low-dimensional datasets.

Inf-sup neural networks for high-dimensional elliptic PDE problems

Solving high dimensional partial differential equations (PDEs) has historically posed a considerable challenge when utilizing conventional numerical methods, such as those involving domain meshes. Recent advancements in the field have seen the emergence of neural PDE solvers, leveraging deep networks to effectively tackle high dimensional PDE problems. This study introduces Inf-SupNet, a model-based unsupervised learning approach designed to acquire solutions for a specific category of elliptic PDEs. The fundamental concept behind Inf-SupNet involves incorporating the inf-sup formulation of the underlying PDE into the loss function. The analysis reveals that the global solution error can be bounded by the sum of three distinct errors: the numerical integration error, the duality gap of the loss function (training error), and the neural network approximation error for functions within Sobolev spaces. To validate the efficacy of the proposed method, numerical experiments conducted in high dimensions demonstrate its stability and accuracy across various boundary conditions, as well as for both semi-linear and nonlinear PDEs.

Towards multi-objective high-dimensional feature selection via evolutionary multitasking

Feature selection (FS) plays a crucial role in high-dimensional classification problems by identifying relevant features that contribute to model performance. Evolutionary multitasking (EMT) has recently shown success in FS problems. However, existing EMT-based FS methods have limitations in terms of diversity in task construction, evolutionary search, and knowledge transfer, leading to inadequate acquisition, exploration, and utilization of knowledge. To this end, this paper develops a novel EMT framework for multi-objective high-dimensional FS problems, namely MO-FSEMT. In particular, multiple auxiliary tasks are constructed by distinct formulation methods to provide diverse search spaces and information representations and then simultaneously addressed with the original task by leveraging multiple evolutionary solvers with different biases and search preferences. A task-specific-based knowledge transfer mechanism is designed to leverage the advantageous information from each task, facilitating the discovery and effective transmission of high-quality solutions during the search process. Comprehensive experimental results on 27 real high-dimensional datasets demonstrate the superiority of MO-FSEMT over state-of-the-art FS methods in terms of effectiveness and efficiency. Ablation studies further confirm the contributions of key components of the proposed MO-FSEMT.

Application of clustering algorithms for dimensionality reduction in infrastructure resilience prediction models

Recent studies increasingly adopt simulation-based machine learning (ML) models to analyse critical infrastructure system resilience. For realistic applications, these ML models consider the component-level characteristics that influence the network response during emergencies. However, such an approach could result in a large number of features and cause ML models to suffer from the ’curse of dimensionality’. A clustering-based method is presented that simultaneously minimises the problem of high-dimensionality and improves the prediction accuracy of ML models developed for resilience analysis in large-scale interdependent infrastructure networks. The methodology has three parts: (a) generation of simulation dataset, (b) network component clustering, and (c) dimensionality reduction and development of prediction models. First, an interdependent infrastructure simulation model simulates the network-wide consequences of various disruptive events. The component-level features are extracted from the simulated data. Next, clustering algorithms are used to derive the cluster-level features by grouping component-level features based on their topological and functional characteristics. Finally, ML algorithms are used to develop models that predict the network-wide impacts of disruptive events using the cluster-level features. The applicability of the method is demonstrated using an interdependent power-water-transport testbed. The proposed method can be used to develop decision-support tools for post-disaster recovery of infrastructure networks.

A Novel Finite Element-Based Method for Predicting the Permeability of Heterogeneous and Anisotropic Porous Microstructures.

Permeability is a fundamental property of porous media. It quantifies the ease with which a fluid can flow under the effect of a pressure gradient in a network of connected pores. Porous materials can be natural, such as soil and rocks, or synthetic, such as a densified network of fibres or open-cell foams. The measurement of permeability is difficult and time-consuming in heterogeneous and anisotropic porous media; thus, a numerical approach based on the calculation of the tensor components on a 3D image of the material can be very advantageous. For this type of microstructure, it is important to perform calculations on large samples using boundary conditions that do not suppress the transverse flows that occur when flow is forced out of the principal directions. Since these are not necessarily known in complex media, the permeability determination method must not introduce bias by generating non-physical flows. A new finite element-based method proposed in this study allows us to solve very high-dimensional flow problems while limiting the biases associated with boundary conditions and the small size of the numerical samples addressed. This method includes a new boundary condition, full permeability tensor identification based on the multiscale homogenization approach, and an optimized solver to handle flow problems with a large number of degrees of freedom. The method is first validated against academic test cases and against the results of a recent permeability benchmark exercise. The results underline the suitability of the proposed approach for heterogeneous and anisotropic microstructures.

A Sampling Strategy for High-Dimensional, Simulation-Based Transportation Optimization Problems

When tackling high-dimensional, continuous simulation-based optimization (SO) problems, it is important to balance exploration and exploitation. Most past SO research focuses on the enhancement of exploitation techniques. The exploration technique of an SO algorithm is often defined as a general-purpose sampling distribution, such as the uniform distribution, which is inefficient at searching high-dimensional spaces. This work is motivated by the formulation of exploration techniques that are suitable for large-scale transportation network problems and high-dimensional optimization problems. We formulate a sampling mechanism that combines inverse cumulative distribution function sampling with problem-specific structural information of the underlying transportation problem. The proposed sampling distribution assigns greater sampling probability to points with better expected performance as defined by an analytical network model. Validation experiments on a toy network illustrate that the proposed sampling distribution has important commonalities with the underlying and typically unknown true sampling distribution of the simulator. We study a high-dimensional traffic signal control case study of Midtown Manhattan in New York City. The results show that the use of the proposed sampling mechanism as part of an SO framework can help to efficiently identify solutions with good performance. Using the analytical information for exploration, regardless of whether it is used for exploitation, outperforms benchmarks that do not use it, including standard Bayesian optimization. Using the analytical information for exploration only yields solutions with similar performance than when the information is used for exploitation only, reducing the total compute times by 65%. This paper sheds light on the importance of developing suitable exploration techniques to enhance both the scalability and the compute efficiency of general-purpose SO algorithms. Funding: T. Tay thanks the Agency for Science, Technology and Research (A*STAR) Singapore for funding his work. Supplemental Material: The online appendix is available at https://doi.org/10.1287/trsc.2023.0110 .