diff --git a/changelog.md b/changelog.md index 6fc380eb7..11cca5502 100644 --- a/changelog.md +++ b/changelog.md @@ -23,6 +23,7 @@ * Added improved printing of calibration items under `channel.calibration` providing a more convenient overview of the items associated with a `Slice`. * Added improved printing of calibrations performed with `Pylake`. * Improved error message that includes the name of the model when trying to access a model that was not added in an [`FdFit`](https://lumicks-pylake.readthedocs.io/en/latest/_api/lumicks.pylake.FdFit.html) using angular brackets. +* Allow customizing the minimum step size during step size determination for [`DwelltimeModel.profile_likelihood()`](https://lumicks-pylake.readthedocs.io/en/latest/_api/lumicks.pylake.DwelltimeModel.html#lumicks.pylake.DwelltimeModel.profile_likelihood) and choose a more sensible default. Also ensured that the warning only gets emitted at most once per direction. ## v1.5.3 | 2024-10-29 diff --git a/lumicks/pylake/fitting/fit.py b/lumicks/pylake/fitting/fit.py index c85d465f2..088240447 100644 --- a/lumicks/pylake/fitting/fit.py +++ b/lumicks/pylake/fitting/fit.py @@ -301,6 +301,14 @@ def profile_likelihood( for the profiled parameter, then fixes that parameter and re-optimizes all the other parameters [4]_ [5]_. + For each profile likelihood point, the profiled parameter is changed such that the fit + error goes up by an amount between `min_chi2_step` and `max_chi2_step`. After step size + determination, the step is applied and the other parameters are re-optimized. + + To prevent step sizes from becoming too small, a minimal step size can be set. When the + step size goes below the minimum step size a warning is issued. When this happens, check + to verify that the profile likelihood curve is smooth. + Parameters ---------- parameter_name: str diff --git a/lumicks/pylake/fitting/profile_likelihood.py b/lumicks/pylake/fitting/profile_likelihood.py index db973e1e3..dff275e65 100644 --- a/lumicks/pylake/fitting/profile_likelihood.py +++ b/lumicks/pylake/fitting/profile_likelihood.py @@ -235,10 +235,7 @@ def do_step( last_step = StepCode.identify_step(new_step_size, current_step_size) - if current_step_size == step_config.min_abs_step: - warn("Warning: Step size set to minimum step size.", RuntimeWarning) - - return current_step_size, p_trial + return current_step_size, p_trial, current_step_size == step_config.min_abs_step def scan_dir_optimisation( @@ -276,11 +273,14 @@ def scan_dir_optimisation( p_next = parameter_vector chi2_list = [] parameter_vectors = [] + minimum_step_warning = False for step in range(num_steps): - current_step_size, p_next = scan_config.step_function( + current_step_size, p_next, minimum_step_hit = scan_config.step_function( chi2_last, p_next, current_step_size, step_sign ) + minimum_step_warning = minimum_step_warning or minimum_step_hit + try: p_next = fit_function( p_next, scan_config.lower_bounds, scan_config.upper_bounds, scan_config.fitted @@ -314,8 +314,14 @@ def scan_dir_optimisation( if chi2_last > scan_config.termination_level: break - return np.array(chi2_list), ( - np.vstack(parameter_vectors) if parameter_vectors else np.empty((0, parameter_vector.size)) + return ( + np.array(chi2_list), + ( + np.vstack(parameter_vectors) + if parameter_vectors + else np.empty((0, parameter_vector.size)) + ), + minimum_step_warning, ) @@ -547,10 +553,19 @@ def _extend_profile(self, chi2_function, fit_function, parameters, num_steps, fo chi2_list = self._chi2[field] if field in self._chi2 else [self.profile_info.minimum_chi2] chi2_last = chi2_list[-1] - chi2_new, parameters_new = scan_direction( + chi2_new, parameters_new, minimum_step_warning = scan_direction( chi2_last, initial_parameters, 1 if forward else -1, num_steps, verbose ) + if minimum_step_warning: + warn( + "Step size was set to minimum step size during step size determination " + f"for parameter {self.parameter_name}. Check whether the profile likelihood " + "result is smooth. If not, consider reducing the minimum step size; " + "otherwise, ignore this warning.", + RuntimeWarning, + ) + self._parameters[field] = np.vstack((parameter_vectors, parameters_new)) self._chi2[field] = np.hstack((chi2_list, chi2_new)) diff --git a/lumicks/pylake/fitting/tests/test_profiles.py b/lumicks/pylake/fitting/tests/test_profiles.py index 21982c006..c751b40f0 100644 --- a/lumicks/pylake/fitting/tests/test_profiles.py +++ b/lumicks/pylake/fitting/tests/test_profiles.py @@ -216,7 +216,7 @@ def linear_with_exception(x, a=0.01, b=1): fit["model/a"].value = 1 fit.fit() with warnings.catch_warnings(): - warnings.filterwarnings("ignore", message="Warning: Step size set to minimum step size") + warnings.filterwarnings("ignore", message="Step size was set to minimum step size") warnings.filterwarnings("ignore", message="Optimization error encountered") profile = fit.profile_likelihood("model/a") diff --git a/lumicks/pylake/fitting/tests/test_stepping.py b/lumicks/pylake/fitting/tests/test_stepping.py index 3e2a382cd..6f8ff53f1 100644 --- a/lumicks/pylake/fitting/tests/test_stepping.py +++ b/lumicks/pylake/fitting/tests/test_stepping.py @@ -82,7 +82,7 @@ def __call__(self, x): ], ) def test_stepper(chi2_function, start_pos, step_sign, target_step_size, target_p_trial): - step_size, p_trial = do_step( + step_size, p_trial, will_warn = do_step( chi2_function=chi2_function, step_direction_function=lambda: np.array([1, 1]), chi2_last=5, @@ -94,6 +94,7 @@ def test_stepper(chi2_function, start_pos, step_sign, target_step_size, target_p assert step_size == target_step_size np.testing.assert_allclose(p_trial, target_p_trial) + assert not will_warn def test_minimum_step(): @@ -102,17 +103,15 @@ def test_minimum_step(): # step size to take is 10 (see max_chi2_step_size in the cfg), this will result in the step size # shrinking to the minimum step size. In a non-pathological case, a smaller step size would # result in a smaller chi2 increase. - with pytest.warns(RuntimeWarning, match="Warning: Step size set to minimum step size."): - step_size, p_trial = do_step( - chi2_function=lambda x: 25, - step_direction_function=lambda: np.array([1, 1]), - chi2_last=5, - parameter_vector=np.array([2, 2]), - current_step_size=1, - step_sign=1, - step_config=cfg, - ) - assert step_size == cfg.min_abs_step - np.testing.assert_allclose( - p_trial, np.array([2.0 + cfg.min_abs_step, 2.0 + cfg.min_abs_step]) - ) + step_size, p_trial, will_warn = do_step( + chi2_function=lambda x: 25, + step_direction_function=lambda: np.array([1, 1]), + chi2_last=5, + parameter_vector=np.array([2, 2]), + current_step_size=1, + step_sign=1, + step_config=cfg, + ) + assert step_size == cfg.min_abs_step + np.testing.assert_allclose(p_trial, np.array([2.0 + cfg.min_abs_step, 2.0 + cfg.min_abs_step])) + assert will_warn diff --git a/lumicks/pylake/population/dwelltime.py b/lumicks/pylake/population/dwelltime.py index 9e6683ee8..7ab742b65 100644 --- a/lumicks/pylake/population/dwelltime.py +++ b/lumicks/pylake/population/dwelltime.py @@ -44,6 +44,8 @@ def _from_dwelltime_model( *, alpha, num_steps, + min_step, + max_step, min_chi2_step, max_chi2_step, verbose, @@ -60,6 +62,11 @@ def _from_dwelltime_model( Significance level. Confidence intervals are calculated as 100*(1-alpha)%. num_steps: integer Number of steps to take. + min_step : float + Minimum step size used in step size determination. Any step size lower will be clamped + to this value and a warning will be issued. + max_step : float + Maximum step size used in step size determination. Any step size larger will be clamped. min_chi2_step: float Minimal desired step in terms of chi squared change prior to re-optimization. When the step results in a fit change smaller than this threshold, the step-size will be @@ -109,6 +116,8 @@ def calculate_profile(param): num_dof=1, min_chi2_step=min_chi2_step, max_chi2_step=max_chi2_step, + min_step=min_step, + max_step=max_step, confidence_level=1.0 - alpha, bound_tolerance=1e-8, # Needed because constraint is not always exactly fulfilled ) @@ -689,6 +698,8 @@ def profile_likelihood( *, alpha=0.05, num_steps=150, + min_step=1e-7, + max_step=1.0, min_chi2_step=0.01, max_chi2_step=0.1, verbose=False, @@ -699,6 +710,14 @@ def profile_likelihood( confidence intervals. It iteratively performs a step for the profiled parameter, then fixes that parameter and re-optimizes all the other parameters [2]_ [3]_. + For each profile likelihood point, the profiled parameter is changed such that the fit + error goes up by an amount between `min_chi2_step` and `max_chi2_step`. After step size + determination, the step is applied and the other parameters are re-optimized. + + To prevent step sizes from becoming too small, a minimal step size can be set. When the + step size goes below the minimum step size a warning is issued. When this happens, check + to verify that the profile likelihood curve is smooth. + Parameters ---------- alpha : float @@ -706,6 +725,12 @@ def profile_likelihood( default: 0.05 num_steps: int Number of steps to take, default: 150. + min_step : float + Minimum step size used during step size determination. Any step size lower will be + clamped to this value and a warning will be issued. Default: 1e-7. + max_step : float + Maximum step size used during step size determination. Any step size larger will be + clamped to this value. Default: 1.0 min_chi2_step: float Minimal desired step in terms of chi squared change prior to re-optimization. When the step results in a fit change smaller than this threshold, the step-size will be @@ -731,6 +756,8 @@ def profile_likelihood( self, alpha=alpha, num_steps=num_steps, + min_step=min_step, + max_step=max_step, min_chi2_step=min_chi2_step, max_chi2_step=max_chi2_step, verbose=verbose,