From 8f19732cdddaac5117354fe456c1109451d26bcf Mon Sep 17 00:00:00 2001 From: Ulpu Remes Date: Wed, 13 Mar 2024 22:23:34 +0200 Subject: [PATCH 1/6] sample randmaxvar batches --- elfi/methods/bo/acquisition.py | 18 +++++++++++++----- 1 file changed, 13 insertions(+), 5 deletions(-) diff --git a/elfi/methods/bo/acquisition.py b/elfi/methods/bo/acquisition.py index f4062465..8b0a645c 100644 --- a/elfi/methods/bo/acquisition.py +++ b/elfi/methods/bo/acquisition.py @@ -492,7 +492,7 @@ class RandMaxVar(MaxVar): """ - def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, + def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, warmup=None, limit_faulty_init=10, sigma_proposals=None, *args, **opts): """Initialise RandMaxVar. @@ -504,6 +504,8 @@ def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, Name of the sampler (options: metropolis, nuts). n_samples : int, optional Length of the sampler's chain for obtaining the acquisitions. + warmup : int, optional + Number of samples discarded as warmup. Defaults to n_samples/2. limit_faulty_init : int, optional Limit for the iterations used to obtain the sampler's initial points. sigma_proposals : dict, optional @@ -515,6 +517,7 @@ def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, self.name = 'rand_max_var' self.name_sampler = sampler self._n_samples = n_samples + self._warmup = warmup or n_samples // 2 self._limit_faulty_init = limit_faulty_init if self.name_sampler == 'metropolis': self._sigma_proposals = resolve_sigmas(self.model.parameter_names, @@ -538,8 +541,8 @@ def acquire(self, n, t=None): """ if n > self._n_samples: - raise ValueError(("The number of acquisitions ({0}) has to be lower " - "than the number of the samples ({1}).").format(n, self._n_samples)) + raise ValueError(("The number of acquisitions ({0}) has to be lower than the number " + "of the samples ({1}).").format(n, self._n_samples - self._warmup)) logger.debug('Acquiring the next batch of %d values', n) gp = self.model @@ -593,8 +596,13 @@ def _evaluate_logpdf(theta): raise ValueError( "Incompatible sampler. Please check the options in the documentation.") - # Using the last n points of the MH chain for the acquisition batch. - batch_theta = samples[-n:, :] + if n > 1: + # Remove warmup samples and return n random points + samples = samples[self._warmup:] + batch_theta = self.random_state.permutation(samples)[:n] + else: + # Return the last point + batch_theta = samples[-1:] break return batch_theta From 53362c49f83900ce2ff83323590835d96209146d Mon Sep 17 00:00:00 2001 From: Ulpu Remes Date: Wed, 13 Mar 2024 22:25:32 +0200 Subject: [PATCH 2/6] sample candidate initial points from prior --- elfi/methods/bo/acquisition.py | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/elfi/methods/bo/acquisition.py b/elfi/methods/bo/acquisition.py index 8b0a645c..c2ff181c 100644 --- a/elfi/methods/bo/acquisition.py +++ b/elfi/methods/bo/acquisition.py @@ -493,7 +493,7 @@ class RandMaxVar(MaxVar): """ def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, warmup=None, - limit_faulty_init=10, sigma_proposals=None, *args, **opts): + limit_faulty_init=1000, sigma_proposals=None, *args, **opts): """Initialise RandMaxVar. Parameters @@ -571,9 +571,14 @@ def _evaluate_logpdf(theta): raise SystemExit("Unable to find a suitable initial point.") # Proposing the initial point. - theta_init = np.zeros(shape=len(gp.bounds)) - for idx_param, range_bound in enumerate(gp.bounds): - theta_init[idx_param] = self.random_state.uniform(range_bound[0], range_bound[1]) + if self.prior is None: + theta_init = np.zeros(shape=len(gp.bounds)) + for idx_param, bound in enumerate(gp.bounds): + theta_init[idx_param] = self.random_state.uniform(bound[0], bound[1]) + else: + theta_init = self.prior.rvs(random_state=self.random_state) + for idx_param, bound in enumerate(gp.bounds): + theta_init[idx_param] = np.clip(theta_init[idx_param], bound[0], bound[1]) # Refusing to accept a faulty initial point. if np.isinf(_evaluate_logpdf(theta_init)): From 20234b72f186db5bf056f10a62535497fa7bacdd Mon Sep 17 00:00:00 2001 From: Ulpu Remes Date: Thu, 14 Mar 2024 00:47:15 +0200 Subject: [PATCH 3/6] update changelog --- CHANGELOG.rst | 1 + 1 file changed, 1 insertion(+) diff --git a/CHANGELOG.rst b/CHANGELOG.rst index 10975429..9313c10e 100644 --- a/CHANGELOG.rst +++ b/CHANGELOG.rst @@ -1,6 +1,7 @@ Changelog ========= +- Improve `randmaxvar` initialisation and batch acquisitions - Enable using `maxiter` in `bo.utils.minimize` - Fix surrogate model copy operation - Fix typo in requirements.txt From 209c1c3c3f1695c6b5cb92569ed027907ccd2527 Mon Sep 17 00:00:00 2001 From: Ulpu Remes Date: Thu, 14 Mar 2024 00:47:59 +0200 Subject: [PATCH 4/6] fix whitespace error --- elfi/methods/bo/utils.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/elfi/methods/bo/utils.py b/elfi/methods/bo/utils.py index c7f9b866..45ec826b 100644 --- a/elfi/methods/bo/utils.py +++ b/elfi/methods/bo/utils.py @@ -97,7 +97,7 @@ def minimize(fun, for i in range(n_start_points): result = scipy.optimize.minimize(fun, start_points[i, :], method=method, jac=grad, - bounds=bounds, constraints=constraints, + bounds=bounds, constraints=constraints, options={'maxiter': maxiter}) locs.append(result['x']) vals[i] = result['fun'] From 432933e73723c8c0e48f06807ece4d1e580b2f7f Mon Sep 17 00:00:00 2001 From: Ulpu Remes Date: Thu, 14 Mar 2024 20:29:55 +0200 Subject: [PATCH 5/6] update sampling candidate points --- elfi/methods/bo/acquisition.py | 18 ++++++++++++------ 1 file changed, 12 insertions(+), 6 deletions(-) diff --git a/elfi/methods/bo/acquisition.py b/elfi/methods/bo/acquisition.py index c2ff181c..1a7c09d2 100644 --- a/elfi/methods/bo/acquisition.py +++ b/elfi/methods/bo/acquisition.py @@ -493,7 +493,8 @@ class RandMaxVar(MaxVar): """ def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, warmup=None, - limit_faulty_init=1000, sigma_proposals=None, *args, **opts): + limit_faulty_init=1000, init_from_prior=False, sigma_proposals=None, + *args, **opts): """Initialise RandMaxVar. Parameters @@ -508,6 +509,9 @@ def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, warmup=None, Number of samples discarded as warmup. Defaults to n_samples/2. limit_faulty_init : int, optional Limit for the iterations used to obtain the sampler's initial points. + init_from_prior : bool, optional + Controls whether the sampler's initial points are sampled from the prior or + a uniform distribution within model bounds. Defaults to model bounds. sigma_proposals : dict, optional Standard deviations for Gaussian proposals of each parameter for Metropolis Markov Chain sampler. Defaults to 1/10 of surrogate model bound lengths. @@ -519,6 +523,7 @@ def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, warmup=None, self._n_samples = n_samples self._warmup = warmup or n_samples // 2 self._limit_faulty_init = limit_faulty_init + self._init_from_prior = init_from_prior if self.name_sampler == 'metropolis': self._sigma_proposals = resolve_sigmas(self.model.parameter_names, sigma_proposals, @@ -571,15 +576,16 @@ def _evaluate_logpdf(theta): raise SystemExit("Unable to find a suitable initial point.") # Proposing the initial point. - if self.prior is None: - theta_init = np.zeros(shape=len(gp.bounds)) - for idx_param, bound in enumerate(gp.bounds): - theta_init[idx_param] = self.random_state.uniform(bound[0], bound[1]) - else: + if self._init_from_prior: theta_init = self.prior.rvs(random_state=self.random_state) for idx_param, bound in enumerate(gp.bounds): theta_init[idx_param] = np.clip(theta_init[idx_param], bound[0], bound[1]) + else: + theta_init = np.zeros(shape=len(gp.bounds)) + for idx_param, bound in enumerate(gp.bounds): + theta_init[idx_param] = self.random_state.uniform(bound[0], bound[1]) + # Refusing to accept a faulty initial point. if np.isinf(_evaluate_logpdf(theta_init)): continue From 6ee5a7e894ffe8f73a036fa190d89763b68a18bd Mon Sep 17 00:00:00 2001 From: Ulpu Remes Date: Thu, 14 Mar 2024 21:29:27 +0200 Subject: [PATCH 6/6] fix maxvar initialisation --- elfi/methods/bo/acquisition.py | 29 ++++++++++++++++++++--------- 1 file changed, 20 insertions(+), 9 deletions(-) diff --git a/elfi/methods/bo/acquisition.py b/elfi/methods/bo/acquisition.py index 1a7c09d2..fa94bba7 100644 --- a/elfi/methods/bo/acquisition.py +++ b/elfi/methods/bo/acquisition.py @@ -327,16 +327,20 @@ class MaxVar(AcquisitionBase): """ - def __init__(self, quantile_eps=.01, *args, **opts): + def __init__(self, model, prior, quantile_eps=.01, **opts): """Initialise MaxVar. Parameters ---------- + model : elfi.GPyRegression + Gaussian process model used to calculate the unnormalised approximate likelihood. + prior : scipy-like distribution + Prior distribution. quantile_eps : int, optional Quantile of the observed discrepancies used in setting the ABC threshold. """ - super(MaxVar, self).__init__(*args, **opts) + super(MaxVar, self).__init__(model, prior=prior, **opts) self.name = 'max_var' self.label_fn = 'Variance of the Unnormalised Approximate Posterior' self.quantile_eps = quantile_eps @@ -492,13 +496,16 @@ class RandMaxVar(MaxVar): """ - def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, warmup=None, - limit_faulty_init=1000, init_from_prior=False, sigma_proposals=None, - *args, **opts): + def __init__(self, model, prior, quantile_eps=.01, sampler='nuts', n_samples=50, warmup=None, + limit_faulty_init=1000, init_from_prior=False, sigma_proposals=None, **opts): """Initialise RandMaxVar. Parameters ---------- + model : elfi.GPyRegression + Gaussian process model used to calculate the unnormalised approximate likelihood. + prior : scipy-like distribution + Prior distribution. quantile_eps : int, optional Quantile of the observed discrepancies used in setting the ABC threshold. sampler : string, optional @@ -517,7 +524,7 @@ def __init__(self, quantile_eps=.01, sampler='nuts', n_samples=50, warmup=None, Markov Chain sampler. Defaults to 1/10 of surrogate model bound lengths. """ - super(RandMaxVar, self).__init__(quantile_eps, *args, **opts) + super(RandMaxVar, self).__init__(model, prior, quantile_eps, **opts) self.name = 'rand_max_var' self.name_sampler = sampler self._n_samples = n_samples @@ -648,13 +655,17 @@ class ExpIntVar(MaxVar): """ - def __init__(self, quantile_eps=.01, integration='grid', d_grid=.2, + def __init__(self, model, prior, quantile_eps=.01, integration='grid', d_grid=.2, n_samples_imp=100, iter_imp=2, sampler='nuts', n_samples=2000, - sigma_proposals=None, *args, **opts): + sigma_proposals=None, **opts): """Initialise ExpIntVar. Parameters ---------- + model : elfi.GPyRegression + Gaussian process model used to calculate the approximate unnormalised likelihood. + prior : scipy-like distribution + Prior distribution. quantile_eps : int, optional Quantile of the observed discrepancies used in setting the discrepancy threshold. integration : str, optional @@ -680,7 +691,7 @@ def __init__(self, quantile_eps=.01, integration='grid', d_grid=.2, Markov Chain sampler. Defaults to 1/10 of surrogate model bound lengths. """ - super(ExpIntVar, self).__init__(quantile_eps, *args, **opts) + super(ExpIntVar, self).__init__(model, prior, quantile_eps, **opts) self.name = 'exp_int_var' self.label_fn = 'Expected Loss' self._integration = integration