Incorrectly reported Func-count in lsqcurvefit. Also, MaxIterations not obeyed

In R2020b and R20201a, I have done the following test of lsqcurvefit. The test function below test(M) uses lsqcurvefit to solve a simple linear equation A*x=y under the assumption that the first M variables are known (enforced by setting lb(1:M)=ub(1:M)) and counts the number of function evaluations. I find, strangely, that for any M>1, an additional N=100 function evaluations are done beyond what is reported by lsqcurvefit (and beyond what should be necessary). Also, even though I have specified MaxIterations=5, six iterations are sometimes done. Does anyone know why this occurs?

 TRUECOUNTS= test(0), 

                                         Norm of      First-order 
 Iteration  Func-count     f(x)          step          optimality
     0        101     6.44539e+08                      1.45e+05
     1        202     6.19343e+08             10       1.42e+05      
     2        303      5.7046e+08             20       1.36e+05      
     3        404     4.78727e+08             40       1.25e+05      
     4        505       3.194e+08             80       1.02e+05      
     5        606     9.73003e+07            160       5.61e+04      
     6        707      5.4583e-07        346.511       0.000622      

Solver stopped prematurely.

lsqcurvefit stopped because it exceeded the iteration limit,
options.MaxIterations = 5.000000e+00.

 

TRUECOUNTS = 707

 
 TRUECOUNTS=test(1), 

                                         Norm of      First-order 
 Iteration  Func-count     f(x)          step          optimality
     0        100     6.54782e+08                      1.45e+05
     1        200     6.29379e+08             10       1.42e+05      
     2        300     5.80081e+08             20       1.36e+05      
     3        400     4.87524e+08             40       1.25e+05      
     4        500     3.26557e+08             80       1.02e+05      
     5        600     1.01211e+08            160       5.69e+04      
     6        700     5.41302e-07        343.045       0.000567      

Solver stopped prematurely.

lsqcurvefit stopped because it exceeded the iteration limit,
options.MaxIterations = 5.000000e+00.

 

TRUECOUNTS = 800

 
 TRUECOUNTS=test(99),

                                         Norm of      First-order 
 Iteration  Func-count     f(x)          step          optimality
     0          2          324460                      3.24e+03
     1          4          262813             10       2.92e+03      
     2          6          158986             20       2.27e+03      
     3          8         29201.4             40            973      
     4         10     1.06417e-10             30       5.88e-05      
     5         12               0    1.81103e-06              0      

Local minimum found.

Optimization completed because the size of the gradient is less than
the value of the optimality tolerance.

 

TRUECOUNTS = 112

 function TRUECOUNTS=test(M)
   
   if nargin<1, M=0; end

    N=100;
    xtrue=(1:N).';  %true parameter values
    Adata=rand(N);

    
    TRUECOUNTS=0;
    ydata=F(xtrue,Adata);

    lb=-inf(N,1); ub=+inf(N,1);
    
   
    known=1:M;               %The first M variables are known      
    lb(known)=xtrue(known); 
    ub(known)=xtrue(known);

    opts=optimoptions(@lsqcurvefit,  'Display','iter','MaxIterations',5);

    
    TRUECOUNTS=0;
    lsqcurvefit(@F,zeros(N,1),Adata,ydata,lb,ub,opts);
    disp ' ' 

    function y=F(x,A)

        y=A*x;
        
        TRUECOUNTS=TRUECOUNTS+1;

        
    end

end

Answer 1

[Alan Weiss]

Thank you for reporting this odd behavior, Matt. I investigated and found that an internal check for stopping in the nonlinear least-squares code is

if iterationNumber > maxiter

Therefore, when you set 5 as the iteration limit, the check is triggered at the sixth iteration, not the fifth.

In contrast, several other solvers have the more natural

if iterationNumber >= maxiter

This stops the solver at iteration 5. So this is at least an inconsistency, if not a bug. I will now report this in our internal development system.

By the way, the reason that your test sometimes found only 5 iterations, not 6, is that the problem was solved, not that the iteration limit was hit. Your test(99) is that case.

As for why an extra 100 function evaluations are done, your code line

ydata=F(xtrue,Adata);

causes those evaluations, not the lsqcurvefit code.

Thanks again for the report,

Alan Weiss

MATLAB mathematical toolbox documentation

[Matt J]

Thanks Alan,

That does explain the number of iterations. However, I don't think the initial call to F() in the line

ydata=F(xtrue,Adata);

could explain the extra 100 function evaluations. For one thing, this line only calls F() once. For another, I reset TRUECOUNTS to zero later in the code.

Is it possible that lsqcurvefit is doing a final jacobian calculation in case the user requests 7 output arguments? If so, it would make sense to suppress this step when nargout<7.

[Alan Weiss]

Oh, right, I missed that line.

I looked through the code again, and there is an additional 100 function evaluations due to the initial gradient approximation. This is different than the initial function evaluation, and it calls your objective function another 100 times. You can see that the internal reported number of function evaluations is correct at the end, though the iterative display for some reason doesn't show that initial 100 evaluations. Change your code like this to see what I mean:

function [sol,a,b,eflag,outpt,TRUECOUNTS] = test(M)
   
   if nargin<1, M=0; end
    N=100;
    xtrue=(1:N).';  %true parameter values
    Adata=rand(N);
    
    TRUECOUNTS=0;
    ydata=F(xtrue,Adata);
    lb=-inf(N,1); ub=+inf(N,1);
    
   
    known=1:M;               %The first M variables are known      
    lb(known)=xtrue(known); 
    ub(known)=xtrue(known);
    opts=optimoptions(@lsqcurvefit,  'Display','iter','MaxIterations',5);
    
    TRUECOUNTS=0;
    [sol,a,b,eflag,outpt] = lsqcurvefit(@F,zeros(N,1),Adata,ydata,lb,ub,opts);
    disp ' ' 
    function y=F(x,A)
        y=A*x;
        
        TRUECOUNTS=TRUECOUNTS+1;
        
    end
end

Then call it like this:

[~,~,~,eflag,outpt,TRUECOUNTS]=test(1)

You will see that outpt.funcCount = 800.

Alan Weiss

MATLAB mathematical toolbox documentation

[Matt J]

OK, but it seems unnecessary isn't it? When there are N-M unknowns with M>1 you really shouldn't need to do a full N-dimensional Jacobian calculation, except maybe after the final iteration, and then only if the user has requested the final jacobian as an output argument.