Unexpected Improvements to Expected Improvement for Bayesian Optimization

Motivation

• Bayesian optimization combines a probabilistic surrogate with an acquisition function to select queries for expensive black-box objectives and is used in fields including aerospace, biology, materials science, civil engineering, and machine learning [26].Its applications include aerospace engineering, biology and medicine, materials science, civil engineering [4], and machine-learning hyperparameter optimization [68].Source: §1 · p. 1

Contribution

• LogEI reformulates expected improvement so its optima are identical or approximately equal to canonical EI optima while being substantially easier to optimize numerically.The analytic LogEI variant yields the same Bayesian-optimization policy as EI but empirically improves optimization performance.Source: Contributions · p. 2

Extensions

• The LogEI approach extends to CEI, EHVI, qEI, and qEHVI by using smooth utility approximations that provide non-vanishing gradients.The resulting methods are available in BoTorch [6].Source: Contributions · p. 2

Results

• The proposed acquisition functions substantially outperform their canonical analogues across broad benchmarks without meaningful additional computational cost, often matching or exceeding recent methods.The reformulations target numerical pathologies that prior work addressed mainly through optimization algorithms and initialization strategies rather than EI computation itself.Source: Contributions · p. 2, §1 · p. 2

Motivation

• Acquisition optimization is generally non-convex with many local maxima, so gradient-based methods are typically preferred on continuous domains, especially in higher dimensions.Multiple starting points can reduce the risk of local-optimum traps, but improvement-based acquisition optimization becomes harder as more data reduces the likelihood of improving over the incumbent.Source: Motivation · pp. 2–3
• EI values and gradients frequently become numerically zero across large parts of the domain, making gradient-based acquisition optimization exceptionally difficult as observations accumulate.Although EI is not mathematically zero under a Gaussian posterior except at points perfectly correlated with previous noiseless observations, floating-point computation can produce exact zeros.Source: Motivation · pp. 2–3

Bayesian Optimization

• Bayesian optimization targets expensive noisy black-box functions, using a surrogate that predicts output distributions for candidate batches and an acquisition function to guide selection.The objective may have multiple outputs for multiple objectives or black-box constraints.Source: §2 · p. 3

Surrogate Models

• Gaussian Processes are widely used Bayesian-optimization surrogates because they provide data-efficient predictions with uncertainty quantification and multivariate Normal output distributions.In the sequential single-outcome case, the posterior is univariate Normal with mean µ(x) and variance σ^2(x).Source: §2.1 · p. 3

Acquisition Optimization

• Acquisition-function optimization is non-convex, and parallel acquisition functions are harder still because batch optimization increases dimensionality from d to qd and complicates the optimization surface.Common approaches include gradient-free, evolutionary, first-order, quasi-second-order, and multi-start methods; greedy sequential selection can simplify parallel optimization.Source: §2.3 · p. 4

Related Work

• Prior work introduced many acquisition functions and initialization strategies, but paid comparatively little attention to effectively implementing and optimizing them or to numerical challenges specific to improvement-based acquisitions.The covariance-matrix issues discussed by McLeod et al. [59] are described as orthogonal to improvement-based numerical pathologies.Source: §2.4 · p. 5

Conditions for Vanishing Gradients

• As the optimality gap f*−y* and posterior uncertainty decrease, EI is exceedingly likely to exhibit numerically vanishing gradients.This identifies the objective and surrogate conditions underlying the pathology shown in Figure 1.Source: §3 · p. 5

Theorem 1

• Theorem 1 shows that, with high probability, the probability of inputs producing a sufficiently small EI argument is bounded by the objective’s probability of falling below f*−ϵ_n.Here, ϵ_n depends on the optimality gap f*−y* and the maximum posterior uncertainty max_x σ_n(x).Source: §3 · p. 5

Non-asymptotic Behavior

• The bound need not be tight in early iterations, when both the optimality gap and maximum posterior standard deviation can be large.Thus, the theorem provides non-asymptotic insight rather than implying that vanishing gradients occur only late in optimization.Source: §3 · p. 5

Asymptotic Behavior

• Under technical kernel and training-data conditions, the maximum posterior variance is guaranteed to vanish as the number of observations increases [52].This asymptotic result helps explain why decreasing surrogate uncertainty is relevant to EI’s numerical pathology.Source: §3 · p. 5

Analytic EI

• LogEI computes the analytically equivalent log-transformed EI stably across inputs where naïve implementations become numerically zero, including highly negative standardized improvement values.The reformulation uses a stable implementation of log_h and adds log(σ(x)); an asymptotic case remedies instabilities as z grows negative.Source: §4.1 · p. 6

Monte Carlo parallel EI

• qLogEI smooths qEI’s per-sample discrete maximum and applies stable log-space computation, avoiding the mathematically zero gradients caused by differentiating Monte Carlo samples.The smoothing temperature must be sufficiently low for good optimization characteristics, which becomes feasible through log-space computation.Source: §4.2 · p. 6
• The relative approximation error between qLogEI and qEI is tightly bounded as a function of the temperature parameters τ0, τmax and batch size q.This guarantee quantifies how closely the smoothed acquisition follows canonical qEI.Source: §4.2 · p. 6

Constraints

• LogEI extends to constrained optimization through stable log-probability terms, while Monte Carlo constraints use smooth indicator approximations and accurate log implementations.For independent constraints fi(x) ≤0, analytic LogCEI combines LogEI with log(P(fi(x) ≤0)).Source: §4.3 · p. 7

Multi-objective optimization

• The same smoothing and log-transformation principles address qEHVI’s numerical difficulties, using a differentiable inclusion-exclusion formulation of hypervolume improvement.qEHVI is more complex because many operations have mathematically zero gradients with respect to some inputs.Source: §4.4 · p. 7

Single-objective sequential BO

• LogEI continues making progress on Sum-of-Squares after EI stalls at about 75 observations because EI develops vanishing gradients.The acquisition value for EI becomes jagged as optimization stalls, whereas LogEI maintains steady progress.Source: §5 · p. 7, Figure 2 · p. 7
• LogEI substantially outperforms EI on Ackley, with the performance gap widening as dimensionality increases, and shows a smaller advantage on Michalewicz.Vanishing gradients make the exploration–exploitation tradeoff especially difficult on multimodal Ackley.Source: §5 · p. 7, Figure 3 · p. 8

BO with Black Box Constraints

• LogCEI outperforms naive CEI and converges faster than SCBO on four constrained engineering-design problems, with gains over CEI increasing with dimensionality [22].All constraints use the same bilog transform as SCBO to facilitate Gaussian-process modeling.Source: Figure 4 · p. 8, §5 · p. 8

Parallel Expected Improvement with qLogEI

• qLogEI improves parallel optimization over qEI on 16-dimensional Ackley, while jointly optimizing batches outperforms sequential greedy optimization.The result holds across the evaluated batch sizes and uses fat-tailed nonlinearities from Appendix A.4.Source: §5 · p. 8, Figure 5 · p. 9

High-dimensional BO with qLogEI

• qLogEI significantly improves the embedded Hartmann problem, but differences on the SVM and Rover problems are not significant, with results varying by GP or SAASBO [21].On Hartmann, canonical-GP BO eventually catches up with the SAAS-prior-equipped model when qLogEI is included.Source: Figure 6 · p. 9, §5 · p. 9

Multi-Objective optimization with qLogEHVI

• qLogEHVI consistently outperforms qEHVI and JES across all tested batch sizes on two multi-objective problems, including cell network design.The comparison covers DTLZ2 and the real-world-inspired cell network design problem.Source: §5 · p. 9

• EI’s gradients vanish when high objective values are highly concentrated or as optimization progresses, but these conditions do not hold for every Bayesian optimization application.Consequently, LogEI’s performance depends on whether the problem exhibits these numerical pathologies and on surrogate quality.Source: §6 · p. 9

On problem dimensionality

• LogEI’s advantages generally grow with dimensionality because high objective values become concentrated in the search space, not because of dimensionality itself.For high ambient but low intrinsic dimensionality, such as the SVM problem in Figure 6, LogEI may not significantly improve over EI.Source: §6 · p. 9

On asymptotic improvements

• LogEI’s performance improvements may be small or emerge only in later iterations, despite generally enabling better optimization and higher acquisition values.Examples include the convex 10D sum of squares and q = 16 on DTLZ2 in Figure 7.Source: §6 · p. 10

On model quality

• Even with concentrated good values and many iterations, poor surrogate predictions or uninformative uncertainties can prevent LogEI from outperforming EI.In such cases, better acquisition values do not necessarily produce better Bayesian optimization performance, as illustrated by Rover in Figure 6.Source: §6 · p. 10

• The authors strongly recommend replacing EI variants with LogEI counterparts, because EI dominance would indicate that the EI family is itself sub-optimal.Any resulting improvement could otherwise be attributed to the exploratory quality of randomly distributed candidates, which could be incorporated explicitly.Source: §6 · p. 10

7 Conclusion

• Careful reformulations and implementations mitigate vanishing-gradient problems, substantially improving optimization across modified EI variants and making joint batch optimization competitive with, and sometimes better than, alternatives.The conclusion reports these gains across a broad range of problems, while noting that the methods may also inform other numerically challenging acquisition functions [78], [53] [68].Source: §7 · p. 10

A Acquisition Function Details

• For increasingly negative z, stable analytic LogEI evaluates Φ(z)/ϕ(z) with erfcx and the final logarithm with log1mexp, while expm1 and log1p provide canonical stable primitives.These reformulations expose the quadratic factor and avoid direct unstable evaluations.Source: §A.1 · p. 17

A.1 Analytic Expected Improvement

• Naïve analytic EI becomes numerically zero for moderately negative z, whereas the proposed log_h evaluation exhibits quadratic asymptotic behavior suited to numerical optimization.This failure mode and the improved behavior are shown in Figure 8.Source: §A.1 · p. 17, Figure 8 · p. 17
• When max_x EI(x) > 0, LogEI and EI have the same maximizers; when max_x EI(x) = 0, EI treats every point as a maximizer while LogEI is undefined everywhere.Thus, the two acquisition functions induce the same Bayesian-optimization policy in the zero-maximum case.Source: §A.1 · p. 17

A.2 Analytical LogEI’s Asymptotics

• A Laurent expansion supplies an approximate asymptotic computation for large negative z, addressing instabilities that can remain in existing analytic LogEI implementations.The appendix derives an asymptotic error bound after the logarithmic transformation.Source: §A.2 · p. 18

A.3 Monte-Carlo Expected Improvement

• Monte-Carlo LogEI first smooths max(0, x) with softplus and the maximum over q candidates with an ℓ_(1/τmax) norm, then transforms computations into log space.The approximation error from both smoothings is tightly bounded by temperature parameters τ0 and τmax.Source: §A.3 · p. 20
• The lower branch of logsoftplus has O(z^2) approximation error and converges exponentially quickly as x → −∞, with the threshold chosen to avoid losing significant digits.For large batches, however, logsoftplus outputs can still have a large dynamic range and produce vanishing gradients, motivating smoother O(1/x^2)-decaying nonlinearities.Source: §A.3 · p. 20

A.4 A Class of Smooth Approximations with Fat Tails for Larger Batches

• Fat-tailed approximations replace exponential decay with polynomial decay while preserving smoothness, positivity, monotonicity, and convexity, moderating dynamic range for larger batches.The proposed functions use Lorentzian/Cauchy-inspired terms with O(1/x^2) decay as x → −∞, unlike canonical softplus decay O(exp(x)).Source: §A.4 · pp. 21–22, Figure 12 · p. 21
• The fat softplus approximates max(x, 0) tightly while remaining monotonic, convex, and positive, and its convexity condition also guarantees monotonicity.For the proposed φ+, 0 ≤ α < 1/2 is sufficient for monotonicity, while the stronger convexity condition subsumes it.Source: Figure 11 · p. 21, §A.4 · p. 22

A.5 Constrained Expected Improvement

• For constrained EI, the fat sigmoid provides a positive, monotone constraint-indicator approximation that approaches 1 at positive infinity and decays as O(1/x^2) at negative infinity.A temperature parameter makes the approximation tighter to the Heaviside step function away from the origin; analytical benchmarks use stable log Gaussian-CDF implementations, while Monte Carlo parallel EI uses the fat-tailed indicator.Source: §A.4 · p. 23, §A.5 · p. 23

A.6 Parallel Expected Hypervolume Improvement

• For qEHVI, the stable Log formulation uses fat-tailed minimum and positive-part approximations, separates positive and negative inclusion-exclusion terms, and applies logsumexp reductions.The construction handles signed sums by computing a stable logarithm of the difference between exponentiated positive and negative components before reducing over hyperrectangles and Monte Carlo samples.Source: §A.6 · p. 24

A.7 Probability of Improvement

• Probability of improvement gains numerical stability by taking the logarithm of its Gaussian-CDF expression and evaluating the positive-tail case with erfcx.The supplied formulation specifically rewrites the logarithm using log(erfcx(x)) − x^2 for x > 0.Source: §A.7 · p. 24

A.8 q-Noisy Expected Improvement

• qLogNEI extends qLogEI’s numerical improvements to noisy EI by smoothing sample-level improvements, approximating the batch maximum, and operating in logarithmic space.Its incumbent is computed with a hard maximum over previously evaluated designs because gradients with respect to those designs are unnecessary; pruning and cached Cholesky updates provide speedups.Source: §A.8 · pp. 24–25

B Strategies for Optimizing Acquisition Functions

• For continuous acquisition optimization, gradient-based methods are standard, whereas grid or random search does not scale to high-dimensional domains and zeroth-order methods generally underperform gradient-based alternatives.Analytic acquisitions commonly use L-BFGS-B, while constrained or Monte Carlo settings use general nonlinear optimizers or stochastic gradient ascent.Source: B Strategies for Optimizing… · p. 25

B.1 Common initialization heuristics for multi-start gradient-descent

• Multi-start gradient descent mitigates local optima through varied initializations, but random restarts become ineffective when acquisition values are numerically zero across large regions, especially in higher dimensions and with more data.Initialization heuristics can help, but the foundational remedy is to ensure acquisition functions provide numerically meaningful gradients.Source: §B.1 · pp. 25–26

C Proofs

• Theorem 1 derives a probability-1−δ bound by decomposing the improvement argument into posterior-mean error and the gap between the true objective and the best observed value.The posterior standardized error is Gaussian, enabling a standard Gaussian-tail bound that yields a sufficient condition for controlling posterior mean minus incumbent.Source: C Proofs · p. 26
• The resulting bound depends on the optimality gap f*−y* and maximal posterior standard deviation, not posterior statistics evaluated at any specific x.The proof introduces ϵ_n from the optimality gap and maximal posterior uncertainty after rearranging the bound.Source: C Proofs · p. 27
• Lemma 2 bounds the approximation error of qLogEI to qEI for τ_0, τ_max > 0.The proof treats the Monte Carlo approximation first and then obtains the mathematical acquisition-value result as the sample size tends to infinity.Source: C Proofs · p. 27
• The qLogEI-to-qEI proof bounds the sample-wise difference using logsumexp, logsoftplus, ReLU, and the triangle inequality.The argument compares the relevant nonlinear terms and then divides by the number of samples to obtain the Monte Carlo acquisition-value bound.Source: C Proofs · p. 27
• The softplus approximation satisfies ||softplus_τ−ReLU||_∞ = log(2)τ, with the maximum discrepancy at x=0 and softplus_τ(x) above ReLU(x) for all x.The difference approaches zero as x tends to either infinity, while its extremum occurs at the origin.Source: C Proofs · pp. 27–28
• The same style of approximation guarantee can be derived for the fat-tailed nonlinearities introduced in Appendix A.4.The passage states that these guarantees follow by an analogous argument.Source: C Proofs · p. 28

• All algorithms are implemented in BoTorch, including analytic EI, qEI, and cEI through standard BoTorch implementations.Source: §D.1 · p. 28
• Single-objective JES [38], GIBBON [63], and multi-objective JES use the original authors’ implementations available in BoTorch’s main repository.Source: §D.1 · p. 28
• Each simulation uses 32 replicates, with error bars representing ±2 times the standard error of the mean.Source: §D.1 · p. 28
• The experiments use a Matérn-5/2 kernel with automatic relevance determination (ARD).Source: §D.1 · p. 28

D.2 Additional Empirical Results on Vanishing Values and Gradients

• The vanishing-gradient simulations use training data with 80% uniformly sampled points and 20% Gaussian points centered at the function maximum, mimicking a Bayesian optimization loop.A representative 60-point fit has decent out-of-sample mean predictions and reasonably calibrated uncertainty, despite the numerical pathology.Source: §D.2 · p. 28, Figure 15 · p. 29
• Figure 16 shows that most test points fall below thresholds where h(z) is small, explaining why EI values and gradients numerically vanish across much of the domain.Because σ(x) is close to 1 for most test points, the threshold behavior of h(z) approximately corresponds to EI threshold behavior.Source: §D.2 · p. 28, Figure 16 · p. 29

D.3 Parallel Expected Improvement

• On 16-dimensional Ackley, qLogEI outperforms all baselines, while joint batch optimization improves over sequential greedy; on Levy, joint qLogEI eventually wins as batch size increases.On Levy, joint qLogEI initially underperforms sequential optimization but becomes superior as parallel batch size grows.Source: Figure 17 · p. 29, §D.3 · p. 29

D.4 Noisy Expected Improvement

• qLogNEI’s advantage over canonical noisy EI grows with dimensionality and noise level, significantly outperforming both canonical counterparts and Gibbon in most cases.The benchmark covers noisy Hartmann 6d, Ackley 8d, and Ackley 16d problems across varying noise levels.Source: §D.4 · p. 29, Figure 19 · p. 30

D.5 Multi-Objective optimization with qLogEHVI

• qLogEHVI significantly outperforms qEHVI and all baselines across six multi-objective problems, with larger gains at higher dimensionality.The problems span 2–3 objectives and 2–30 dimensions, including real-world-inspired applications in cell network, laser plasma, and vehicle design [40] [70].Source: §D.5 · p. 30, Figure 20 · p. 31

D.6 Combining LogEI with TuRBO for High-Dimensional Bayesian Optimization

• TuRBO-1 combined with LogEI substantially outperforms TuRBO-1 + EI, LogEI, and Thompson sampling at small evaluation budgets on 50d Ackley.Thompson sampling performs better after 10,000 evaluations, while joint q = 50 optimization with qLogEI needs no random restarts for good performance.Source: §D.6 · p. 31, Figure 21 · p. 31

Constrained Problems

• For pressure-vessel design, LogEI reached feasible objectives of 5659.1108 worst-case and 5651.8862 best-case after 110 evaluations, versus 6059.946341 from Coello and Montes.The LogEI results came from 16 runs and used almost three orders of magnitude fewer evaluations.Source: §D.7 · p. 32
• For welded-beam design, LogEI found 1.7496 after 110 evaluations versus 1.728226 from Coello and Montes, slightly worse but using three orders of magnitude fewer evaluations.Source: §D.7 · p. 32
• For tension-compression design, LogEI found a feasible value of 0.0129 after 110 evaluations versus 0.012681 reported in.Source: §D.7 · p. 32
• Across these constrained problems, Bayesian optimization was competitive with and could outperform genetic algorithms while using only a tiny fraction of the sample budget.Bayesian optimization emphasizes sample efficiency and small-data regimes, whereas genetic algorithms generally use substantially more function evaluations.Source: §D.7 · p. 32

D.8 Parallel Bayesian Optimization with cross-batch constraints

• The constrained experiment models a five-dimensional silver nanoparticle synthesis problem with a total-flow throughput limit that couples batch elements.The study uses a continuous RBF surrogate, batch size q = 3, five random initial points, and 75 replicated runs.Source: §D.8 · pp. 32–33
• qLogEI with proper joint batch constraints outperforms qEI, whereas the advantage disappears under the random-flow heuristic.Both constrained approaches outperform purely random sampling, but additional gains require combining LogEI with explicit constraints and jointly optimizing the batch.Source: Figure 22 · p. 33, §D.8 · p. 33

D.9 Details on Multi-Objective Problems

• Performance is evaluated on synthetic biobjective Branin-Currin, ZDT1, and DTLZ2 benchmarks alongside real-world-inspired problems.The laser plasma acceleration problem uses independent GP surrogates fitted to public data at the highest queried fidelity.Source: §D.9 · p. 33

D.10 Effect of Temperature Parameter

• Smaller fixed softplus temperatures typically produce better optimization performance.The temperature study examines the effect of fixed τ values on optimization performance.Source: §D.10 · p. 33

D.11 Effect of the initialization strategy

• The initialization study compares uniformly random restart points with BoTorch’s acquisition-weighted Boltzmann strategy.The comparison is motivated by common smart-initialization heuristics for acquisition-function optimization.Source: §D.11 · p. 33, Figure 24 · p. 35
• BoTorch’s default Boltzmann initialization improves regret across all evaluated cases compared with uniformly random initialization.The default evaluates 1024 scrambled-Sobol candidates and samples restart points using probabilities derived from acquisition values.Source: §D.11 · p. 34