Profile Bayesian Optimization for Expensive Computer Experiments

Technical Summary: Profile Bayesian Optimization for Expensive Computer Experiments

Problem Statement
The paper addresses the challenge of optimizing expensive, deterministic black-box computer simulations where the objective is to identify "profile optima" rather than a single global optimum. In these scenarios, inputs are divided into a single control parameter ($x^\star$) and multiple nuisance parameters ($x_{-\star}$). The goal is to estimate the function $T(x^\star) = \min_{x_{-\star}} f(x^\star, x_{-\star})$ across the entire plausible range of $x^\star$.

Standard Bayesian Optimization (BO) is ill-suited for this task because it typically targets a single global minimum, often neglecting the exploration of the full control parameter space. Existing profile optimization methods, such as Profile Expected Improvement (PEI) proposed by Ginsbourger et al. (2014), suffer from two main limitations:

1. Computational Intractability: They rely on dense grids for candidate selection, which becomes infeasible in dimensions higher than two.
2. Exploration-Exploitation Imbalance: Like traditional Expected Improvement (EI), PEI tends to over-exploit regions near the global minimum, failing to adequately explore the entire range of the control parameter.

The authors are motivated by a specific engineering application: designing a diffuser for a Rotating Detonation Combustor (RDC). The simulation is computationally expensive (10 minutes to several hours per run) and nonstationary due to combustion dynamics. The objective is to map the minimum energy loss (pressure loss) as a function of diffuser length ($x^\star$), optimizing over five nuisance parameters (curvature points and bypass location) to inform design trade-offs between efficiency and physical constraints.

Methodology
The authors propose a novel Profile Bayesian Optimization (PBO) procedure that integrates surrogate modeling with a two-stage acquisition strategy.

1. Surrogate Modeling:

   • The method utilizes either Gaussian Process (GP) or Deep Gaussian Process (DGP) surrogates. DGPs are highlighted for their ability to handle nonstationary response surfaces common in physical simulations.
   • To estimate $T(x^\star)$, the authors generate joint posterior samples from the surrogate. For each sample, they identify the minimum response value for every unique $x^\star$. The empirical distribution of these minimums provides an estimate $\hat{T}(x^\star)$ along with uncertainty quantification (UQ) via credible intervals.
   • To manage computational costs in high dimensions, Vecchia approximation is employed for efficient joint posterior sampling.

2. Modified Tricands (Candidate Generation):

   • To avoid the computational bottleneck of dense grids, the authors adapt the "tricands" (triangulation candidates) approach by Gramacy et al. (2022).
   • Instead of a full grid, candidates are generated using a Delaunay triangulation of the observed nuisance parameters ($X_{-\star}$). "Internal candidates" are placed at the geometric centers of triangulation components, and "fringe candidates" are placed near the boundaries.
   • A semi-dense 1D grid is created for the control parameter $x^\star$. The final candidate set ($X_{tri+}$) is the Cartesian product of the $x^\star$ grid and the nuisance parameter tricands. This allows for efficient minimum estimation without evaluating billions of points.

3. Two-Stage Acquisition Scheme:
   The selection of the next evaluation point $x_{n+1} = \{x^\star_{n+1}, x_{-\star, n+1}\}$ occurs in two distinct steps:

   • Stage 1 (Exploration of $x^\star$): The algorithm selects the control parameter value $x^\star_{n+1}$ that exhibits the largest uncertainty (widest credible interval width) in the current estimate of $\hat{T}(x^\star)$. This is inspired by the Active Learning MacKay criterion and ensures the entire range of the control parameter is explored.
   • Stage 2 (Exploitation of $x_{-\star}$): Once $x^\star_{n+1}$ is fixed, the algorithm selects the nuisance parameters $x_{-\star, n+1}$ by maximizing the Profile Expected Improvement (PEI) along that specific slice. This focuses the search on finding the optimal nuisance configuration for the chosen control parameter.

Key Results
The paper validates the PBO procedure through synthetic benchmarks and the motivating RDC application:

• Synthetic Benchmarks:

  • Stationary Functions (2D Branin, 3D Kyger): PBO outperformed standard BO, random Latin Hypercube Sampling (LHS), and the original PEI method. PBO achieved lower Root Mean Squared Error (RMSE) and narrower credible intervals while maintaining coverage near 95%. Crucially, PBO successfully explored the full range of $x^\star$, whereas PEI and standard BO clustered acquisitions near global minima, leaving parts of the profile unexplored.
  • Nonstationary Function (4D Squiggle): Using DGP surrogates, PBO demonstrated superior accuracy compared to GP surrogates. While DGP predictions showed slightly lower coverage (indicating overconfidence), they yielded significantly lower RMSE and Maximum Absolute Deviation (MaxAD) than GP-based methods.

• Motivating Application (RDC Diffuser):

  • Applied to a 6D simulation (1 control, 5 nuisance parameters), PBO with DGP surrogates was run for 150 evaluations.
  • The resulting profile estimate provided uncertainty quantification across the entire range of diffuser lengths.
  • When compared to a state-of-the-art multi-objective evolutionary optimizer (CADO) run on the same budget (113 evaluations), PBO identified designs with lower pressure losses than the Pareto front found by CADO. Furthermore, PBO provided a continuous profile of optima, whereas CADO only provided a discrete set of trade-off points.

Significance and Claims
The paper claims that Profile Bayesian Optimization offers a robust framework for estimating profile optima in expensive computer experiments where the control parameter range is of primary interest.

• Novelty: The primary contribution is the two-stage acquisition scheme that explicitly separates the exploration of the control parameter from the exploitation of nuisance parameters. This prevents the "over-exploitation" common in traditional BO and PEI.
• Scalability: By combining modified tricands with Vecchia approximation, the method extends profile optimization to higher dimensions (demonstrated up to 6D) where previous grid-based methods failed.
• Flexibility: The procedure is compatible with both stationary (GP) and nonstationary (DGP) surrogates, making it applicable to complex physical simulations with regime shifts.
• Practical Impact: In the RDC application, the method enabled informed decision-making by quantifying the trade-off between diffuser length and energy loss, identifying designs that outperformed existing multi-objective optimization approaches.

The authors modestly note limitations, including the assumption of unique profile optima for each $x^\star$ slice and the computational cost of generating sufficient posterior samples for DGP coverage in high dimensions. They suggest future work could involve gradient-enhanced surrogates and Bayesian calibration with physical data.