How do I use the Levenberg-Marquardt algorithm for lsqcurvefit?

15 次查看(过去 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 个评论
Warren Boschen
Warren Boschen 2023-3-2
I agree with you, but unfortunately this is the most readily available data at the moment. I'm hopeful that the documentation specifically cites the Levenberg-Marquardt algorithm as a way to handle such a system, but I can't get it to work in the first place.
Do you think there is any way that I can use this data? Or would it require a complete reformulation of the problem?
Torsten
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
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

Matt J
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);

类别

Help CenterFile Exchange 中查找有关 Interpolation 的更多信息

产品


版本

R2019a

Community Treasure Hunt

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

Start Hunting!

Translated by