[Bug]: Poor performance of `compute_best_feasible_objective` with no feasible points

[TobyBoyne]

What happened?

I am performing an optimization with qLogNEI and an outcome constraint. When generating candidates with no feasible values in the training data, the acqf is zero almost everywhere, leading to candidates that ignore the constraint. This leads to degenerate behaviour, where the optimizer will keep proposing the same (infeasible) point throughout the BO loop.

I think that the acqf is zero everywhere because the compute_best_feasible_objective does not give a good lower bound for infeasible values. I have found some issues with the current compute_best_feasible_objective, each of which contributes some part to this problem. Sorry for the very long ticket!!

1 Normalization. The convex weights are normalized in the wrong dimension below - the weights are (32, N) and should be normalized along N so that the weights assigned to each point sum to 1. The current implementation essentially scales down all points to search near the origin. Solution is to just fix the dimension here.
https://github.com/pytorch/botorch/blob/1518b304f47f5cdbaf9c175e808c90b3a0a6b86d/botorch/acquisition/utils.py#L192

2 Convex search excludes boundaries. By using convex weights, we don't explore the edges of the space which is where points satisfying the threshold are likely to lie. I don't think there is a good solution to this, since we don't have access to bounds when the class is initialized.

3 Random convex weights favour middle. By using random convex weights, we tend to explore very close to the middle. I think this is due to the Central Limit Theorem, but I haven't thought about this too deeply. Try simulating:

w = torch.rand(32, N)
(w / w.sum(dim=1, keepdim=True)) @ torch.linspace(0, 1, N)

As N increases, the above samples tend closer and closer to 0.5, so we don't actually get to explore the space very much. A solution would be to pass torch.cat((X, convex_weights @ X)) to get_infeasible_cost, to at least ensure that we evaluate the minima at the points on the convex hull.
https://github.com/pytorch/botorch/blob/1518b304f47f5cdbaf9c175e808c90b3a0a6b86d/botorch/acquisition/utils.py#L198

4 Prune baseline makes this worse When we pass prune_baseline=True (the default), we only keep the infeasible point that maximizes the objective. However, we should also keep the point that is closest to being feasible, as we otherwise only evaluate the lower bound at an area of the search space that is very far away from being feasible. Solution is to either change this default, or be less aggressive in pruning the baseline.

5 Account for objective direction This is a bit of an aside, but also an issue here. The code below assumes that the objective is increasing in Y, and that evaluating objective(mean - 6*std) will give a lower bound on the objective. However if the objective is decreasing, then the lower bound should be given by objective(mean + 6*std). A simple solution is to compute both - and +, then just take the minimum of each of those.
https://github.com/pytorch/botorch/blob/1518b304f47f5cdbaf9c175e808c90b3a0a6b86d/botorch/acquisition/utils.py#L238

Please provide a minimal, reproducible example of the unexpected behavior.

We train a model on very simple linear problem of maximizing Y_1 = -X, with the constraint that Y_2 = X > 0.81 where the greatest value in the training data is 0.8. We observe that the recommended point is at X=0., which is very far from satisfying the constraint. We would expect the recommended point to satisfy the constraint (since we have enough data to model the constraint function well).

Features of this problem that lead to the buggy behaviour:

• The acqf optimum is at X=0, which is the area searched due to problem (1) above
• The constraint and the objective are competing (the problem would happen if the constraint were Y_2<0.19)

import torch
from botorch.models import SingleTaskGP
from botorch.acquisition import qLogNoisyExpectedImprovement
from botorch.acquisition.objective import GenericMCObjective
from botorch.optim import optimize_acqf
from gpytorch.mlls import ExactMarginalLogLikelihood
from botorch import fit_gpytorch_mll

torch.set_default_dtype(torch.float64)

bounds = torch.tensor([[0.], [1.]])
train_X = torch.linspace(0.2, 0.8, steps=20).unsqueeze(-1)
train_Y = torch.cat((-train_X, train_X), dim=-1)

def constraint(Y: torch.Tensor):
    # Y1 > 0.81
    return -(Y[..., 1] - 0.81)

objective = GenericMCObjective(objective=lambda samples, X: samples[..., 0])

model = SingleTaskGP(train_X, train_Y)
mll = ExactMarginalLogLikelihood(model.likelihood, model)
fit_gpytorch_mll(mll)

acq_function = qLogNoisyExpectedImprovement(
    model=model, 
    X_baseline=train_X,
    constraints=[constraint],
    objective=objective,
    prune_baseline=False, # issue exists for True/False
)

candidates, acq_value = optimize_acqf(
    acq_function,
    bounds,
    q=1,
    num_restarts=20,
    raw_samples=512,
)
print(candidates) # tensor([[0.]]), we expect to be > 0.81
print(acq_function(torch.tensor([[0.90]]))) # tensor([-42.0]), very small!

Please paste any relevant traceback/logs produced by the example provided.

BoTorch Version

0.14.1.dev103+g1518b304f (almost latest on GitHub)

Python Version

3.12.9

Operating System

WSL 2

(Optional) Describe any potential fixes you've considered to the issue outlined above.

1. Fix the dim=0 bug listed in (1) above.
2. Add a constant value to the unconstrained acqf/increase the minimum value. This doesn't change anything when no constraints are present, but biases the optimizer to favour satisfying the constraints when the acqf is zero everywhere. This seems to work okay, but not entirely sure of the downstream effects of doing so.
3. Expose the infeasible_obj argument in compute_best_feasible_objective to the qLogNEI class so that the user can manually provide a value.

Pull Request

Yes

Code of Conduct

• I agree to follow BoTorch's Code of Conduct

[esantorella]

Thanks so much for reporting this!

1. Yeah this looks like a bug. And when N is large, we can wind up with values in X much larger than 1 when those values should like between 0 and 1.

2, 3: Yes, the weights wind up near the middle and not at the edges, especially as the number of points gets large. With the fix from (1), the following gives me a minimum value of 0.4840 and a maximum value of 0.5138:

X = torch.rand(1000, 1)

convex_weights = torch.rand(
    32,
    X.shape[-2],
    dtype=X.dtype,
    device=X.device,
)
weights_sum = convex_weights.sum(dim=1, keepdim=True)
convex_weights = convex_weights / weights_sum
convex_weights @ X

I'm not sure what would be ideal here, but I think almost any way of finding samples from within the search space would be better.

4. I don't have a strong opinion on this, considering prune_baseline can be helpful at least in unconstrained settings -- maybe others can chime in on this.

5. This does look like a bug.

A pull request to address any of these would be appreciated!

[TobyBoyne]

Thanks for the reply! I'm happy to create a PR.

I've been thinking a bit more, and I think this problem of computing _estimate_objective_lower_bound is not easily solved without knowing the bounds of the problem. For example, if all of X_baseline is near the edges of the search space, then there is no need to search beyond the convex hull to find the lowest bound. However, if they are far away, then you would want to search near the boundaries since that is where you'd have the highest variance.

considering prune_baseline can be helpful at least in unconstrained settings

I agree, and I wouldn't want to change any BoTorch defaults without very good reason. Perhaps it would be better to change the behaviour of prune_baseline when all of the points are infeasible (ie. don't prune the "most feasible" point)?

Proposed changes in the PR:

• Instead of random weights, find the bounding box of the points, and uniformly sample within that.
• Check mean +- 6*std to support increasing and decreasing objectives
• Change the behaviour of prune_baseline to keep the "most feasible" point as well as the best objective.
• Add an infeasible_obj argument to the qLogNEI constructor, which allows the user to implicitly solve the "not knowing the bounds of the problem" issue mentioned above. Providing infeasible_obj would also be necessary for when the objective is neither strictly increasing nor decreasing (eg. abs(samples - target_value)), or for mixed discrete spaces where interpolating in baseline_X no longer makes sense.

I was considering somehow using the GP lengthscale to expand the sampling box, or maybe even using some heavier method to bound the value of GP function draws, but I think the acqf is model-agnostic which prevents these approaches.

[esantorella]

When nothing is feasible, qLogProbabilityOfFeasibility should work much better than qLogNEI.

[TobyBoyne]

Thank you for pointing that out, I will start using that! Do you think there should be any warning raised to the user when they use EI with no feasible points, recommending qLogProbabilityOfFeasibility? I couldn't see anything in the docs/issues/discussions making that suggestion.

[Balandat]

Do you think there should be any warning raised to the user when they use EI with no feasible points, recommending qLogProbabilityOfFeasibility?

Yes, I think this would make sense (in Ax we actually dispatch to qLogProbabilityOfFeasibility from a higher level if we observe infeasibility of all points before dropping to botorch).

Maybe you can add this warning as part of #3011?

[TobyBoyne]

TL;DR: The downstream BO impact of this change is massively problem dependent, and most important change is just raising the warning to use ProbabilityOfFeasibility.

Sorry again for so much text!

I've investigated the impact of these changes a bit further. I still recommend the changes in the PR (especially raising a warning to switch to ProbabilityOfFeasibility). However, the downstream impact of the changes on BO is highly variable.

First, as I noted in the original issue report, this behaviour arises when:

1. The objective is large at X=(0, ..., 0), which is the area searched by the buggy convex weight calculation.
2. The constraint and objective are competing.
   If the objective doesn't satisfy these properties, then the poor performance may not arise.

Also, since the buggy behaviour in the weighted sum effectively evaluates the posterior at (0, ..., 0), the lower bound (=mean - 6*std) massively depends whether a point near zero is in train_X. This highlights that maybe the best way to compute the lower bound is to evaluate the model very far from all data, however we don't have access to the domain so we do not know what values are within the problem bounds. If the domain of the problem does not contain (0, ..., 0), then the current behaviour isn't well defined (to me, at least).

One last option could be to change the behaviour of prune_inferior_points. Currently, there is a strong correlation between the quality of the acqf optimization, and the number of points retained after pruning. This is because the bounding box of the pruned baseline is larger, and so can be better sampled. Maybe we could guarantee that at least N points are retained after pruning if all points are infeasible?

Experiments

I have evaluated the change on the Hartmann function, translated such that the optimum is at (0, ..., 0), and constrained such that points are only feasible in a small radius around (0.9, ..., 0.9). This means that the objective and the constraint are competing.

Computing lower bounds

Below is a table avged across ten repeats, of the output of qLogNEI.compute_best_f. The new method more reliably computes a better lower bound. The old method is 0.0 since the function clamps the lower bound to 0.

	Old method	New Method
Prune baseline	0.0	-30.4
No pruning	0.0	-4.412

Downstream BO performance

Here, we show that the improved method finds a feasible value faster than the old method (although both are slower than ProbabilityOfFeasibility). This is also dependent on the point (0, ..., 0) being in the training set. Plot shows median and quartiles across 10 repeats.

For more details, see the notebook here:
improve_best_feasible_objective.ipynb.txt