High-dimensional Bayesian optimization (BO) often fits Gaussian process (GP) surrogates from far fewer observations than input dimensions. Modern Vanilla BO can perform well in this regime with dimension-aware priors, initialization, and acquisition optimization, but it typically retains automatic relevance determination (ARD), fitting one lengthscale per input coordinate. We study this modeling choice and propose Iso-BO, a controlled modification that replaces the ARD GP with an isotropic GP using one shared lengthscale while keeping the surrounding BO pipeline matched. For radial kernels, we show that the marginal log likelihood (MLL) depends on the inverse-squared ARD lengthscales only through weighted pairwise distances among the observed inputs. The current design can therefore leave some ARD directions exactly invisible or only weakly constrained by the MLL. Iso-BO removes coordinatewise reweighting and fits a single shared scale instead. Lengthscale-fitting and predictive-density diagnostics show that this finite-data effect appears in practice, including when the data-generating process is anisotropic. Across GP-prior, synthetic, and real-world benchmarks, Iso-BO often improves over matched modern Vanilla BO and remains competitive with the included high-dimensional BO baselines under the tested budgets. Stress tests also show the expected boundary wherein sufficiently strong, learnable anisotropy can favor the more flexible ARD model.
Figures & tables
Figure 1: Controlled lengthscale fitting on isotropic SE GP-prior realizations with ground-truth lengthscale ℓ=1 . Isotropic and ARD GPs are fit by MAP using L-BFGS-B or GD, with Sobol training inputs. Curves show the mean and 95% confidence interval over 10 independent realizations; the dotted line is the ground truth and the green line is the prior mode. For ARD fits, the plotted value is the arithmetic mean across coordinatewise lengthscales and is used only as a scalar summary of the fitted ARD vector.
Figure 2: Held-out NLPD for isotropic and ARD SE GP fits on functions sampled from anisotropic SE GP priors with γ=0.3 . Both models are fit by MAP with L-BFGS-B using the same Sobol training inputs; held-out responses are standardized using training-set statistics. Curves show the mean and 95% confidence interval over 20 independent GP-prior realizations. Lower NLPD is better, and the isotropic fit has lower NLPD throughout the reported finite-data regimes.
Figure 3: Best-observed values on 100D anisotropic SE GP-prior objectives with γ∈{0.3,0.5,0.7} . Curves show the median and IQR over 10 runs; lower is better. The generating GP has coordinatewise lengthscales, so these objectives are deliberately misspecified for Iso-MAP. Despite this, Iso-MAP reaches lower median values than Vanilla-MAP in all settings.
Figure 4: Best-observed values on 100D Schwefel, Styblinski–Tang, and Powell. Curves show the median and IQR over 10 runs; lower is better. Iso-MAP reaches lower final median values than matched Vanilla-MAP on all three benchmarks. TuRBO attains the lowest median values on Powell.
Figure 5: Best-observed values on six high-dimensional non-synthetic benchmarks. Curves show the median and IQR over 10 runs; lower is better. Iso-MAP shows its clearest gains over matched Vanilla-MAP on Rover, HalfCheetah-v2, and HumanoidStandup-v5, while remaining competitive across the broader benchmark suite. SAASBO uses 300 evaluations on HumanoidStandup-v5 and is not run on Humanoid-v2. BAxUS is omitted from HalfCheetah-v2 because its standard subspace initialization produce substantially lower initial function value compared with other baseline.
Appendix figures & tables33 assets
Supplementary material from the paper’s appendix.
Appendix
Figure 6: Training MAP objective values for SE-kernel ARD GP fits using L-BFGS-B and GD. GD uses learning rates 0.01 , 0.1 , and 1 . The plotted quantity is the average unnormalized log posterior per training point, i.e., MLL plus the log-prior contribution; higher is better. Top row: isotropic GP-prior data with ground-truth lengthscale ℓ=1 ; middle row: Schwefel; bottom row: Styblinski–Tang. Training inputs are Sobol points. Curves show the mean and 95% confidence interval over 10 independent replicates. L-BFGS-B generally attains higher training MAP objective values than the tested GD fits. GD is run for at most 1000 steps.
Figure 7: Learning-rate ablation for SE-kernel ARD GP fits using GD. The figure reports the average fitted ARD lengthscale across coordinates for learning rates 0.01 , 0.1 , and 1 ; the black line shows the prior mode. Top row: isotropic GP-prior data with ground-truth ℓ=1 ; middle row: Schwefel; bottom row: Styblinski–Tang. Training inputs are Sobol points. Curves show the mean and 95% confidence interval over 10 independent replicates. GD is run for at most 1000 steps.
Figure 8: Held-out NLPD for isotropic and ARD GP fits along Vanilla-MAP trajectories on anisotropic SE GP-prior objectives with γ=0.3 . At each reported training size, both models are fit to the first n observations and evaluated on later observations from the same completed trajectory. MAP estimation uses L-BFGS-B. Curves show the mean and 95% confidence interval over 20 independent BO trajectories; lower NLPD is better.
Figure 9: Held-out NLPD for isotropic and ARD GP fits on Schwefel and Styblinski–Tang in 30D, 50D, and 100D. MAP estimation uses L-BFGS-B and training inputs are Sobol points. Curves show the median and 95% confidence interval over 10 independent replicates; lower NLPD is better. The isotropic fits have lower NLPD throughout the reported low-data regimes.
Figure 10: Controlled lengthscale fitting on isotropic GP-prior data using Matérn 5/2 fitted kernels. The data-generating shared lengthscale is ℓ=1 . MAP estimation uses L-BFGS-B or GD. Dashed curves denote isotropic GPs and solid curves denote ARD GPs; the black dotted line shows the ground truth and the green line shows the prior mode. Training inputs are Sobol points. Results show the mean and 95% confidence interval over 10 independent GP-prior realizations. For ARD fits, the plotted value is the average fitted lengthscale across coordinates.
Figure 11: Lengthscale fitting on anisotropic GP-prior data using SE fitted kernels. The data-generating coordinatewise lengthscales are sampled from Eq. 3 with γ=0.5 . MAP estimation uses L-BFGS-B or GD. Dashed curves denote isotropic GPs and solid curves denote ARD GPs. The black dotted line shows the replicate-averaged mean of the generating coordinatewise lengthscales and the green line shows the prior mode. Results show the mean and 95% confidence interval over 10 independent GP-prior realizations. For ARD fits, the plotted value is the average fitted lengthscale across coordinates.
Figure 12: Lengthscale fitting on the same anisotropic GP-prior setting as Figure 11 , now using Matérn 5/2 fitted kernels. The data-generating coordinatewise lengthscales are sampled from Eq. 3 with γ=0.5 . The black dotted line shows their replicate-averaged mean and the green line shows the prior mode. Results show the mean and 95% confidence interval over 10 independent GP-prior realizations. For ARD fits, the plotted value is the average fitted lengthscale across coordinates. GD is run for at most 1000 steps.
Figure 13: Fitted lengthscales on Styblinski–Tang ( top row ) and Schwefel ( bottom row ) using SE fitted kernels. MAP estimation uses L-BFGS-B or GD with Sobol training inputs. The green line shows the prior mode. Results show the mean and 95% confidence interval over 10 independent replicates; for ARD fits, the plotted value is the average fitted lengthscale across coordinates. These deterministic objectives do not have a GP-generating lengthscale. GD is run for at most 1000 steps.
Method
Lengthscale model
Hyperparameter fit
Main kernel
Vanilla-MAP
ARD ( D parameters)
MAP
SE
Iso-MAP
Shared (1 parameter)
MAP
SE
MSR
ARD ( D parameters)
MLE
Matérn 5/2
Iso-MLE
Shared (1 parameter)
MLE
SE
Appendix
Table 1: GP variants used in the main BO experiments. Vanilla-MAP and Iso-MAP form the primary matched comparison: they use the same kernel family and fitting procedure and differ in the lengthscale parameterization. The current MSR and Iso-MLE main runs use different kernel families and are therefore treated as a secondary comparison.
Figure 14: Final best-observed value distributions for the 100D Schwefel, Styblinski–Tang, and Powell experiments. Each violin summarizes the final best-observed value over 10 independent runs. Lower values are better. This figure complements the optimization trajectories in Figure 4 by showing the run-to-run variation at the final budget.
Figure 15: Final best-observed value distributions for the 100D anisotropic GP-prior experiments with γ∈{0.3,0.5,0.7} . Each violin summarizes the final best-observed value over 10 independent runs. Lower values are better. The objectives and baseline implementations are the same as in Figure 3 .
Figure 16: Final best-observed value distributions for the real-world benchmarks. Each violin summarizes the final value over 10 independent runs. The final budget is 1000 evaluations for the main methods; SAASBO uses 500 evaluations on Rover, HalfCheetah-v2, Walker2d-v5, and SVM and 300 evaluations on HumanoidStandup-v5, and is not run on Humanoid-v2. Lower values are better; reward-maximization tasks are reported as negative reward.
Problem
Median Δ
95% paired CI
Iso-MAP wins
p
Synthetic objectives
100D Schwefel
-9451.71
[-14169.33, -4001.42]
9/10
0.002
100D Styblinski–Tang
-1150.74
[-1339.72, -1049.55]
10/10
0.001
100D Powell
-3898.00
[-5283.70, -2611.16]
9/10
0.002
Anisotropic GP-prior objectives
100D, γ=0.3
-0.88
[-1.33, -0.09]
8/10
0.010
Appendix
Table 2: Paired final-value comparison between Iso-MAP and Vanilla-MAP. Here Δ=fIso∗−fVanilla∗ , so negative values favor Iso-MAP. The 95% interval is the bootstrap percentile confidence interval for the median paired difference. “Iso wins” is the number of paired replicates for which Iso-MAP has a lower final best-observed value. The final column gives the one-sided Wilcoxon signed-rank p -value for a negative paired shift.
Figure 17: Best-observed values on highly anisotropic GP-prior objectives in 50D and 100D with γ∈{0.9,0.99} . The budgets are 500 evaluations in 50D and 1000 in 100D. Iso-MAP remains close to matched Vanilla-MAP at γ=0.9 , while Vanilla-MAP reaches lower final median values at γ=0.99 .
Figure 18: Best-observed values on 100D anisotropic GP-prior objectives with 50 inactive dimensions. The budget is 500 evaluations. Iso-MAP remains competitive with matched Vanilla-MAP for γ=0.3 and 0.7 and reaches a lower median trace for γ=0.5 .
Figure 19: Best-observed values on anisotropic GP-prior objectives in 50D with γ∈{0.3,0.5,0.7} . The budget is 500 evaluations. Iso-MAP reaches lower final median values than matched Vanilla-MAP for γ=0.3 and 0.5 and remains close for γ=0.7 .
Figure 20: Best-observed values on the 50D Schwefel, Styblinski–Tang, and Powell benchmarks with a budget of 500 evaluations. Iso-MAP finishes below matched Vanilla-MAP in all three reported settings.
Figure 21: Best-observed values on the 60D Rover and 102D Walker2d-v2 benchmarks. The budget is 1000 evaluations for the main methods and 500 for SAASBO. Iso-MAP remains close to Vanilla-MAP on 60D Rover and reaches substantially lower median values on Walker2d-v2.
Figure 22: Best-observed values on the 124D MOPTA08 vehicle-design benchmark with 68 constraints and N=3000 evaluations. The constrained benchmark is converted to an unconstrained objective using the soft penalty specified in Appendix C.2 . The comparison is restricted to Iso-MAP, Linear BO, and CMA-ES; lower values are better.
Figure 23: SE versus Matérn 5/2 kernels on the 50D and 100D Schwefel, Styblinski–Tang, and Powell benchmarks. Budgets are 500 evaluations in 50D and 1000 in 100D. The two kernel choices give similar Iso-MAP traces across these settings, while Vanilla-MAP is more sensitive to the kernel choice in several panels.
Figure 24: SE versus Matérn 5/2 kernels on anisotropic GP-prior objectives in 50D and 100D with γ∈{0.3,0.5,0.7} . Budgets are 500 evaluations in 50D and 1000 in 100D. Iso-MAP shows similar behavior under the two kernel families across the reported settings.
Figure 25: SE versus Matérn 5/2 kernels on the real-world benchmarks. The budget is 1000 evaluations for each task. Kernel choice affects performance differently across problems; the Iso-MAP traces are comparatively similar under the two kernel families in most of the reported settings.
Figure 26: Query-spread metrics for SE-kernel Vanilla-MAP and Iso-MAP on the 50D anisotropic GP-prior objectives with γ=0.3 . Left: OTSD; middle: normalized OTSD; right: observation entropy (OE). Curves show the mean and ±1 standard deviation over 10 runs. Larger values indicate more spatially dispersed query locations.
Figure 27: Query-spread metrics for SE-kernel Vanilla-MAP and Iso-MAP on the 100D anisotropic GP-prior objectives with γ=0.3 . Left: OTSD; middle: normalized OTSD; right: observation entropy (OE). Curves show the mean and ±1 standard deviation over 10 runs. Larger values indicate more spatially dispersed query locations.
Figure 28: Wall-clock MAP hyperparameter fitting time for isotropic and ARD GPs. Objective functions are sampled from isotropic SE GP-prior realizations, while both fitted models use Matérn 5/2 kernels and L-BFGS-B. Training inputs are Sobol points. Curves show the mean and 95% confidence interval over 20 independent replicates. Experiments use 4-core CPU jobs on nodes with Intel Xeon CPU Max 9470 processors and 512 GB DDR5 RAM. Both models retain the same leading O(n3) exact-GP inference cost.
Method
Evaluations
GPU time (s)
Vanilla-MAP
1000
3255.2
Iso-MAP
1000
1783.0
MSR
1000
7334.5
Iso-MLE
1000
1431.1
Linear BO
1000
1982.0
TuRBO
1000
3556.6
Appendix
Table 3: End-to-end GPU runtime on the 100D anisotropic GP-prior objective with γ=0.3 . GP hyperparameters are refit at every BO iteration for this timing experiment. All runs use an NVIDIA L40 GPU with 45 GB of memory. SAASBO uses 500 evaluations; all other methods use 1000.
Figure 29: Batch BO on anisotropic GP-prior objectives in 50D and 100D with γ∈{0.3,0.5,0.7} . Each BO iteration selects q=10 points with qLogEI. The total budgets are 500 evaluations in 50D and 1000 in 100D. Lower best-observed values are better. Iso-MAP reaches lower median values over most of the budget in all six settings, with the smallest final separation in the 50D γ=0.7 case.
Figure 30: Batch BO on the same anisotropic GP-prior objectives as Figure 29 , now selecting q=50 points per iteration. The total budgets are 500 evaluations in 50D and 1000 in 100D. Lower best-observed values are better. Iso-MAP reaches lower median values across all six reported settings.
Figure 31: Batch BO on real-world benchmarks with q=10 and a total budget of 1000 evaluations per problem. Lower best-observed values are better; reward-maximization tasks are reported as negative reward. Iso-MAP reaches lower final median values on Rover, HalfCheetah-v2, Walker2d-v5, and HumanoidStandup-v5, while Vanilla-MAP finishes lower on Walker2d-v2 and SVM.
Figure 32: Batch BO on the same real-world benchmarks as Figure 31 , now with q=50 . The total budget is 1000 evaluations per problem. Iso-MAP reaches lower final median values on Walker2d-v2, HalfCheetah-v2, Walker2d-v5, and SVM; the final median values are similar for the two methods on Rover and HumanoidStandup-v5.
Figure 33: Batch BO on four 256D GuacaMol latent-space molecular optimization tasks with q=10 and a total budget of 3000 evaluations. The plotted objective is the negative GuacaMol task score, so lower values are better. Iso-MAP reaches lower final median values on fexo and med2 , Vanilla-MAP finishes lower on rano , and the two methods reach similar final values on med1 .
Figure 34: Batch BO on the same 256D molecular optimization tasks as Figure 33 , now with q=50 . The total budget is 3000 evaluations. Iso-MAP reaches lower final median values on fexo , med1 , and med2 ; the two methods finish at similar values on rano .
Figure 35: BO-FIM on strongly anisotropic GP-prior objectives in 50D and 100D with γ∈{0.9,0.99} . “Eigen Product” uses the product of the FIM eigenvalues with threshold 0.1 , while “Eigen Minimum” uses the minimum eigenvalue with threshold 0.01 . Budgets are 500 evaluations in 50D and 1000 in 100D. Curves show the median over 10 independent runs and shaded regions show the interquartile range. The adaptive variants remain competitive at γ=0.9 and recover substantial performance relative to fixed Iso-MAP at γ=0.99 .
Gaussian Process (GP) kernels are central to Bayesian optimization (BO), yet designing effective kernels for high-dimensional problems still relies on extensive manual engineering. Existing automated approaches struggle in high dimensions for two bottlenecks: their kernel search space is limited to additions and multiplications of base kernels, and LLM-based approaches require conditioning on raw observations, which becomes infeasible due to context-length limits and the difficulty of extracting meaningful patterns. We introduce \textbf{Kernel Discovery}, a LLM-driven evolutionary framework for high-dimensional BO that searches a broader kernel space beyond predefined composition rules and does not require conditioning on observations. Motivated by the observation that directly prompting an LLM to generate kernel code yields syntactically varied but functionally identical kernels, we adopt a two-stage approach: an LLM first proposes novel mathematical forms, then a second LLM call converts each form into validated, executable code. We also propose a leave-one-out continuous ranked probability score (LOO-CRPS) as a selection criterion that penalizes overfitted kernels. On five high-dimensional BO benchmarks, our method achieves an average rank of \textbf{1.2 out of 17}, outperforming competitive baselines. We further analyze the discovered kernels to identify which kernels lead to improvements in high-dimensional BO.
Taeyoung Yun, Woocheol Shin, Inhyuck Song +2
1Korea Advanced Institute of Science and Technology (KAIST)
Trust Region Bayesian Optimization (TuRBO) is an effective strategy for alleviating the curse of dimensionality in high-dimensional black-box optimization. However, inappropriate lengthscale design can cause the local Gaussian process (GP) model within the trust region to degenerate, leading to suboptimal performance in high dimensions. In this work, we show that TuRBO's local GP may remain either excessively complex or overly simple as the dimension D and trust region side length L vary. To address this issue, we propose a straightforward variant, AdaScale-TuRBO, which scales the GP lengthscale with both the problem dimension and trust region size, thereby preserving kernel geometry and maintaining consistent prior complexity. Empirically, we show that AdaScale-TuRBO can robustly outperform standard TuRBO and other popular high-dimensional BO methods on synthetic benchmarks and real-world trajectory planning tasks.
Gaussian-process Bayesian optimization (GP-BO) excels at black-box optimization of costly functions, e.g., hyperparameter optimization (HPO) and multi-agent system (MAS) design. Convergence-rate guarantees exist for select methods, notably GP upper confidence bound (GP-UCB), but require a fixed kernel. Critically, the kernel encodes how input proximity affects objective value similarity. When raw coordinates poorly match this geometry - as with log-scaled hyperparameters or localized peaks - input warping can greatly improve sample efficiency, yet known GP-UCB proofs require a fixed kernel. We propose Finite-Library Input-Warped Bayesian Optimization (FLIWBO), which selects warps from a finite library of smooth input maps by any history-dependent rule. It adapts the input geometry to accelerate learning while retaining high-probability convergence guarantees under mild hypotheses, with an explicit (Nε) library-size cost. Controlled diagnostics show that finite-library warping repairs planted geometry mismatches and identify FLIWBO failure cases. Across four repeated benchmarks - warped synthetic objectives, a confidence-fence trap, and Fashion-MNIST HPO - FLIWBO-UCB beats raw-coordinate GP-UCB under misspecified geometry, escapes traps that defeat even oracle-warp expected improvement, and recovers much of the gain from manual log scaling, while leading the tested methods that admit a matching regret guarantee. A 20-dimensional MAS design study further shows feasibility under costly noisy evaluations. Code for experiments is available: https://github.com/edvin-ketabati/bogp-paper-experiments.