From 44d9b1f429b58eb5425a540aab25edd6bcfffb23 Mon Sep 17 00:00:00 2001 From: John Jakeman Date: Fri, 2 Aug 2024 17:09:56 -0600 Subject: [PATCH] Add Hyperparameter.set_active_indices() instead of determining active indices by looking for np.nan in bounds --- pyapprox/surrogates/autogp/exactgp.py | 2 +- .../autogp/tests/test_gaussian_process.py | 14 ++-- pyapprox/surrogates/bases/functiontrain.py | 15 ++--- .../bases/tests/test_functiontrain.py | 4 +- pyapprox/surrogates/kernels/kernels.py | 9 ++- pyapprox/util/hyperparameter.py | 65 ++++++++++++++----- pyapprox/util/tests/test_hyperparameter.py | 16 ++++- 7 files changed, 88 insertions(+), 37 deletions(-) diff --git a/pyapprox/surrogates/autogp/exactgp.py b/pyapprox/surrogates/autogp/exactgp.py index c3ff909a..49079cef 100644 --- a/pyapprox/surrogates/autogp/exactgp.py +++ b/pyapprox/surrogates/autogp/exactgp.py @@ -1,4 +1,4 @@ -bfrom abc import ABC, abstractmethod +from abc import ABC, abstractmethod from typing import Tuple import warnings diff --git a/pyapprox/surrogates/autogp/tests/test_gaussian_process.py b/pyapprox/surrogates/autogp/tests/test_gaussian_process.py index bee66cbd..875842bd 100644 --- a/pyapprox/surrogates/autogp/tests/test_gaussian_process.py +++ b/pyapprox/surrogates/autogp/tests/test_gaussian_process.py @@ -257,11 +257,11 @@ def test_compare_with_deprecated_gp(self): sigma = 1 lenscale = 0.5 kernel = ConstantKernel( - sigma, [np.nan, np.nan], backend=bkd + sigma, [0.01, 1], fixed=True, backend=bkd ) * MaternKernel( - np.inf, lenscale, [np.nan, np.nan], nvars, backend=bkd + np.inf, lenscale, [lenscale, lenscale], nvars, backend=bkd ) + GaussianNoiseKernel( - noise, [np.nan, np.nan], backend=bkd + noise, [0.02, 2], fixed=True, backend=bkd ) gp = ExactGaussianProcess(nvars, kernel) @@ -414,7 +414,7 @@ def fun(xx): exact_gp = ExactGaussianProcess( nvars, kernel - + GaussianNoiseKernel(noise_var, [np.nan, np.nan], backend=bkd), + + GaussianNoiseKernel(noise_var, [0.1, 1], fixed=True, backend=bkd), trend=None, values_trans=values_trans, kernel_reg=0, @@ -431,17 +431,19 @@ def fun(xx): "noise_std", 1, np.sqrt(noise_var), - [np.nan, np.nan], + [0.1, 1], LogHyperParameterTransform(backend=bkd), + fixed=True, ) inducing_samples = InducingSamples( nvars, ninducing_samples, inducing_samples=inducing_samples, - inducing_sample_bounds=bkd._la_atleast1d([np.nan, np.nan]), + inducing_sample_bounds=bkd._la_atleast1d([-1, 1]), noise=noise, backend=bkd, ) + inducing_samples.hyp_list.set_all_inactive() values_trans = IdentityTransform(backend=bkd) # use correlation length learnt by exact gp vi_kernel = kernel diff --git a/pyapprox/surrogates/bases/functiontrain.py b/pyapprox/surrogates/bases/functiontrain.py index 8ed04fb9..2fab09e4 100644 --- a/pyapprox/surrogates/bases/functiontrain.py +++ b/pyapprox/surrogates/bases/functiontrain.py @@ -1,9 +1,6 @@ import copy from abc import ABC, abstractmethod -# import numpy for np.nan -import numpy as np - from pyapprox.surrogates.bases.basis import MonomialBasis from pyapprox.surrogates.bases.basisexp import Regressor, BasisExpansion, MonomialExpansion @@ -72,14 +69,14 @@ def _core_params_eval(self, active_opt_params): def _core_jacobian_ad(self, samples, core_id): """Compute Jacobian with automatic differentiation.""" self._samples = samples - hyplist_bounds = self._bkd._la_copy(self.hyp_list.get_bounds()) + active_indices = self.hyp_list.get_active_indices() for ii in range(self.nvars()): if ii != core_id: - self._cores[ii].hyp_list.set_bounds([np.nan, np.nan]) + self._cores[ii].hyp_list.set_all_inactive() jac = self._bkd._la_jacobian( self._core_params_eval, self.hyp_list.get_active_opt_params() ) - self.hyp_list.set_bounds(hyplist_bounds) + self.hyp_list.set_active_indices(active_indices) # create list of jacobians for each qoi # jac will be of shape (nsamples, nqoi, ntotal_core_active_parameters) # so if using basis expansion with same expansion for each core and nqoi = 1 @@ -213,11 +210,11 @@ def __init__(self, univariate_funs : list[BasisExpansion], nqoi: int): def _init_constant_fun(self, fill_value, nqoi): constant_basis = MonomialBasis(backend=self._bkd) constant_basis.set_hyperbolic_indices(1, 0, 1.) - # set coef_bounds to [np.nan, np.nan] so these coefficients are set as - # inactive fun = MonomialExpansion( - constant_basis, coef_bounds=[np.nan, np.nan], nqoi=nqoi + constant_basis, coef_bounds=[fill_value, fill_value], nqoi=nqoi ) + # set as inactive so they cannot be changed during optimization + fun.hyp_list.set_all_inactive() fun.set_coefficients(self._bkd._la_full((1, nqoi), fill_value)) return fun diff --git a/pyapprox/surrogates/bases/tests/test_functiontrain.py b/pyapprox/surrogates/bases/tests/test_functiontrain.py index a59b2b9d..0d096ccf 100644 --- a/pyapprox/surrogates/bases/tests/test_functiontrain.py +++ b/pyapprox/surrogates/bases/tests/test_functiontrain.py @@ -38,13 +38,14 @@ def test_additive_function_train(self): train_values = sum(univariate_vals) assert np.allclose(train_values, ft_vals) - ft.hyp_list.set_bounds([-np.inf,np.inf]) + ft.hyp_list.set_all_active() jac_ans = [ ft._core_jacobian(train_samples, ii) for ii in range(ft.nvars()) ] if bkd._la_jacobian_implemented(): for ii in range(ft.nvars()): for qq in range(ft.nqoi()): + print(ii) assert np.allclose( ft._core_jacobian_ad(train_samples, ii)[qq], jac_ans[ii][qq], atol=1e-14 ) @@ -56,6 +57,7 @@ def test_additive_function_train(self): params = bkd._la_asarray(np.random.normal(0, 1, (ft.hyp_list.nactive_vars(),))) ft.hyp_list.set_values(params) + ft.hyp_list.set_all_active() solver = AlternatingLeastSquaresSolver(verbosity=2) solver.solve(ft, train_samples, train_values) ft_vals = ft(train_samples) diff --git a/pyapprox/surrogates/kernels/kernels.py b/pyapprox/surrogates/kernels/kernels.py index 3c5eb8e7..d6b5e85e 100644 --- a/pyapprox/surrogates/kernels/kernels.py +++ b/pyapprox/surrogates/kernels/kernels.py @@ -121,6 +121,7 @@ def __init__( lenscale, lenscale_bounds, nvars: int, + fixed: bool = False, backend: LinAlgMixin = None, ): """The matern kernel for varying levels of smoothness.""" @@ -134,6 +135,7 @@ def __init__( lenscale, lenscale_bounds, transform, + fixed=fixed, backend=self._bkd, ) self.hyp_list = HyperParameterList([self._lenscale]) @@ -169,7 +171,7 @@ def nvars(self): class ConstantKernel(Kernel): def __init__( - self, constant, constant_bounds=None, transform=None, backend=None + self, constant, constant_bounds=None, transform=None, fixed=False, backend=None ): if backend is None and transform is not None: backend = transform._bkd @@ -179,7 +181,7 @@ def __init__( if constant_bounds is None: constant_bounds = [-self._bkd._la_inf(), self._bkd._la_inf()] self._const = HyperParameter( - "const", 1, constant, constant_bounds, transform, backend=self._bkd + "const", 1, constant, constant_bounds, transform, fixed=fixed, backend=self._bkd ) self.hyp_list = HyperParameterList([self._const]) @@ -205,7 +207,7 @@ def __call__(self, X1, X2=None): class GaussianNoiseKernel(Kernel): - def __init__(self, constant, constant_bounds=None, backend=None): + def __init__(self, constant, constant_bounds=None, fixed=False, backend=None): super().__init__(backend) self._const = HyperParameter( "const", @@ -213,6 +215,7 @@ def __init__(self, constant, constant_bounds=None, backend=None): constant, constant_bounds, LogHyperParameterTransform(backend=self._bkd), + fixed=fixed, backend=self._bkd, ) self.hyp_list = HyperParameterList([self._const]) diff --git a/pyapprox/util/hyperparameter.py b/pyapprox/util/hyperparameter.py index 1fb2e701..ca45eb93 100644 --- a/pyapprox/util/hyperparameter.py +++ b/pyapprox/util/hyperparameter.py @@ -48,6 +48,7 @@ def __init__( values, bounds, transform: HyperParameterTransform = None, + fixed: bool = False, backend: LinAlgMixin = None, ): """A possibly vector-valued hyper-parameter to be used with @@ -79,8 +80,32 @@ def __init__( self._values.shape, self.nvars() ) ) + if fixed: + self.set_all_inactive() + else: + self.set_all_active() self.set_bounds(bounds) + def set_active_indices(self, indices): + if indices.shape[0] == 0: + self._active_indices = indices + return + + if max(indices) >= self.nvars(): + raise ValueError("indices exceed nvars") + if min(indices) < 0: + raise ValueError("Ensure indices >= 0") + self._active_indices = indices + + def get_active_indices(self): + return self._active_indices + + def set_all_inactive(self): + self.set_active_indices(self._bkd._la_zeros((0, ), dtype=int)) + + def set_all_active(self): + self.set_active_indices(self._bkd._la_arange(self.nvars(), dtype=int)) + def set_bounds(self, bounds): self.bounds = self._bkd._la_atleast1d(bounds) if self.bounds.shape[0] == 2: @@ -93,21 +118,6 @@ def set_bounds(self, bounds): self.bounds = self._bkd._la_reshape( self.bounds, (self.bounds.shape[0] // 2, 2) ) - if ( - self._bkd._la_where( - (self._values < self.bounds[:, 0]) - | (self._values > self.bounds[:, 1]) - )[0].shape[0] - > 0 - ): - raise ValueError("values outside bounds") - self._active_indices = self._bkd._la_tointeger( - self._bkd._la_atleast1d( - self._bkd._la_arange(self.nvars())[ - ~self._bkd._la_isnan(self.bounds[:, 0]) - ] - ) - ) def nvars(self): """Return the number of hyperparameters.""" @@ -236,13 +246,36 @@ def get_bounds(self): [hyp.get_bounds() for hyp in self.hyper_params] ) - def get_values(self): """Get the values of the parameters in the user space.""" return self._bkd._la_hstack( [hyp.get_values() for hyp in self.hyper_params] ) + def get_active_indices(self): + cnt = 0 + active_indices = [] + for hyp in self.hyper_params: + active_indices.append(hyp.get_active_indices()+cnt) + cnt += hyp.nvars() + return self._bkd._la_hstack(active_indices) + + def set_active_indices(self, active_indices): + cnt = 0 + for hyp in self.hyper_params: + hyp_indices = self._bkd._la_where( + (active_indices>=cnt)&(active_indices