From 4a3335e54330e18b03bc90323e6822ee417f6f90 Mon Sep 17 00:00:00 2001 From: Tim Mensinger Date: Sun, 21 Jan 2024 16:19:19 +0100 Subject: [PATCH 1/6] Add method get_best to history class --- src/tranquilo/history.py | 19 +++++++++++++++++++ tests/test_history.py | 7 +++++++ 2 files changed, 26 insertions(+) diff --git a/src/tranquilo/history.py b/src/tranquilo/history.py index 7f81fefc..ddef245c 100644 --- a/src/tranquilo/history.py +++ b/src/tranquilo/history.py @@ -1,4 +1,5 @@ import numpy as np +import pandas as pd class History: @@ -173,6 +174,24 @@ def get_fvals(self, x_indices): ) return out + def get_best(self): + """Retrieve best fval and corresponding x and index. + + If there are multiple evaluations per x, the best fval is computed as the mean + of all evaluations per x. + + Returns: + tuple: (x, fval, index), the current minimizer x, the corresponding fval, + which is a mean if there are multiple evaluations per x, and the index. + + """ + average_fvals = { + x_index: np.mean(self.fvals[f_indices]) + for x_index, f_indices in self.index_mapper.items() + } + index = pd.Series(average_fvals).idxmin() + return self.get_xs(index), average_fvals[index], index + def get_n_evals(self, x_indices): fvals = self.get_fvals(x_indices) n_evals = {k: len(v) for k, v in fvals.items()} diff --git a/tests/test_history.py b/tests/test_history.py index efc0e298..13177a6a 100644 --- a/tests/test_history.py +++ b/tests/test_history.py @@ -228,3 +228,10 @@ def test_get_model_data_with_repeated_evaluations(noisy_history, average): else: aaae(got_xs, np.arange(6).reshape(2, 3).repeat([2, 3], axis=0)) aaae(got_fvecs, np.arange(25).reshape(5, 5)) + + +def test_get_best(noisy_history): + x, fval, index = noisy_history.get_best() + assert index == 0 + assert fval == 142.5 + aaae(x, np.array([0, 1, 2])) From 567433500533d655766a4fd06859eb6224842207 Mon Sep 17 00:00:00 2001 From: Tim Mensinger Date: Sun, 21 Jan 2024 16:36:53 +0100 Subject: [PATCH 2/6] Add method get_best to history class --- src/tranquilo/acceptance_decision.py | 103 +++++++++++++++++++++++++++ tests/test_acceptance_decision.py | 37 ++++++++++ 2 files changed, 140 insertions(+) diff --git a/src/tranquilo/acceptance_decision.py b/src/tranquilo/acceptance_decision.py index 534b3397..1904a5e9 100644 --- a/src/tranquilo/acceptance_decision.py +++ b/src/tranquilo/acceptance_decision.py @@ -25,6 +25,7 @@ def get_acceptance_decider( "naive_noisy": accept_naive_noisy, "noisy": accept_noisy, "classic_line_search": accept_classic_line_search, + "greedy": accept_greedy, } out = get_component( @@ -38,6 +39,38 @@ def get_acceptance_decider( return out +def accept_greedy( + subproblem_solution, + state, + history, + *, + wrapped_criterion, + min_improvement, +): + """Do a greedy acceptance step for a trustregion algorithm. + + Args: + subproblem_solution (SubproblemResult): Result of the subproblem solution. + state (State): Namedtuple containing the trustregion, criterion value of + previously accepted point, indices of model points, etc. + wrapped_criterion (callable): The criterion function. + min_improvement (float): Minimum improvement required to accept a point. + + Returns: + AcceptanceResult + + """ + out = _accept_greedy( + subproblem_solution=subproblem_solution, + state=state, + history=history, + wrapped_criterion=wrapped_criterion, + min_improvement=min_improvement, + n_evals=1, + ) + return out + + def _accept_classic( subproblem_solution, state, @@ -239,6 +272,76 @@ def accept_classic_line_search( return res +def _accept_greedy( + subproblem_solution, + state, + history, + *, + wrapped_criterion, + min_improvement, + n_evals, +): + """Do a simple greedy acceptance step for a trustregion algorithm. + + Args: + subproblem_solution (SubproblemResult): Result of the subproblem solution. + state (State): Namedtuple containing the trustregion, criterion value of + previously accepted point, indices of model points, etc. + wrapped_criterion (callable): The criterion function. + min_improvement (float): Minimum improvement required to accept a point. + + Returns: + AcceptanceResult + + """ + candidate_x = subproblem_solution.x + candidate_index = history.add_xs(candidate_x) + wrapped_criterion({candidate_index: n_evals}) + + candidate_fval = np.mean(history.get_fvals(candidate_index)) + actual_improvement = -(candidate_fval - state.fval) + + rho = calculate_rho( + actual_improvement=actual_improvement, + expected_improvement=subproblem_solution.expected_improvement, + ) + + best_x, best_fval, best_index = history.get_best() + + if best_fval < candidate_fval: + candidate_x = best_x + candidate_fval = best_fval + candidate_index = best_index + overall_improvement = -(candidate_fval - state.fval) + else: + overall_improvement = actual_improvement + + is_accepted = overall_improvement >= min_improvement + + if np.isfinite(candidate_fval): + res = _get_acceptance_result( + candidate_x=candidate_x, + candidate_fval=candidate_fval, + candidate_index=candidate_index, + rho=rho, + is_accepted=is_accepted, + old_state=state, + n_evals=n_evals, + ) + else: + res = _get_acceptance_result( + candidate_x=state.x, + candidate_fval=state.fval, + candidate_index=state.index, + rho=-np.inf, + is_accepted=False, + old_state=state, + n_evals=n_evals, + ) + + return res + + def _accept_simple( subproblem_solution, state, diff --git a/tests/test_acceptance_decision.py b/tests/test_acceptance_decision.py index 2be2bab9..71ac3ae6 100644 --- a/tests/test_acceptance_decision.py +++ b/tests/test_acceptance_decision.py @@ -4,6 +4,7 @@ import pytest from tranquilo.sample_points import get_sampler from tranquilo.acceptance_decision import ( + _accept_greedy, _accept_simple, _get_acceptance_result, calculate_rho, @@ -82,6 +83,42 @@ def wrapped_criterion(eval_info): assert_array_equal(res_got.candidate_x, 1.0 + np.arange(2)) +# ====================================================================================== +# Test accept_greedy +# ====================================================================================== + + +@pytest.mark.parametrize("state", states) +def test_accept_greedy( + state, + subproblem_solution, +): + history = History(functype="scalar") + + idxs = history.add_xs(np.arange(10).reshape(5, 2)) + + history.add_evals(idxs.repeat(2), np.arange(10)) + + def wrapped_criterion(eval_info): + indices = np.array(list(eval_info)).repeat(np.array(list(eval_info.values()))) + history.add_evals(indices, -indices) + + res_got = _accept_greedy( + subproblem_solution=subproblem_solution, + state=state, + history=history, + wrapped_criterion=wrapped_criterion, + min_improvement=0.0, + n_evals=2, + ) + + assert res_got.accepted + assert res_got.index == 5 + assert res_got.candidate_index == 5 + assert_array_equal(res_got.x, subproblem_solution.x) + assert_array_equal(res_got.candidate_x, 1.0 + np.arange(2)) + + # ====================================================================================== # Test _get_acceptance_result # ====================================================================================== From e3c1ac2d8dcdd06333ede612505e487ef4a1721f Mon Sep 17 00:00:00 2001 From: Tim Mensinger Date: Sun, 21 Jan 2024 17:09:32 +0100 Subject: [PATCH 3/6] Use history.get_fvals method --- src/tranquilo/history.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/src/tranquilo/history.py b/src/tranquilo/history.py index ddef245c..96841a5b 100644 --- a/src/tranquilo/history.py +++ b/src/tranquilo/history.py @@ -185,10 +185,8 @@ def get_best(self): which is a mean if there are multiple evaluations per x, and the index. """ - average_fvals = { - x_index: np.mean(self.fvals[f_indices]) - for x_index, f_indices in self.index_mapper.items() - } + fvals = self.get_fvals(np.arange(self.n_xs)) + average_fvals = {key: np.mean(val) for key, val in fvals.items()} index = pd.Series(average_fvals).idxmin() return self.get_xs(index), average_fvals[index], index From f63d502a267f7cd8d083cb7c6df7e83a08e8797e Mon Sep 17 00:00:00 2001 From: Tim Mensinger Date: Sun, 21 Jan 2024 17:18:07 +0100 Subject: [PATCH 4/6] Remove private function _accept_greedy --- src/tranquilo/acceptance_decision.py | 123 ++++++++++----------------- 1 file changed, 45 insertions(+), 78 deletions(-) diff --git a/src/tranquilo/acceptance_decision.py b/src/tranquilo/acceptance_decision.py index 1904a5e9..6d62e8ad 100644 --- a/src/tranquilo/acceptance_decision.py +++ b/src/tranquilo/acceptance_decision.py @@ -60,15 +60,52 @@ def accept_greedy( AcceptanceResult """ - out = _accept_greedy( - subproblem_solution=subproblem_solution, - state=state, - history=history, - wrapped_criterion=wrapped_criterion, - min_improvement=min_improvement, - n_evals=1, + candidate_x = subproblem_solution.x + candidate_index = history.add_xs(candidate_x) + wrapped_criterion({candidate_index: 1}) + + candidate_fval = np.mean(history.get_fvals(candidate_index)) + actual_improvement = -(candidate_fval - state.fval) + + rho = calculate_rho( + actual_improvement=actual_improvement, + expected_improvement=subproblem_solution.expected_improvement, ) - return out + + best_x, best_fval, best_index = history.get_best() + + if best_fval < candidate_fval: + candidate_x = best_x + candidate_fval = best_fval + candidate_index = best_index + overall_improvement = -(candidate_fval - state.fval) + else: + overall_improvement = actual_improvement + + is_accepted = overall_improvement >= min_improvement + + if np.isfinite(candidate_fval): + res = _get_acceptance_result( + candidate_x=candidate_x, + candidate_fval=candidate_fval, + candidate_index=candidate_index, + rho=rho, + is_accepted=is_accepted, + old_state=state, + n_evals=1, + ) + else: + res = _get_acceptance_result( + candidate_x=state.x, + candidate_fval=state.fval, + candidate_index=state.index, + rho=-np.inf, + is_accepted=False, + old_state=state, + n_evals=1, + ) + + return res def _accept_classic( @@ -272,76 +309,6 @@ def accept_classic_line_search( return res -def _accept_greedy( - subproblem_solution, - state, - history, - *, - wrapped_criterion, - min_improvement, - n_evals, -): - """Do a simple greedy acceptance step for a trustregion algorithm. - - Args: - subproblem_solution (SubproblemResult): Result of the subproblem solution. - state (State): Namedtuple containing the trustregion, criterion value of - previously accepted point, indices of model points, etc. - wrapped_criterion (callable): The criterion function. - min_improvement (float): Minimum improvement required to accept a point. - - Returns: - AcceptanceResult - - """ - candidate_x = subproblem_solution.x - candidate_index = history.add_xs(candidate_x) - wrapped_criterion({candidate_index: n_evals}) - - candidate_fval = np.mean(history.get_fvals(candidate_index)) - actual_improvement = -(candidate_fval - state.fval) - - rho = calculate_rho( - actual_improvement=actual_improvement, - expected_improvement=subproblem_solution.expected_improvement, - ) - - best_x, best_fval, best_index = history.get_best() - - if best_fval < candidate_fval: - candidate_x = best_x - candidate_fval = best_fval - candidate_index = best_index - overall_improvement = -(candidate_fval - state.fval) - else: - overall_improvement = actual_improvement - - is_accepted = overall_improvement >= min_improvement - - if np.isfinite(candidate_fval): - res = _get_acceptance_result( - candidate_x=candidate_x, - candidate_fval=candidate_fval, - candidate_index=candidate_index, - rho=rho, - is_accepted=is_accepted, - old_state=state, - n_evals=n_evals, - ) - else: - res = _get_acceptance_result( - candidate_x=state.x, - candidate_fval=state.fval, - candidate_index=state.index, - rho=-np.inf, - is_accepted=False, - old_state=state, - n_evals=n_evals, - ) - - return res - - def _accept_simple( subproblem_solution, state, From cdbd65e84ef4a5fe28b6a575dffb692fd14c3223 Mon Sep 17 00:00:00 2001 From: Tim Mensinger Date: Sun, 21 Jan 2024 17:44:06 +0100 Subject: [PATCH 5/6] Add test for accept_greedy --- tests/test_acceptance_decision.py | 41 ++++++++++++++++++++++--------- 1 file changed, 29 insertions(+), 12 deletions(-) diff --git a/tests/test_acceptance_decision.py b/tests/test_acceptance_decision.py index 71ac3ae6..8227fa21 100644 --- a/tests/test_acceptance_decision.py +++ b/tests/test_acceptance_decision.py @@ -1,10 +1,11 @@ from collections import namedtuple +from functools import partial import numpy as np import pytest from tranquilo.sample_points import get_sampler from tranquilo.acceptance_decision import ( - _accept_greedy, + accept_greedy, _accept_simple, _get_acceptance_result, calculate_rho, @@ -93,30 +94,46 @@ def test_accept_greedy( state, subproblem_solution, ): + """Test accept greedy. + + Tests that the best point is chosen in the acceptance step, even though it is added + to the history before the acceptance step. + + """ history = History(functype="scalar") - idxs = history.add_xs(np.arange(10).reshape(5, 2)) + def criterion(x): + return np.sum(x**2) - history.add_evals(idxs.repeat(2), np.arange(10)) + def _wrapped_criterion(eval_info, history): + for x_index, _ in eval_info.items(): + xs = history.get_xs(x_index) + crit_value = criterion(xs) + history.add_evals(np.array([x_index]), crit_value) - def wrapped_criterion(eval_info): - indices = np.array(list(eval_info)).repeat(np.array(list(eval_info.values()))) - history.add_evals(indices, -indices) + wrapped_criterion = partial(_wrapped_criterion, history=history) - res_got = _accept_greedy( + # Add existing xs to history and evaluate wrapped criterion + existing_xs = np.zeros((1, 2)) + existing_xs_indices = history.add_xs(existing_xs) + + eval_info = {x_index: 1 for x_index in existing_xs_indices} + wrapped_criterion(eval_info) + + res_got = accept_greedy( subproblem_solution=subproblem_solution, state=state, history=history, wrapped_criterion=wrapped_criterion, min_improvement=0.0, - n_evals=2, ) assert res_got.accepted - assert res_got.index == 5 - assert res_got.candidate_index == 5 - assert_array_equal(res_got.x, subproblem_solution.x) - assert_array_equal(res_got.candidate_x, 1.0 + np.arange(2)) + assert res_got.index == 0 + assert res_got.candidate_index == 0 + assert res_got.fval == 0.0 + assert_array_equal(res_got.x, np.zeros(2)) + assert_array_equal(res_got.candidate_x, np.zeros(2)) # ====================================================================================== From 9c2d376ca40c41d5fca0b3f69fa486ed9ff6e5a5 Mon Sep 17 00:00:00 2001 From: Tim Mensinger Date: Sun, 21 Jan 2024 18:30:40 +0100 Subject: [PATCH 6/6] Debugging --- src/tranquilo/acceptance_decision.py | 19 +++++++++++++------ src/tranquilo/history.py | 2 +- tests/test_history.py | 1 + 3 files changed, 15 insertions(+), 7 deletions(-) diff --git a/src/tranquilo/acceptance_decision.py b/src/tranquilo/acceptance_decision.py index 6d62e8ad..2384958a 100644 --- a/src/tranquilo/acceptance_decision.py +++ b/src/tranquilo/acceptance_decision.py @@ -65,22 +65,29 @@ def accept_greedy( wrapped_criterion({candidate_index: 1}) candidate_fval = np.mean(history.get_fvals(candidate_index)) - actual_improvement = -(candidate_fval - state.fval) + candidate_improvement = -(candidate_fval - state.fval) rho = calculate_rho( - actual_improvement=actual_improvement, + actual_improvement=candidate_improvement, expected_improvement=subproblem_solution.expected_improvement, ) best_x, best_fval, best_index = history.get_best() - - if best_fval < candidate_fval: + + assert np.isfinite(best_fval) + assert isinstance(best_x, np.ndarray) + assert isinstance(best_index, int) + assert isinstance(best_fval, float) + assert best_x.ndim == 1 + assert np.mean(history.get_fvals(best_index)) == best_fval + + if best_fval < candidate_fval and best_fval < state.fval: candidate_x = best_x candidate_fval = best_fval candidate_index = best_index - overall_improvement = -(candidate_fval - state.fval) + overall_improvement = -(best_fval - state.fval) else: - overall_improvement = actual_improvement + overall_improvement = candidate_improvement is_accepted = overall_improvement >= min_improvement diff --git a/src/tranquilo/history.py b/src/tranquilo/history.py index 96841a5b..39f63406 100644 --- a/src/tranquilo/history.py +++ b/src/tranquilo/history.py @@ -187,7 +187,7 @@ def get_best(self): """ fvals = self.get_fvals(np.arange(self.n_xs)) average_fvals = {key: np.mean(val) for key, val in fvals.items()} - index = pd.Series(average_fvals).idxmin() + index = int(pd.Series(average_fvals).idxmin()) return self.get_xs(index), average_fvals[index], index def get_n_evals(self, x_indices): diff --git a/tests/test_history.py b/tests/test_history.py index 13177a6a..ab9a5484 100644 --- a/tests/test_history.py +++ b/tests/test_history.py @@ -233,5 +233,6 @@ def test_get_model_data_with_repeated_evaluations(noisy_history, average): def test_get_best(noisy_history): x, fval, index = noisy_history.get_best() assert index == 0 + assert isinstance(index, int) assert fval == 142.5 aaae(x, np.array([0, 1, 2]))