Improving LSQCcurvefit solutions through improving x0 through loops

Question

Sam Mahdi 2019-3-12

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/449628-improving-lsqccurvefit-solutions-through-improving-x0-through-loops

编辑： Sam Mahdi 2019-3-13

Hello everyone,

I had a quick question in regards to attempting to improve my lsqccurvefit values by improving x0 through multiple datasets.

Data1=[0.596421471	0.06
0.5859375	0.119284294
0.566037736	0.29296875
0.530035336	0.622641509
0.418994413	1.219081272
0.388619855	3.058659218
];
Data2=[5.00E+04	3.82E+04	3.45E+04	3.42E+04	3.74E+04	3.21E+04	2.81E+04	2.88E+04
5.00E+04	3.82E+04	3.45E+04	3.42E+04	3.74E+04	3.21E+04	2.81E+04	2.88E+04
3.08E+04	3.07E+04	2.19E+04	2.23E+04	2.53E+04	2.05E+04	1.98E+04	1.89E+04
2.05E+04	1.87E+04	1.30E+04	1.10E+04	1.62E+04	1.31E+04	1.05E+04	1.05E+04
8963	1.11E+04	6243	3504	6454	4346	4448	4357
0.00E+00	0.00E+00	0.00E+00	0.00E+00	0.00E+00	0.00E+00	0.00E+00	0.00E+00
];
Pt=Data1(:,1);
Lt=Data1(:,2);
for n = 1:size( Data2, 2 )
   Int=Data2(:,n);
  ...
plot(Lt,Int, 'ro');
xdata=[Pt,Lt]
F=@(x,xdata)((xdata(:,1)+xdata(:,2)+x(2))-((xdata(:,1)+xdata(:,2)+x(2)).^2-4.*xdata(:,1).*xdata(:,2).^1/2).*x(1)/2.*xdata(:,1))
x0=[1 1]
options=optimoptions('lsqcurvefit','MaxFunEvals',5e2)
lb=[];
up=[]
[x,resnorm,~,exitflag,output]=lsqcurvefit(F,x0,xdata,Int,lb,up,options)
options = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt');
lb = [];
ub = [];
[x2,resnorm2,~,exitflag2,output2]=lsqcurvefit(F,x0,xdata,Int,lb,up,options)
M=mean(x) 
M=mean(x2)
hold on 
plot (Lt, F(x,xdata))
plot (Lt, F(x2,xdata))
hold off
end

In a nutshell, I have a function F that is trying to fit my Data 2 using 2 different algorithms. My end plot currently has a poor fit with my initial plot. I've seen that the initial starting values [x0] can have a signifigant impact on your results. Is there a way to have the later iterations in my loop change the x0 value and sample various x0 values and try and get closer to the real value that fits my data well (this can be done by doing an RMSD or something like that with my initial plot and vlues). The simplest way I was thinking of having this done is using having x0 change to the results of the previous run. I.E. The first run using the x0 I gave it, the 2nd x0 that is used in the 2nd run uses the results of the previous run as it's own x0.

I don't quite know how to do this, whether this is a modification through x0 (like my example), or through using a loop within my current loop that runs the same dataset n=1 repeatedly through multiple iterations using the RMSD example I gave above to try and improve it. Then do the same thing for the 2nd dataset (n=2), although I can see this method taking quite a while and maybe potentially creating an endless loop.

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Answer 1

Matt J 2019-3-12

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/449628-improving-lsqccurvefit-solutions-through-improving-x0-through-loops#answer_365016

在 MATLAB Online 中打开

For such a small problem, it is a simple matter to do a vectorized exhaustive search, rather than looping

[X1,X2]=meshgrid(linspace(-50000,0,500), linspace(-60,60,100));
[mm,nn]=size(X1);
XX1=reshape(X1,1,1,[]);
XX2=reshape(X2,1,1,[]);
Y=reshape(vecnorm(F([XX1,XX2],xdata)-Int,inf,1),mm,nn);
[minval,idx]=min(Y(:)),
x0=[X1(idx),X2(idx)]

17 个评论
显示 15更早的评论隐藏 15更早的评论

Matt J 2019-3-12

编辑：Matt J 2019-3-12

在 MATLAB Online 中打开

Exhaustive search just means to search for the best x by doing a sweep over all parameter combinations. In the code I've shown you the objective function is calculated at all lattice points on some 2D grid search grid (except I do it without using explicit loops - instead I just use Matlab's matrix arithmetic) and I can just take the minimum over all the sample points. Exhaustive is a bit of a misnomer here - it's not really exhaustive because you have a continuum of feasible solutions, and this search sweeps through a discretization of that. Also, you still have to guess the grid box size big enough to contain the optimal solution. However, the search can give you an x0 as close to the global optimum as the fineness of the search grid allows, and also tells you approximately the best fitting accuracy (minval) that you can achieve.

You can plot the objective function at all the search points by adding a few lines to the code,

[X1,X2]=meshgrid(linspace(-1e5,1e5,100), linspace(-600,600,100));
[mm,nn]=size(X1);
XX1=reshape(X1,1,1,[]);
XX2=reshape(X2,1,1,[]);
Y=reshape(vecnorm(F([XX1,XX2],xdata)-Int,inf,1),mm,nn);
[minval,idx]=min(Y(:)),
x0=[X1(idx),X2(idx)]
surf(X1,X2,Y);scf
xlabel x1, ylabel x2

Sam Mahdi 2019-3-12

在 MATLAB Online 中打开

Thank you for the explanation. So from my understanding, your code defines x0 as being 2 variables that are searching a grid with the size you provided for the local minimum. So I included this into my loop, since I'd desirably want this done to every single datapoint that I have

for n = 1:size( Data2, 2 )
   Int=Data2(:,n);
  ...
hold on
plot(Lt,Int, 'ro');
xdata=[Pt,Lt]
F=@(x,xdata)((xdata(:,1)+xdata(:,2)+x(2))-((xdata(:,1)+xdata(:,2)+x(2)).^2-4.*xdata(:,1).*xdata(:,2).^1/2).*x(1)/2.*xdata(:,1))
options=optimoptions('lsqcurvefit','MaxFunEvals',5e2)
lb=[];
up=[]
[x,resnorm,~,exitflag,output]=lsqcurvefit(F,x0,xdata,Int,lb,up,options)
options = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt');
lb = [];
ub = [];
[x2,resnorm2,~,exitflag2,output2]=lsqcurvefit(F,x0,xdata,Int,lb,up,options)
[x1,x2]=meshgrid(linspace(-50000,0,500), linspace(-60,60,100));
[mm,nn]=size(x1);
xx1=reshape(x1,1,1,[]);
xx2=reshape(x2,1,1,[]);
Y=reshape(vecnorm(F([xx1,xx2],xdata)-Int,inf,1),mm,nn);
[minval,idx]=min(Y(:)),
x0=[x1(idx),x2(idx)]
M=mean(x) 
M=mean(x2) 
plot (Lt, F(x,xdata))
plot (Lt, F(x2,xdata))
hold off
end

And I got this error.

Error using reshape
To RESHAPE the number of elements must not
change.
Error in Testscript4loop (line 36)
Y=reshape(vecnorm(F([xx1,xx2],xdata)-Int,inf,1),mm,nn); 
>> 

Initially you had everything as capital x's, so I changed them all to lower case (since in my function, it was a lower case x, although I don't know if it's cap sensitive or not. Then I thought it might be an issue since you have it written as x1 and x2, but I don't technically have an x1 in y script (unless x is automatically assumed to be x1 by the program), but I still got the exact same error. I did the reverse as well, changing all the x1s to just x, and that also did not fix it. So it doesn't appear to be an issue regarding syntax.

Sam Mahdi 2019-3-13

在 MATLAB Online 中打开

Thank you, this definitely helped my function fit my data better. However, a few questions again.

1) I noticed the main difference here was the change in the my function (so xdata 2nd column [xdata :,2] is now xdata 1,2, ) and my x(1) variable is now x(1,1,:). I'm curious what this is doing exactly, and why you only changed one of the xdata (:,2) and not the other ones. These changes seem to make distinct changes to my output, but I don't exactly understand what's going on. As well as why xdata (:,1) isn't touched.

F=@(x,xdata)((xdata(:,1)+xdata(:,2) [why this isn't changed]+x(2))-((xdata(:,1)+xdata(1,2,:)[but this is changed]+x(2)).^2-4.*xdata(:,1).*xdata(:,2).^1/2).*x(1,1,:)/2.*xdata(:,1))

2) While my fit is much better now, my results are still very poor. The issue with this is I get negative values; however these results are actually indicating concentration, so negative values don't actually make sense. I have attempted to resolve this by inputing a lower bound of 0

lb=[0];

for my functions but when I do, I get this error.

Warning: Length of lower bounds is <
length(x); filling in missing lower bounds
with -Inf. 

And I can't use lower bounds for the levenberg algorithm since it doesn't seem to accept it (which is a shame, since this is the one that's actually fitting my data decently). I believe this is the main other issue I'm now at, simply that the function should approach zero (and this will give it the nice logarithmic curve), but as it stands, the function goes straight through zero, thus giving me the negative values for x.

Matt J 2019-3-13

编辑：Matt J 2019-3-13

I'm curious what this is doing exactly, and why you only changed one of the xdata (:,2) and not the other ones.

I made a mistake. The redefinition should look like

F=@(x,xdata)((xdata(:,1)+xdata(:,2)+x(1,2,:))-((xdata(:,1)+xdata(:,2)+x(1,2,:)).^2-4.*xdata(:,1).*xdata(:,2).^1/2).*x(1,1,:)/2.*xdata(:,1));

In other words, only the indexing of x should change. This makes it so that you can evaluate a stack of pairs [x1,x2] of dimension 1x2xN simultaneously.

I have attempted to resolve this by inputing a lower bound of 0

You have 2 unknowns, so it should be lb=[0,0]. You should also choose a search grid that excludes negative values, perhaps something like

[x1,x2]=meshgrid(linspace(0,1e4,500), linspace(0,600,100));

While my fit is much better now, my results are still very poor.

That doesn't necessarily indicate a failure in lsqcurvefit. Your model equation and your data could simply be badly matched, i.e., this might be the best fit you can hope for. The minval variable computed by my code will give you the approximate magnitude of the fitting errors you should expect, assuming the search grid is chosen large enough to contain the optimal solution.

It might also be worth mentioning that the surface plot that I presented to you for your first data set shows a very suspicious shape. It is flat looking over a very large region. Essentially, every [x1,x2] in that region looks like it would give nearly the same curve. This makes me suspicious even moreso about the validity of the model equation.

And I can't use lower bounds for the levenberg algorithm since it doesn't seem to accept it (which is a shame, since this is the one that's actually fitting my data decently)

Your other technique limits the number of iterations greatly. This is usually a bad idea - it stopes the solver sooner than it wants to and may be the reason the resulting fit is poorer.

Sam Mahdi 2019-3-13

在 MATLAB Online 中打开

Here is my script so far.

Data1=[0.596421471	0.06
0.5859375	0.119284294
0.566037736	0.29296875
0.530035336	0.622641509
0.418994413	1.219081272
0.388619855	3.058659218
];
Data2=[5.00E+04	3.82E+04	3.45E+04	3.42E+04	3.74E+04	3.21E+04	2.81E+04	2.88E+04
5.00E+04	3.82E+04	3.45E+04	3.42E+04	3.74E+04	3.21E+04	2.81E+04	2.88E+04
3.08E+04	3.07E+04	2.19E+04	2.23E+04	2.53E+04	2.05E+04	1.98E+04	1.89E+04
2.05E+04	1.87E+04	1.30E+04	1.10E+04	1.62E+04	1.31E+04	1.05E+04	1.05E+04
8963	1.11E+04	6243	3504	6454	4346	4448	4357
0.00E+00	0.00E+00	0.00E+00	0.00E+00	0.00E+00	0.00E+00	0.00E+00	0.00E+00
];
Pt=Data1(:,1);
Lt=Data1(:,2);
for n = 1:size( Data2, 2 )
   Int=Data2(:,n);
  ...
hold on
plot(Lt,Int, 'ro');
xdata=[Pt,Lt]
F=@(x,xdata)((xdata(:,1)+xdata(:,2)+x(1,2,:))-((xdata(:,1)+xdata(:,2)+x(1,2,:)).^2-4.*xdata(:,1).*xdata(:,2).^1/2).*x(1,1,:)/2.*xdata(:,1))
[x1,x2]=meshgrid(linspace(0,1e4,5e4), linspace(0,60,100));
[mm,nn]=size(x1);
xx1=reshape(x1,1,1,[]);
xx2=reshape(x2,1,1,[]);
Y=reshape(vecnorm(F([xx1,xx2],xdata)-Int,inf,1),mm,nn);
[minval,idx]=min(Y(:)),
x0=[x1(idx),x2(idx)]
options=optimoptions('lsqcurvefit','MaxFunEvals',5e2)
lb=[0,0];
up=[]
[x,resnorm,~,exitflag,output]=lsqcurvefit(F,x0,xdata,Int,lb,up,options)
options = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt');
lb = [];
ub = [];
[x2,resnorm2,~,exitflag2,output2]=lsqcurvefit(F,x0,xdata,Int,lb,up,options)
M=mean(x) 
M=mean(x2)
plot (Lt, F(x,xdata))
plot (Lt, F(x2,xdata))
hold off
end

And yes, I do expect the x2 values to be close to zero (looking at hopefully somewhere in the 0.1-0.001 range.

Matt J 2019-3-13

在 MATLAB Online 中打开

Note that x.^1/2 is not the square root of x. Compare:

>> 2.^1/2
ans =
     1
>> sqrt(2)
ans =
    1.4142

Sam Mahdi 2019-3-13

编辑：Sam Mahdi 2019-3-13

在 MATLAB Online 中打开

Oh, well that makes a huge difference. With that modification, my graph is no longer in the negative. Although, now it looks even weirder! I should note, I had forgotten to normalize my data prior, it won't change the plot, but does make the values fit better to the quadratic function you have below.

Data2=[0.914458052	0.698775361	0.630597697	0.625479803	0.684335588	0.587461159	0.512703345	0.526777554
914458052	0.698775361	0.630597697	0.625479803	0.684335588	0.587461159	0.512703345	0.526777554
535527691	0.534656914	0.38174852	0.388714734	0.440961338	0.356844305	0.34430512	0.328631139
401724476	0.36664707	0.253968254	0.215755438	0.316872428	0.257103665	0.206153243	0.204977464
182397232	0.224867725	0.127045177	0.071306471	0.131339031	0.088441188	0.090516891	0.088665039
0	0	0	0	0	0	0
];

请先登录，再进行评论。

Answer 2

Matt J 2019-3-13

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/449628-improving-lsqccurvefit-solutions-through-improving-x0-through-loops#answer_365225

编辑：Matt J 2019-3-13

在 MATLAB Online 中打开

I can see that your model function derives in some way from the quadratic formula, but I can't quite see where all the pieces are supposed to go. In any case, using the implementation below, I get a fit that approaches sensible. Maybe you can tweak it to perfection.

for n = 1:size( Data2, 2 )
   Int=Data2(:,n);
 
xdata=[Pt,Lt]
%F=@(x,xdata)((xdata(:,1)+xdata(:,2)+x(1,2,:))-((xdata(:,1)+xdata(:,2)+x(1,2,:)).^2-4.*xdata(:,1).*xdata(:,2).^1/2).*x(1,1,:)/2.*xdata(:,1))
F=@modelfun;
[X1,X2]=meshgrid(linspace(0,1e4,500), linspace(0,60,100));
[mm,nn]=size(X1);
xx1=reshape(X1,1,1,[]);
xx2=reshape(X2,1,1,[]);
Y=reshape(vecnorm(F([xx1,xx2],xdata)-Int,inf,1),mm,nn);
[minval,idx]=min(Y(:)),
x0=[X1(idx),X2(idx)];
options = optimoptions('lsqcurvefit','Algorithm','levenberg-marquardt');
[x2,resnorm2,~,exitflag2,output2]=lsqcurvefit(F,x0,xdata,Int,[],[],options)
plot (Lt,Int,'ro',Lt, F(x2,xdata))
xlabel 'Lt'
end
function out=modelfun(x,xdata)
     Pt=xdata(:,1);
     Lt=xdata(:,2);
     x1=x(1,1,:);
     x2=x(1,2,:);
     
     A=Lt.*x1;
     B=-(Pt+Lt+x2);
     C=Pt;
     out = ( -B + sqrt(B.^2-4.*A.*C) )./(2.*A);
 
end

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Sam Mahdi 2019-3-13

It is the quadratic equation. The only addition is the quadratic is being multiplied by a variable (x1). You have the right set up in your script in terms of what everything stands for, except x1 is in the numerator, and you currently have it in the denominator. This is why I'm confused as to why my graphs look funky, because if you plot it out using it the way you have it written, you get a nice curve. While my function is essentially the same thing, but provides the above funky results.

请先登录，再进行评论。

Improving LSQCcurvefit solutions through improving x0 through loops

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

回答（2 个）

17 个评论
显示 15更早的评论隐藏 15更早的评论

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

另请参阅

类别

标签

Community Treasure Hunt

Improving LSQCcurvefit solutions through improving x0 through loops

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

回答（2 个）

17 个评论 显示 15更早的评论隐藏 15更早的评论

1 个评论 显示 -1更早的评论隐藏 -1更早的评论

另请参阅

类别

标签

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

17 个评论
显示 15更早的评论隐藏 15更早的评论

1 个评论
显示 -1更早的评论隐藏 -1更早的评论