gpr#

This module contains the Gaussia Process Regressor class. For now, the only implemented one is a modified version of the sklearn implementation, with some tweaks to allow disabling input validation and implementing more control at the hyperparameter optimization stage.

At the moment, the GPR kernel is a product of a constant kernel and an anisotropic length-correlation one. Noise is treated as uncorrelated standard deviations for inputs, and thus simply added to the kernel matrix diagonal.

class gpry.gpr.GaussianProcessRegressor(kernel='RBF', output_scale_prior=[0.01, 1000.0], length_scale_prior=[0.001, 10.0], noise_level=0.01, noise_fixed=True, optimizer='fmin_l_bfgs_b', n_restarts_optimizer=0, random_state=None)[source]#

Bases: GaussianProcessRegressor

Modified version of the GaussianProcessRegressor of sklearn.

The implementation is based on Algorithm 2.1 of Gaussian Processes for Machine Learning (GPML) by Rasmussen and Williams.

This modified interface provides, in addition to the sklearn-GPR:

  • Re-implements the fit method to allow for more control in the noise (alpha) update and the hyperparameter optimization.

  • Implements derivative return values in the predict method, as well as a predict_std method to return the standard deviation of the target only (useful for acquisition).

  • In the relevant methods, exposes flags to disable input data validation, for an additional speed boost.

Parameters:
  • kernel (kernel object, string, dict, optional (default: "RBF")) –

    The kernel specifying the covariance function of the GP. If “RBF”/”Matern” is passed, the asymmetric kernel:

    ConstantKernel() * RBF/Matern()
    

    is used as default, where the length correlation kernel is assumed anisotropic if a list of bounds is passed using the length_scale_prior argument. To pass different arguments to the kernel, e.g. nu=5/2 for Matern, pass a single-key dict as {"Matern": {"nu": 2.5}}". Note that the kernel’s hyperparameters are optimized during fitting.

  • output_scale_prior (tuple as (min, max), optional (default: [1e-2, 1e3])) – Prior for the (non-squared) scale parameter, in normalised logp units.

  • length_scale_prior (tuple as (min, max), optional (default: [1e-3, 1e1])) – Prior for the length parameters, as a fraction of the parameter priors sizes.

  • noise_level (float or array-like, optional (default: 1e-2)) – Square-root of the value added to the diagonal of the kernel matrix during fitting. Larger values correspond to increased noise level in the observations and reduce potential numerical issue during fitting. If an array is passed, it must have the same number of entries as the data used for fitting and is used as datapoint-dependent noise level. Note that this is equivalent to adding a WhiteKernel with c=noise_level.

  • noise_fixed (bool (default: True)) – Whether the noise is treated as a fixed diagonal offset in the kernel matrix, or an additive white noise term in the kernel. In the latter case, the value given as noise_level is set as the upper bound.

  • optimizer (str or callable, optional (default: "fmin_l_bfgs_b")) –

    Can either be one of the internally supported optimizers for optimizing the kernel’s parameters, specified by a string, or an externally defined optimizer passed as a callable. If a callable is passed, it must have the signature:

    def optimizer(obj_func, initial_theta, bounds):
        # * 'obj_func' is the objective function to be maximized, which
        #   takes the hyperparameters theta as parameter and an
        #   optional flag eval_gradient, which determines if the
        #   gradient is returned additionally to the function value
        # * 'initial_theta': the initial value for theta, which can be
        #   used by local optimizers
        # * 'bounds': the bounds on the values of theta
        ....
        # Returned are the best found hyperparameters theta and
        # the corresponding value of the target function.
        return theta_opt, func_min
    

    Per default, the ‘fmin_l_bfgs_b’ algorithm from scipy.optimize is used. If None is passed, the kernel’s parameters are kept fixed. Available internal optimizers are:

    'fmin_l_bfgs_b'
    

  • n_restarts_optimizer (int, optional (default: 0)) – The number of restarts of the optimizer for finding the kernel’s parameters which maximize the log-marginal likelihood. The first run of the optimizer is performed from the kernel’s initial parameters, the remaining ones (if any) from thetas sampled log-uniform randomly from the space of allowed theta-values. If greater than 0, all bounds must be finite. Note that n_restarts_optimizer == 0 implies that one run is performed.

  • random_state (int or numpy.random.Generator, optional) – The generator used to perform random operations of the GPR. If an integer is given, it is used as a seed for the default global numpy random number generator.

Attributes:
  • X_train_ (array-like, shape = (n_samples, n_features)) – (Possibly transformed) feature values in training data of the GPR (also required for prediction). Mostly intended for internal use.

  • y_train_ (array-like, shape = (n_samples, [n_output_dims])) – (Possibly transformed) target values in training data of the GPR (also required for prediction). Mostly intended for internal use.

  • alpha (array-like, shape = (n_samples, [n_output_dims]) or scalar) – The value which is added to the diagonal of the kernel. This is the square of noise_level.

  • kernel_ (kernels object) – The kernel used for prediction. The structure of the kernel is the same as the one passed as parameter but with optimized hyperparameters.

  • alpha_ (array-like, shape = (n_samples, n_samples)) – Not to be confused with alpha! The inverse Kernel matrix of the training points multiplied with y_train_ (Dual coefficients of training data points in kernel space). Needed at prediction.

  • V_ (array-like, shape = (n_samples, n_samples)) – Lower-triangular Cholesky decomposition of the inverse kernel in X_train_

  • log_marginal_likelihood_value_ (float) – The log-marginal-likelihood of self.kernel_.theta

  • scales (tuple) – Kernel scales as (output_scale, (length_scale_1, ...))

property scales#

Kernel scales as (output_scale, (length_scale_1, ...)).

property noise_level#

Kernel noise level (not squared).

fit(X, y, noise_level=None, fit_hyperparameters=True, validate=True)[source]#

Re-implementation of the sk GPR fit method, that allows for updating the noise level (as alpha), and exposes flags for input validation and hyperparameter fitting.

If hyperparameters are kept constant, fitting here refers to the re-calculation of the GPR inverse matrix \((K(X,X)+\sigma_n^2 I)^{-1}\) which is needed for predictions.

The highest cost incurred by this method is the refitting of the GPR kernel hyperparameters \(\theta\). It can be useful to disable it (fit_hyperparameters=False) in cases where it is worth saving the computational expense in exchange for a loss of information, such as when performing parallelized active sampling (NB: this is only possible when the GPR hyperparameters have been fit at least once).

If called with X=None, y=None, it re-fits the model without adding new points.

Parameters:
  • X (array-like, shape = (n_samples, n_features), or None) – Training data to append to the model.

  • y (array-like, shape = (n_samples, [n_output_dims]), or None) – Target values to append to the data

  • noise_level (number, array-like, shape = (n_samples, [n_output_dims])) – Uncorrelated standard deviation(s) to add to the diagonal part of the covariance matrix. If an array is passed, it needs to have the same number of entries as y. If None, the noise_level set in the instance is used (notice that passing None is not recommended if a preprocessor for the GPR training samples may have been refit after initialization). If a non-None value is passed, it is advisable to refit the hyperparameters of the kernel.

  • fit_hyperparameters (Bool or dict (default: True)) – Whether the GPR \(\theta\)-parameters should be optimised. It can be passed a dict, which is always interpreted as True, containing arguments for the _fit_hyperparameters method, namely n_restarts of the optimizer, start_from_current if the first optimizer run should start from the last optimum, or different hyperparameters_bounds. E.g. to perform a single GPR optimization run starting at the last hyperparameter optimum: fit_hyperparameters={'start_from_current': True, 'n_restarts': 1}.

  • validate (bool, default: True) – If False, X and y are assumed to be correctly formatted, and no checks are performed on them. Reduces overhead. Use only for repeated calls when the input is programmatically generated to be correct at each stage.

Returns:

self – GaussianProcessRegressor class instance.

Return type:

object

log_marginal_likelihood(*args, **kwargs)[source]#

Log-marginal likelihood of the kernel hyperparameters given the training data.

predict(X, return_std=False, return_mean_grad=False, return_std_grad=False, validate=True)[source]#

Predict output for X.

Reimplementation of the sk-learn GPR method: in addition to the mean of the predictive distribution, also its standard deviation (return_std=True), the gradient of the mean and the standard-deviation with respect to X can be optionally provided.

Parameters:
  • X (array-like, shape = (n_samples, n_features)) – Query points where the GP is evaluated.

  • return_std (bool, default: False) – If True, the standard-deviation of the predictive distribution at the query points is returned along with the mean.

  • return_mean_grad (bool, default: False) – Whether or not to return the gradient of the mean. Only valid when X is a single point.

  • return_std_grad (bool, default: False) – Whether or not to return the gradient of the std. Only valid when X is a single point.

  • validate (bool, default: True) – If False, X and y are assumed to be correctly formatted, and no checks are performed on them. Reduces overhead. Use only for repeated calls when the input is programmatically generated to be correct at each stage.

  • note:: (..) – Note that in contrast to the sklearn GP Regressor our implementation cannot return the full covariance matrix. This is to save on some complexity and since the full covariance cannot be calculated if either values are infinite or a y-preprocessor is used.

Returns:

  • y_mean (array, shape = (n_samples, [n_output_dims])) – Mean of predictive distribution a query points

  • y_std (array, shape = (n_samples,), optional) – Standard deviation of predictive distribution at query points. Only returned when return_std is True.

  • y_mean_grad (shape = (n_samples, n_features), optional) – The gradient of the predicted mean.

  • y_std_grad (shape = (n_samples, n_features), optional) – The gradient of the predicted std.

predict_std(X, validate=True)[source]#

Predict output standart deviation for X.

Parameters:
  • X (array-like, shape = (n_samples, n_features)) – Query points where the GP is evaluated.

  • validate (bool, default: True) – If False, X and y are assumed to be correctly formatted, and no checks are performed on them. Reduces overhead. Use only for repeated calls when the input is programmatically generated to be correct at each stage.

Returns:

y_std – Standard deviation of predictive distribution at query points. Only returned when return_std is True.

Return type:

array, shape = (n_samples,), optional