Reputation: 73
I want to make this "by hand" rather than using a surface fitting tool, because depending on the data I have, the surface fitting may vary. So, I first read the data in an excel sheet, then initialize some coefficients, calculate a 3D surface (f(x,y)) and then calculate the total least squares sum, which I'd like to minimise. Every time I run the script it tells me that I'm at the local minimum, even when I change the initial values. Changing the tolerance doesn't affect the result either.
This is the code:
% flow function in a separate .m file (approximation, it’s a negative paraboloid, maybe if required, this function may vary):
function Q = flow(P1,P2,a,b,c,d,e,f)
Q1 = a-b.*P1-c.*P1.^2;
Q2 = d-e.*P2-f.*P2.^2;
Q = Q1 + Q2;
% Variable read, I use a xlsread instead
p1a = [-5, -5, -5, -5, -5, -5, -5, -5, -5, -5];
p2a = [0, 1, 2, 3, 4, 5, 6, 7, 8, 9];
qa = [10, 9, 8, 7, 6, 5, 4, 3, 2, 1];
p1b = [-6, -6, -6, -6, -6, -6, -6, -6, -6, -6];
p2b = [0, 1, 2, 3, 4, 5, 6, 7, 8, 9];
qb = [12, 11, 10, 9, 8, 7, 6, 5, 4, 3];
% Variable initialization
coef = [50, 1, 1, 10, 1, 1];
% Function calculation
q1a = flow(p1a,p2a,coef(1),coef(2),coef(3),coef(4),coef(5),coef(6));
q1b = flow(p1b,p2b,coef(1),coef(2),coef(3),coef(4),coef(5),coef(6));
% Least squares
LQa = (qa-q1a).^2;
LQb = (qb-q1b).^2;
Sa = sum(LQa);
Sb = sum(LQb);
St = Sa+Sb;
% Optimization (minimize the least squares sum)
func = @(coef)(St);
init = coef;
opt = optimoptions('fminunc', 'Algorithm', 'quasi-newton', 'Display', 'iter','TolX', 1e-35, 'TolFun', 1e-30);
[coefmin, Stmin] = fminunc(func, init, opt);
If you run this, you should get a result of 15546
for Stmin
, but if you change the coefficients, you'll get another result, and it will also be considered as a local minimum.
What am I doing wrong?
Upvotes: 1
Views: 179
Reputation: 45741
The problem is that your func
is just a constant. It simply returns a pre-calculated value, St
, which is constant no matter what input you pass to func
. Try calling func
with various different inputs to test this.
Your objective function needs to contain all the calculations that got you to St
. So I suggest you replace your func
with a function saved in an m-file looking something like this:
function St = objectiveFunction(coef, p1a, p2a, p1b, p2b, qa, qb, q1a, q1b)
% Function calculation
q1a = flow(p1a,p2a,coef(1),coef(2),coef(3),coef(4),coef(5),coef(6));
q1b = flow(p1b,p2b,coef(1),coef(2),coef(3),coef(4),coef(5),coef(6));
% Least squares
LQa = (qa-q1a).^2;
LQb = (qb-q1b).^2;
Sa = sum(LQa);
Sb = sum(LQb);
St = Sa+Sb;
end
And then in your script call objectiveFunction
using an anonymous function like this:
[coefmin, Stmin] = fminunc(@(coef)(objectiveFunction(coef, p1a, p2a, p1b, p2b, qa, qb, q1a, q1b)), init, opt);
The idea is to create an anonymous function that only takes a single parameter, coef
, which is the variable that fminunc
will peturb and pass back to your objective function. The other parameters that your objectiveFunction
needs (i.e. p1a, p2a, p1b,...
) are now considered to be pre-calculated by your anonymous function and thus by fminunc
.
The rest of your code can stay the same.
Upvotes: 1