Error in using lsqcurvefit with constrained parameters.

5 次查看（过去 30 天）
Warren Boschen
Warren Boschen 2023-1-19
评论： Matt J 2023-1-20
Hi all,
I'm attempting to use the function lsqcurvefit to fit a single parameter, however I am providing the function with two other parameters that have been fixed using this previous post as a guide. Here is my code:
opts2 = optimset('Display', 'off', 'FinDiffType', 'central');
free_water = zeros(X, Y, Z); % The fraction of non-free water (unitless).
for z = 1:length(slices) % For every particular voxel...
for y = 1:Y
for x = 1:X
if mask(x, y, slices(z)) == 0 % Set all values outside of the mask to be 0.
free_water(x, y, z) = 0;
else
S_low_avg_formatted = zeros(low_num_unique, 1); % Reformat the average signal within a voxel.
for b = 1:low_num_unique % Vectorize eventually!!!
S_low_avg_formatted(b, 1) = S_avg_low(x, y, z, b);
end
params0 = [0; lambda_perp(x, y, z); deltaLambda(x, y, z)];
fixedset = [0 1 1];
fixedvals = [lambda_perp(x, y, z), deltaLambda(x, y, z)];
params0_vary = params0(~fixedset);
free_water(x, y, z) = lsqcurvefit(@(fw, data) low_wrapper(fw, data, fixedset, fixedvals),...
params0_vary, double(low_unique), S_low_avg_formatted, 0, 1, opts2);
end
end
end
end
% --------------------------------------------------------------------------------------------------------------
function [signal] = sm_signal_low(f,perp,delta_lambda,bvalues)
D = 3.00e-3; % D is the diffusivity of free water.
signal = (1-f)*exp(-bvalues*D)+ f*(exp(-bvalues*perp).*sqrt(pi./(4*bvalues*(delta_lambda)))....
*erf(sqrt(bvalues*(delta_lambda))));
end
function pred = low_wrapper(xvary, data, fixedset, fixedvals)
x = zeros(size(fixedset));
x(fixedset) = fixedvals;
x(~fixedset) = xvary;
pred = sm_signal_low(x,data);
end
And here is the error that I am receiving:
Array indices must be positive integers or logical values.
Error in smt_estimate>low_wrapper (line 343)
x(fixedset) = fixedvals;
Error in smt_estimate>@(fw,data)low_wrapper(fw,data,fixedset,fixedvals)
Error in lsqcurvefit (line 213)
initVals.F = feval(funfcn_x_xdata{3},xCurrent,XDATA,varargin{:});
Error in smt_estimate (line 279)
free_water(x, y, z) = lsqcurvefit(@(fw, data) low_wrapper(fw, data, fixedset, fixedvals),...
Caused by:
Failure in initial objective function evaluation. LSQCURVEFIT cannot continue.
I have a similar loop in a prior portion of the code which works fine, however I am not fixing any of the parameters in that fit so I believe this issue comes from an incorrect use of the low_wrapper() function. What am I doing wrong? For reference, lambda_perp and deltaLambda are 128x128x6 array, low_unique is a 5x1 array, and S_low_avg_formatted is also a 5x1 array. Neither lambda_perp nor deltaLambda have any values below 0, and I believe that any 0's would be excluded by the if/else statement.
Thank you!
-Warren

请先登录，再回答此问题。

采纳的回答

Torsten
Torsten 2023-1-19
移动：Torsten 2023-1-19
Ran in:
Maybe
fixedset = logical([0 1 1]);
fixedvals = [3 4];
xvary = 33;
x = zeros(size(fixedset));
x(fixedset) = fixedvals;
x(~fixedset) = xvary;
x
x = 1×3
33 3 4

更多回答（1 个）

Matt J
Matt J 2023-1-19
编辑：Matt J 2023-1-19
I think you can set the upper bounds equal to the lower bounds to fix the parameters as the OP in the previous post originally wanted to do. lsqcurvefit in recent Matlab versions will now do the reformulation you are attempting to do automatically. It probably won't be as computationally efficient, however, as reformulating it yourself.
  6 个评论
显示 4更早的评论
Matt J
Matt J 2023-1-20
Ran in:
We would have to see what code you ran, but it seems more likely that the sqrt() operations in your objective are generating complex values.
erf(1i)
Error using erf
Input must be real and full.
Your objective function should probably have a check in it anyway to verify that bvalues*delta_lambda is real and non-negative, since that is clearly what your model assumes.
Matt J
Matt J 2023-1-20
Do you know when this became true? I'm using MATLAB R2019a so I'm not sure.
It was definitely there in R2019a:
https://www.mathworks.com/help/releases/R2019a/optim/ug/lsqcurvefit.html

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Programming 的更多信息

标签

产品


版本

R2019a

Community Treasure Hunt

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

Start Hunting!

Translated by