How do I use the Levenberg-Marquardt algorithm for lsqcurvefit?
13 次查看(过去 30 天)
显示 更早的评论
I'm trying to run a fit on some data using lsqcurvefit. Here is the basic structure of the code:
opts = optimset('Display', 'off');
params0 = [0; 0.60e-4; 1.30e-3];
lb = [0, 0.40e-4, 0.92e-4];
ub = [1, 0.80e-4, 1.76e-4];
[params, sigval] = lsqcurvefit(@sm_signal, params0, double(b_unique),...
S_avg_formatted, lb, ub, opts);
function [signal] = sm_signal(params,bvalues)
% Assign parameters from params for use within the function sm_signal.
[f, perp, delta_lambda] = deal(params(1),params(2),params(3));
D = 3.00e-3;
v = 0.01;
signal = (1-f)*exp(-bvalues*D) + ...
f*(exp(-bvalues*perp).*sqrt(pi./(4*bvalues*(delta_lambda))).*erf(sqrt(bvalues*(delta_lambda)))) + ...
v*(perp/delta_lambda);
end
In general, this works well. However, I have a specific data set that I want to fit which has fewer b-values than parameters to be fit. That is to say, size(b_unique) = [2 1] and size(params0) = [3 1]. Trying to fit the data gives the following error:
Error using lsqncommon (line 64)
The Levenberg-Marquardt algorithm does not handle bound constraints and the trust-region-reflective algorithm requires at least as many equations as
variables; aborting.
Error in lsqcurvefit (line 271)
lsqncommon(funfcn,xCurrent,lb,ub,options,defaultopt,allDefaultOpts,caller,...
I've been reading the documentation and it says that for an underdetermined system, the Levenberg-Marquardt algorithm could be used instead. But whenever I run the code using that algorithm, I receive the same error. What am I doing wrong? Is it the fact that I am imposing bounds on my parameters with lb and ub? If so, how do I fix it?
Thanks,
Warren
3 个评论
Torsten
2023-3-2
It simply doesn't make sense what you try to do. The parameters you will find are almost arbitrary and will only be useful to reproduce the two measurement data. Transfering the model to different conditions - which is the reason to build a mathematical model - will not be possible.
回答(2 个)
Matt J
2023-3-2
编辑:Matt J
2023-3-2
If the signal is supposed to be smooth, it would make sense to add some sort of roughness penalty term to the reisduals, e.g.,
opts = optimoptions('lsqnonlin','Display', 'off','Algorithm','levenberg-Marquardt');
params0 = [0; 0.60e-4; 1.30e-3];
lb = [0, 0.40e-4, 0.92e-4];
ub = [1, 0.80e-4, 1.76e-4];
bvalues=double(b_unique);
fun=@(params) sm_signal(params,bvalues,S_avg_formatted);
[params, sigval] = lsqnonlin(@sm_signal, params0, lb, ub, opts);
function residuals = sm_signal(params,bvalues,S)
% Assign parameters from params for use within the function sm_signal.
[f, perp, delta_lambda] = deal(params(1),params(2),params(3));
D = 3.00e-3;
v = 0.01;
signal = (1-f)*exp(-bvalues*D) + ...
f*(exp(-bvalues*perp).*sqrt(pi./(4*bvalues*(delta_lambda))).*erf(sqrt(bvalues*(delta_lambda)))) + ...
v*(perp/delta_lambda);
penalty=0.01*diff(signal);
residuals=[signal(:)-S(:) ; penalty(:)];
end
0 个评论
Matt J
2023-3-2
编辑:Matt J
2023-3-2
Here's a silly solution to the problem. The fit will be meaningless, however.
opts = optimset('Display', 'off');
params0 = [0; 0.60e-4; 1.30e-3];
lb = [0, 0.40e-4, 0.92e-4];
ub = [1, 0.80e-4, 1.76e-4];
fun=@(p,b) repmat(sm_signal(p,b),1,2);
[params, sigval] = lsqcurvefit(@sm_signal, params0, double(b_unique),...
[S_avg_formatted, S_avg_formatted], lb, ub, opts);
2 个评论
Matt J
2023-3-2
In the sense that the solution is still under-determined. So, you won't be able to make predictions with it, see also Torsten's comment.
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Interpolation 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!