Add new lm minimiser #208

RichardWaiteSTFC · 2024-11-01T16:17:17Z

New features:

New (simpler) LM minimiser lm4 that support minimising scalars correctly and with optional bounds on parameters
Add support to minimise parameters using function handle returning flattened array of residuals
Add support to pass residual array rather than scalar in sw_fitpowder (powder fitting class)
Add option vary to simplex to specify which parameters to vary in the fit (provided this option in lm4 as well)

Changes to behaviour:

cost_function_wrapper (and minimisers which use it - e.g. simplex and lm4) will reset parameters out of bounds provided using same method as lsqnonlin (previously would reset to be nearest boundary, but this didn't work well for the LM method).
Also cost_function_wrapper will now throw a warning when parameters have been reset inside the bounds.

2 tests omitted for lm4 - seems to be issue if parameter exactly at bound

codecov-commenter · 2024-11-01T16:23:58Z

⚠️ Please install the to ensure uploads and comments are reliably processed by Codecov.

Codecov Report

Attention: Patch coverage is 92.35294% with 13 lines in your changes missing coverage. Please review.

Project coverage is 43.29%. Comparing base (73f604a) to head (1fc202a).

Files with missing lines	Patch %	Lines
swfiles/+ndbase/lm4.m	93.51%	7 Missing ⚠️
swfiles/sw_fitpowder.m	73.33%	4 Missing ⚠️
swfiles/+ndbase/cost_function_wrapper.m	94.11%	2 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##           master     #208      +/-   ##
==========================================
+ Coverage   42.86%   43.29%   +0.42%     
==========================================
  Files         242      243       +1     
  Lines       16240    16383     +143     
==========================================
+ Hits         6962     7093     +131     
- Misses       9278     9290      +12

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

github-actions · 2024-11-01T16:26:15Z

Test Results

4 files ± 0 124 suites ±0 4m 47s ⏱️ -2s
831 tests +40 813 ✅ +40 18 💤 ±0 0 ❌ ±0
2 324 runs +84 2 288 ✅ +84 36 💤 ±0 0 ❌ ±0

Results for commit 1fc202a. ± Comparison against base commit 73f604a.

This pull request removes 2 and adds 42 tests. Note that renamed tests count towards both.

sw_tests.unit_tests.unittest_ndbase_cost_function_wrapper ‑ test_init_with_fcost_both_bounds_with_fixed_param_using_ifix
sw_tests.unit_tests.unittest_sw_fitpowder ‑ test_fit_background

sw_tests.unit_tests.unittest_ndbase_cost_function_wrapper ‑ test_init_with_fcost_all_params_fixed
sw_tests.unit_tests.unittest_ndbase_cost_function_wrapper ‑ test_init_with_fcost_both_bounds_fixed_invalid_param_using_ifix
sw_tests.unit_tests.unittest_ndbase_cost_function_wrapper ‑ test_init_with_fcost_both_bounds_fixed_param_using_ifix
sw_tests.unit_tests.unittest_ndbase_cost_function_wrapper ‑ test_init_with_resid_handle
sw_tests.unit_tests.unittest_ndbase_optimisers ‑ test_optimise_data_struct(optimiser=@ndbase.lm4,poly_func=char_@_x__p__polyval_p__x_)
sw_tests.unit_tests.unittest_ndbase_optimisers ‑ test_optimise_data_struct(optimiser=@ndbase.lm4,poly_func=function_handle)
sw_tests.unit_tests.unittest_ndbase_optimisers ‑ test_optimise_data_struct(optimiser=value2,poly_func=value1)
sw_tests.unit_tests.unittest_ndbase_optimisers ‑ test_optimise_data_struct(optimiser=value2,poly_func=value2)
sw_tests.unit_tests.unittest_ndbase_optimisers ‑ test_optimise_residual_array_lm(optimiser=@ndbase.lm4)
sw_tests.unit_tests.unittest_ndbase_optimisers ‑ test_optimise_residual_array_lm(optimiser=@ndbase.simplex)
…

♻️ This comment has been updated with latest results.

In ndbase.lm4 and simplex optimisers. Also added tests This allows support for ND fitting problems.

And add test

mducle

Thanks for this. It should probably be two different PR, one for the residual return function and one for lm4 but I think it's fine since the residual stuff is small; but I'll split my comments on the two parts:

The residual stuff is good. I don't have a problem with using the resid_handle flag. The only other way I can think of to check it is to evaluate the model function and check if it returns a scalar or not but this could be costly just on construction to test what kind of input it is. For lm4 where we evaluate to get the starting cost value then this might be a way but I'm not sure if the simplex minimiser does this and I don't know how to implement it since the wrapper is separate from the mimisers...
The docstring for lm4 needs a lot of updating. I guess it was taken from another function?
Some of the unit tests (especially for bounded problems) are disabled for lm4 but work for simplex - is this a weakness of lm4? Should we note in the docstring that bounded fitting with lm4 is a bit flaky and maybe suggest users to specify larger bounds?

mducle · 2024-11-12T14:29:27Z

+sw_tests/+unit_tests/unittest_ndbase_cost_function_wrapper.m

+        end
+
+        function test_init_with_fcost_all_params_fixed(testCase)
+            % note second param outside bounds


I think this comment is wrong? It doesn't look like you've set any bounds here, just fixed the parameter values?
Looks like this was copied from test_init_with_fcost_no_bounds_with_fixed_param_using_ifix which also has the same (erroneous?) comment in line 96.

mducle · 2024-11-12T14:32:40Z

+sw_tests/+unit_tests/unittest_ndbase_optimisers.m

+            testCase.verify_val(cost_val, 0, 'abs_tol', 1e-6);
+        end
+
+        function test_optimise_residual_array_lm(testCase, optimiser)


This test should also run for ndbase.simplex too, right? Should the _lm suffix be removed?

+sw_tests/+unit_tests/unittest_ndbase_optimisers.m

mducle · 2024-11-12T16:11:46Z

+sw_tests/+unit_tests/unittest_ndbase_optimisers.m

            testCase.verify_val(pars_fit, testCase.rosenbrock_minimum, 'abs_tol', 1e-3);
            testCase.verify_val(cost_val, 0, 'abs_tol', 1e-6);
        end

-        function test_optimise_rosen_both_bounds_minimum_not_accessible(testCase, optimiser)
+        function test_optimise_rosen_both_bounds_minimum_not_accessible(testCase)


So in this test the bounds given are all outside the optimum parameter values of (1, 1) - in this case would we expect all minimisers to try to find the optimum value within the bounds? If so, do you know why lm4 isn't able to do this? Is it something to do with how we've implemented the bounds? (Is it because the gradient becomes too steep near the boundaries?)

(Just running it manually with these bounds and the Rosenbrock function I get optimised parameters (-0.5, -0.5) which is just the lower bound... likewise with the next test (fixed_minimum_not_accessible) where lm4 just returns the lower bound.)

Just checked and lsqnonlin does find the lowest cost val here at (0,0)
I think the problem could be solved by moving slightly away from the boundary edge when we reset invalid parameters outside bounds?

Yeh if you change to have the initial guess 0.1 inside the lb like so
'lb', [-1.1, -1.1]
This test passes for lm4 minimiser

I just had a play and it still seems flaky. For some bounds it works and for others it doesn't... I wonder if we have to transform the diff_step too?

I think the diff step has to be taken in the transformed coordinate so as not to go over the bounds - the diff step is relative - i.e. scaled by the magnitude of the scaled/free parameter rather than the bound one which I also think it correct. We can discuss tomorrow!

So with implementing the box start parameters in this gist I could get all the tests to pass with for both simplex and lm4, even with the original bounds in both_bounds_minimum_accessible and with a tighter tolerance for upper_bound_minimum_not_accessible (1e-4 compared to 1e-6 for simplex), but I had to change the default nu_up from 10 to 5 in order to get those two tests to fit with the original bounds and tighter tolerance... not sure if this is ok?

Thanks @mducle - I'll look at this today - in theory I don't have a problem with changing nu_up=5 as long as it still passes the NIST data sets.
Here is the advice given on wikipedia...

An effective strategy for the control of the damping parameter, called delayed gratification, consists of increasing the parameter by a small amount for each uphill step, and decreasing by a large amount for each downhill step. The idea behind this strategy is to avoid moving downhill too fast in the beginning of optimization, therefore restricting the steps available in future iterations and therefore slowing down convergence.[7] An increase by a factor of 2 and a decrease by a factor of 3 has been shown to be effective in most cases, while for large problems more extreme values can work better, with an increase by a factor of 1.5 and a decrease by a factor of 5.[8]

Which I think (if I'm interpreting it correctly) would correspond to nu_up=2, nu_dn=0.33 (similar value to existing default value for nu_dn). I will try nu_up = 10 (existing default), 5, 2 and see how they compare in terms of the number of iterations as well as accuracy. FYI I got my value of nu_up from numerical recipes (or something based on numerical recipes...).

So I quickly tried this on the Gauss3 NIST dataset (so-called average difficulty with 8 parameters) and get these results (note a negative dcost indicates the chi-squared value we found is larger than the NIST provided value - but as you can see only by a very small amount, this could probably be reduced by changing the default convergence criteria!)

% minimising residual array % nu_up nu_dn iterations nhess_evals dcost(%) % 10 0.3 11 10 -1.1e-9 % 5 0.3 9 8 -1.1e-9 % 2 0.3 11 8 -1.1e-9 % 1.5 0.2 15 8 -1.1e-9 % minimising scalar % nu_up nu_dn iterations nhess_evals dcost(%) % 10 0.3 12 10 -1.2e-7 % 5 0.3 17 8 -6.5e-8 % 2 0.3 16 10 -1.2e-7 % 1.5 0.2 21 9 -1.2e-7

I will try on a more difficult dataset...

Tried fitting the Eckerle4 NIST dataset (high difficulty - 3 parameters) - note lm3 can't fit this
Here are the results

% minimising residual array % nu_up nu_dn iterations nhess_evals dcost(%) % 10 0.3 121 81 -5.7e-5 % 5 0.3 64 39 -1.6e-5 % 2 0.3 136 54 -2.3e-7 % 1.5 0.2 166 38 -3.1e-5 % minimising scalar % nu_up nu_dn iterations nhess_evals dcost(%) % 10 0.3 55 38 -0.7 % 5 0.3 47 29 -0.7 % 2 0.3 76 17 -2e4 % very bad fit! % 1.5 0.2 74 19 -0.7

Seems like nu_up=5, nu_dn=0.3 is a bit of a sweet spot! Albeit having not given it the rigorous fitbenchmarking treatment...

swfiles/+ndbase/lm4.m

mducle · 2024-11-12T16:46:50Z

swfiles/+ndbase/lm4.m

+stat.func = cost_func_wrap.cost_func;
+stat.param = param;
+stat.param.Np = cost_func_wrap.get_num_free_parameters();
+stat.msg        = message;


In the docstring this is message - it might be good to change it here to match the docstring...

swfiles/+ndbase/lm4.m

RichardWaiteSTFC · 2024-11-12T17:38:39Z

Thanks for the review @mducle - sorry for the long PR, it probably should have been 2 as you say!

I've been on unscheduled leave most of the last two days so as didn't get a chance to give final polish - in particular I didn't get a chance to change docstring (taken from lm3) as I had planned as discussed at stand-up!

Some of the unit tests (especially for bounded problems) are disabled for lm4 but work for simplex - is this a weakness of lm4? Should we note in the docstring that bounded fitting with lm4 is a bit flaky and maybe suggest users to specify larger bounds?

Yes there are 2 tests disabled for the lm4 minimiser - it fails to find the minimum within the constraints (note for these tests the true minima of the cost function is not accessible)! The problem is indeed due to the bounds, if the initial guess is too close to one of the bounds (which happens automatically if the initial parameters are outside the bounds) the lm4 minimiser seems to 'converge' very prematurely (sometimes after only 1 iteration) for one of the change in cost function or parameter step criteria the algorithm. These can be fiddled with (along with the step size used) but it is still flaky - I think this problem is perhaps inherent to any method using finite differences - it would be interesting to compare to lmfit.

I will add a warning when the parameter is reset at the boundary (as lm4 users should want to avoid this).

I also wondered whether a different parameter transformation may be better? But it is worth noting that it is enabled for the majority of tests with bound parameters!

RichardWaiteSTFC · 2024-11-15T19:38:50Z

Since the last review:

Have changed default nu_dn=5 in lm4
cost_function_wrapper (and minimisers which use it - i.e. simplex and lm4) will reset parameters out of bounds provided using same method as lsqnonlin (previously would reset to be nearest boundary, but this didn't work well for the LM method as discussed). A warning is now printed when this has been done.
- Had to use a slightly different method to your gist sorry, was some tricky behaviour with setting parameters correctly if fixed (as these are cached in wrapper class)
Enabled all the unit tests for the lm4 optimiser (and tightened some tolerances).
Add option vary to simplex and lm4 (lm,lm2 and lm3 all have this option) to specify which parameters to vary in the fit.
Updated doc-strings (hopefully giving user some idea of what the LM parameters do...at least as I understand them!)

mducle

Looks good! Just a couple of typos in the docstring.

swfiles/+ndbase/lm4.m

RichardWaiteSTFC added 11 commits September 27, 2024 18:30

Optional inital cost_val and return jacobuan in estimate_hessian

959c38f

Add function to evaluate weighted residuals in cost func wrapper

65d0cfc

Add LM minimiser (currently only supports scalar functions)

2474a3a

Support least_squares residuals in lm4

d8b14e9

Improve inversion for ill-conditioned matrices

7bf3529

Add missing factors of 2

b94c4c7

Fix bug in finite difference not resetting parameter vector

d1e8605

Fix bug with diff step length if parameters fixed

e284a58

Fix bug where absolute and fracatinal step confused if 1 parameter

40c9095

Add lm4 to optimiser to unit test

b00c3d3

2 tests omitted for lm4 - seems to be issue if parameter exactly at bound

Merge branch 'master' into add_new_LM_minimiser

dd978e2

RichardWaiteSTFC added 3 commits November 4, 2024 15:04

Add support for function handles returning residal arrays

44a6290

In ndbase.lm4 and simplex optimisers. Also added tests This allows support for ND fitting problems.

Support resid_handle option in sw_fitpowder

3dc3f2f

Support resid_handle in bg fitting of pwader class

520e53c

And add test

mducle requested changes Nov 12, 2024

View reviewed changes

RichardWaiteSTFC added 4 commits November 15, 2024 15:11

Correct doc-strings for simplex and lm4

e579fe2

Change nu_dn default to 5 after benchmarking with NIST datasets

5c3f460

Use lsqnonlin method to reset invalid parameters within bounds

8cb0a1a

Add vary arg to simplex and lm4 (and tests)

5b72585

mducle requested changes Nov 15, 2024

View reviewed changes

swfiles/+ndbase/lm4.m Outdated Show resolved Hide resolved

swfiles/+ndbase/lm4.m Outdated Show resolved Hide resolved

Fix typos in doc string

1fc202a

mducle approved these changes Nov 18, 2024

View reviewed changes

RichardWaiteSTFC merged commit 7e87dee into master Nov 18, 2024
12 checks passed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Add new lm minimiser #208

Add new lm minimiser #208

RichardWaiteSTFC commented Nov 1, 2024 •

edited

Loading

codecov-commenter commented Nov 1, 2024 •

edited

Loading

github-actions bot commented Nov 1, 2024 •

edited

Loading

mducle left a comment

mducle Nov 12, 2024

mducle Nov 12, 2024

mducle Nov 12, 2024

RichardWaiteSTFC Nov 12, 2024

RichardWaiteSTFC Nov 12, 2024 •

edited

Loading

mducle Nov 12, 2024

RichardWaiteSTFC Nov 12, 2024

mducle Nov 14, 2024

RichardWaiteSTFC Nov 15, 2024

RichardWaiteSTFC Nov 15, 2024 •

edited

Loading

RichardWaiteSTFC Nov 15, 2024 •

edited

Loading

RichardWaiteSTFC Nov 15, 2024

mducle Nov 12, 2024

RichardWaiteSTFC commented Nov 12, 2024 •

edited

Loading

RichardWaiteSTFC commented Nov 15, 2024 •

edited

Loading

mducle left a comment

Add new lm minimiser #208

Add new lm minimiser #208

Conversation

RichardWaiteSTFC commented Nov 1, 2024 • edited Loading

codecov-commenter commented Nov 1, 2024 • edited Loading

Codecov Report

github-actions bot commented Nov 1, 2024 • edited Loading

Test Results

mducle left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

RichardWaiteSTFC Nov 12, 2024 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

RichardWaiteSTFC Nov 15, 2024 • edited Loading

Choose a reason for hiding this comment

RichardWaiteSTFC Nov 15, 2024 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

RichardWaiteSTFC commented Nov 12, 2024 • edited Loading

RichardWaiteSTFC commented Nov 15, 2024 • edited Loading

mducle left a comment

Choose a reason for hiding this comment

RichardWaiteSTFC commented Nov 1, 2024 •

edited

Loading

codecov-commenter commented Nov 1, 2024 •

edited

Loading

github-actions bot commented Nov 1, 2024 •

edited

Loading

RichardWaiteSTFC Nov 12, 2024 •

edited

Loading

RichardWaiteSTFC Nov 15, 2024 •

edited

Loading

RichardWaiteSTFC Nov 15, 2024 •

edited

Loading

RichardWaiteSTFC commented Nov 12, 2024 •

edited

Loading

RichardWaiteSTFC commented Nov 15, 2024 •

edited

Loading