profile_likelihood: silence unnecessary warnings

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
 

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

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))
 

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")
 

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

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,13 +710,27 @@ 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
             Significance level. Confidence intervals are calculated as 100 * (1 - alpha)%,
             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,