How can I use lsqcurvefit with a look-up table search as the objective function? Currently, the initial guess always is identified as a local minimum.

3 views (last 30 days)
Samuel Spink
Samuel Spink on 20 Oct 2021
Edited: Alan Weiss on 26 Oct 2021
I have tried lowering the step, function, and optimality tolerance to 1e-12, and the initial guess is still always a local minimum, regardless of what the initial guess is (at least for the few initial guesses that I have tried). Below is some information about the code I have written, with commentary in bold and actual code not in bold.
The "LUT" (lookup table) variable that is loaded has the following properties:
LUT =
struct with fields:
Mua: [1000×1000 double]
Musp: [1000×1000 double]
M7: [1000×1000 double]
M10: [1000×1000 double]
M13: [1000×1000 double]
M16: [1000×1000 double]
Mueff: [1000×1000 double]
Below is the function that I am attempting to use as the objective function with lsqcurvefit, which calls the LUT:
function [R] = p1seminf_Reflectance_spec_MCLUT(constants,vars)
f_water=constants(1);
f_lipid=constants(2);
a_scat=constants(3);
b_scat=constants(4);
lambda0=constants(5);
SD_all=vars(:,1);
lambda_all=vars(:,2);
x=load('chromophores_SWIR.txt');
wv=x(:,1);
waterExt=x(:,2);
lipidExt=x(:,3);
load LUT_CW_multiDist_g0p9_n1p4_sd_7_10_13_16mm_dense
for i=1:size(vars,1)
SD=SD_all(i);
lambda=lambda_all(i);
if SD==7
R_mat=LUT.M7;
elseif SD==10
R_mat=LUT.M10;
elseif SD==13
R_mat=LUT.M13;
elseif SD==16
R_mat=LUT.M16;
end
wvInd_SWIR=find(wv==lambda);
E_water=waterExt(wvInd_SWIR);
E_lipid=lipidExt(wvInd_SWIR);
mua_temp=f_water*E_water+f_lipid*E_lipid;
musp_temp=a_scat*(lambda/lambda0)^(-b_scat);
muaInd=find(abs(mua_temp-LUT.Mua(1,:))==min(abs(mua_temp-LUT.Mua(1,:))));
muspInd=find(abs(musp_temp-LUT.Musp(:,1))==min(abs(musp_temp-LUT.Musp(:,1))));
R(i,1)=R_mat(muspInd,muaInd);
end
end
Lastly, here is the portion of code that uses lsqcurvefit:
fitfunc=@(p,pdata) p1seminf_Reflectance_spec_MCLUT([p(1:4) lambda0],pdata);
[fitparams_SWIR,resnorm_SWIR,residual_SWIR,exitflag_SWIR,output_SWIR,...
lambda_SWIR,jacobian_SWIR]=lsqcurvefit(fitfunc,[water_guess lipid_guess...
A_guess b_guess],vars_SWIR,R_SWIR,[],[],options);
Any thoughts on why this does not seem to work and potential alternatives would be appreciated! Thanks!

Answers (1)

Alan Weiss
Alan Weiss on 22 Oct 2021
To use lsqcurvefit your function needs to be defined and differentiable for continuous values of p. It looks to me as if your function is locally constant, a stair-step sort of function. When lsqcurvefit tries to optimize, it takes tiny steps delta and evaluates [fitfunc(p+delta,pdata) - fitfunc(p,pdata)]/delta, an approximation to the gradient. Your function returns exactly zero for this approximate gradient, so lsqcurvefit geets stuck at the start.
What can you do? Perhaps you can interpolate your function so that it is defined on continuous values, and has nonzero gradient. Or, instead of using a derivative-based solver, try using a Global Optimization Toolbox solver such as patternsearch that does not care about gradients. You will have to rewrite your objective function as the sum of squares to get a scalar objective function. For an example that uses both interpolation and the patternsearch solver, see Pattern Search Climbs Mount Washington.
Good luck,
Alan Weiss
MATLAB mathematical toolbox documentation
  2 Comments
Alan Weiss
Alan Weiss on 26 Oct 2021
All I meant was that for lsqcurvefit you are supposed to return a vector function value. This value is implicitly squared and summed. If you use a general optimization routine such as patternsearch, you must square and sum the values yourself. For an example showing both kinds of objectives, see Nonlinear Data-Fitting.
Alan Weiss
MATLAB mathematical toolbox documentation

Sign in to comment.

Products


Release

R2018a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!