acerbilab / pybads

PyBADS: Bayesian Adaptive Direct Search optimization algorithm for model fitting in Python
https://acerbilab.github.io/pybads/
BSD 3-Clause "New" or "Revised" License
67 stars 2 forks source link

JOSS review comments #50

Closed vankesteren closed 7 months ago

vankesteren commented 7 months ago

Dear authors,

When reviewing your package here I came up with a few comments / suggestions.

In general, I think your package looks great and has a clear, user-friendly interface. In my testing, it has worked very well. I will recommend acceptance; you may consider the following points as feature requests.

General information about when this method will NOT work

In your documentation and paper, you comment on when you can use PyBADS. In general, for new users it would also be nice to have an idea of when the method is expected to NOT work. An example is this weird step function:

from pybads import BADS
import numpy as np
from pprint import pprint

def crazy_step_fcn(x):
    """Step function which is hard to optimize."""
    if x[0] < 0.5:
        return 1
    elif x[0] > 1.4:
        return 1.1
    return 0.1

plb = np.array([-5, -5])
pub = np.array([5, 5])
x0  = np.array([2, 0])
bads = BADS(crazy_step_fcn, x0, plausible_lower_bounds=plb, plausible_upper_bounds=pub)
pprint(bads.optimize())
{
 'algorithm': 'Bayesian adaptive direct search',
 'fsd': 0.0,
 'fun': <function crazy_step_fcn at 0x000001C383557F60>,
 'func_count': 89,
 'fval': 1.0,
 'iterations': 5,
 'mesh_size': 0.000244140625,
 'message': 'Optimization terminated: change in the function value less than '
            "options['tol_fun'].",
 'non_box_cons': None,
 'overhead': 6985816.8440025635,
 'problem_type': 'unconstrained',
 'random_seed': None,
 'success': True,
 'target_type': 'deterministic',
 'total_time': 0.6077661524282226,
 'version': '1.0.3',
 'x': array([-4.82910156, -4.50195312]),
 'x0': array([[2, 0]]),
 'ysd_vec': None,
 'yval_vec': None
}

Better univariate interface

For univariate problems (as in the example in the previous section, where as a workaround I just added a random second dimension) it would be great to have a simple interface where you don't need to specify a row vector, but rather a single variable, like so:

rfn = lambda x: x^2 + 3.2 + np.random.normal(scale=0.1)
plb = np.array([-5, -5])
pub = np.array([5, 5])
x0 = np.array([10])
opt = BADS(rfn, x0, plausible_lower_bounds=plb, plausible_upper_bounds=pub)

currently, this gives the following error:

ValueError: All input vectors (x0, lower_bounds, upper_bounds,
                 plausible_lower_bounds, plausible_upper_bounds), if specified,
                 need to be row vectors with D elements.

but x0 IS a row vector! Would be even nicer to just allow a float value rather than an np.array.

Access to GP surrogate model in OptimizeResult

It would be great to be able to inspect the optimization result in a bit more detail by looking at the last / best GP surrogate model. The dependency gpyreg seems to have a nice plot method, would be great to have access to this method by passing the last gpyreg.GP object through.

lacerbi commented 7 months ago

Hi @vankesteren, thanks for the review and the helpful comments!

General information about when this method will NOT work

Good point, we will add a section in the docs. It is true that we do assume some degree of nicety in the target. Importantly, we do not claim that this is a global optimization algorithm (when classifying optimization algorithms for my students, I put BADS in the 'semi-local' family together with e.g. CMA-ES).

Better univariate interface

Excellent catch. @GurjeetSinghSangra is having a look, we should get this to work.

Access to GP surrogate model in OptimizeResult

Actually, this is omitted on purpose. The GP surrogate model built by (Py)BADS is a local model and can often be a terrible representation of the target. Point is, this is not what matters for BADS: the role of the surrogate model is to guide the optimization in the right direction at each iteration along the optimization trajectory. There is absolutely no requirement for it to be a good model of the target (even locally, let alone globally), it just needs to be "good enough". On top of that, the final (local) GP model tends to be particularly not representative of the target since it is based on a set of local points in the very close neighborhood of the final point. The final GP can be anything from a constant to something sensible close to the point but crazy outside that neighborhood, etc.

In sum, for all these reasons, returning the GP can be useless or even misleading for the user; so it was a deliberate decision at the time of the original bads - and now PyBADS - not to return the surrogate. You can probably still go and grab it if you are an advanced user, but the barrier is on purpose - the typical user should not even think about looking at the GP...

GurjeetSinghSangra commented 7 months ago

Hello @vankesteren, thanks again for concluding the review procedure, and highlighting this issue.

Better univariate interface

By looking at the problem I found that the error output message was a bit misleading in this case. For the given example, PyBADS should instead tell the user about the incoherent dimension between bounds and the given starting point $x_0$.

More precisely, in your example:

rfn = lambda x: x**2 + 3.2 + np.random.normal(scale=0.1)
plb = np.array([-5, -5])
pub = np.array([5, 5])
x0 = np.array([10])
opt = BADS(rfn, x0, plausible_lower_bounds=plb, plausible_upper_bounds=pub)

the provided bounds are 2-D vectors, and the starting is a scalar. Therefore, it is unclear to PyBADS the dimension of the optimization problem. The correct code for the univariate should use one element vector or scalar values, i.e

plb = -5 # or np.array([-5])
pub = 5 # or np.array([5])
x0 = 10 # or np.array([10])
opt = BADS(rfn, x0, plausible_lower_bounds=plb, plausible_upper_bounds=pub)

The current version of PyBADS is already handling the univariate case and even scalar value. I tested it again and it passed all the tests. See for example the following screenshot:

image

However, I would like to thank you for your review and for pointing out the misleading error message. I just turned it into a more useful and friendly message, allowing the user to inspect better the problem configuration. More precisely with the latest change in the code, if the dimensions between the bounds and starting point $x_0$ do not match, PyBADS will return the following message:

All input vectors (lower_bounds, upper_bounds,
                 plausible_lower_bounds, plausible_upper_bounds), if specified,
                 need to be of dimension D=1 as the starting point x_0=[[10]].
vankesteren commented 7 months ago

Ah, great, yes my mistake. Nice that the error message is more informative now!

Feel free to close this issue as my comments have been addressed.