insert two criteria function to lsqcurvefit

Hi, I have some data points as (X,Y) and I want to use lsqcurvefit function to find the best fit by least square method. My guidance function is a non-linear function Gamma2(h) as bellow:

As above equation shows, it's a two criteria function and because the value of a0 is not clearly known, I have to use lsqcurvefit function simultaneously to reach the a0, but I don't know how to do it.

Please help me to solve problem.

Thanks, Mani

Answer 1

[Star Strider]

This works with simulated data:

%  PARAMETERS: b(1) = C0,    b(2) = a0
gamma2 = @(b,h) b(1).*[(1.5*(h/b(2)) - 0.5*(h./b(2)).^3).*(h<=b(2)) +  (h>b(2))];

% SIMULATE FUNCTION
b = [10; 5];                                    % Created Parameters
h = linspace(0,7,25);                           % ‘h’ (Independent Variable)
gam2data = gamma2(b,h)+0.5*randn(1,length(h));  % Created Data (Dependent Variable)

B0 = [1; 1];
B = lsqcurvefit(gamma2, B0, h, gam2data);       % Estimate Parameters

figure(1)
plot(h, gam2data,'bp')
hold on
plot(h, gamma2(B,h), '-r')
hold off
grid

[Mani Ahmadian]

Dear Star Strider

Thanks a lot for your great guidance, you solved my problem. Now I'm able to insert all criteria at a time and it was great for me.

But I have another problem, too. About least square, the first points are more important than the last ones. (importance: First Point>Second Point>Third Point> ....>Last point)

How can I add suitable weights to my points and then use lsqcurvefit function in this way?

My real data are here:

h=
   1.1180340e+00
   2.2360680e+00
   3.3541020e+00
   4.4721360e+00
   5.5901699e+00
   6.7082039e+00

Gamma2=
   9.8000000e+01
   1.2500000e+01
   3.1100000e+01
   5.3000000e+01
   2.1250000e+01
   3.0500000e+01

I'm waiting for your great answer.

Thanks a lot, Mani

[Star Strider]

You have to provide a numerical vector of weights the same length as your data vectors. I would use nlinfit for this if you have it (Statistics Toolbox) since it can do weighted nonlinear parameter estimation. It is also relatively easy with fminsearch.

Your function does not describe your data. It gives a non-significant fit that is essentially a horizontal line through the mean of the y-values. You need more data over a wider range of ‘h’ values, ideally starting near zero. I doubt if a weight vector would improve the fit.

EDIT — For an excellent discussion on weighted least squares, see Weighted Least Squares.

[Mani Ahmadian]

Dear Star Strider

Thanks a lot again. I google about Weighted Least Square method(WLS) to reach the best fit and compare the results by previous one. But equations are so complicated for me to understand. Would you please help me more to solve the problem?

My data points are as above(h, Gamma2), about weights: As I understand the first points(lower h) are more important than the last ones, so I think consideration of weights as

w=ones(size(h)./h);

will be helpful.

Thanks

Mani

[Star Strider]

That weight vector doesn’t change anything other than the value of ‘C0’. At this point, you simply don’t have enough data over a wide enough span to estimate your function parameters correctly.

I used fminsearch because it handles functions such as yours more efficiently and it also is easier to program a weight vector into the cost function:

W = ones(size(h))./h;

WLSCF = @(b) sum(W.*(Gamma2-gamma2(b,h)).^2); % Weighted Sum-Of-Squares Cost Function
B0 = rand(2,1)*10;
B = fminsearch(WLSCF, B0);      % Estimate Parameters

figure(1)
plot(h, Gamma2,'bp')
hold on
plot(h, gamma2(B,h), '-r')
hold off
grid

[Mani Ahmadian]

Dear Star Strider

Thanks for your great cooperation and nice code. I run it and the result was so bad! I'm sure it was just for my bad weights. I'm so puzzled to solve my problem. I have to find C0 and a0 so that a0>h(1) or a0>h(2) (I used 'while-loop' to reach this clue, the result was an infinit loop...)

Every method I tried before not works and I decided to use weighted least square method to find the best results.

Now, I do not know what I can do to reach my goal.

[Star Strider]

My pleasure.

It’s not the weights. You need more data, preferably 5 or more data points at ‘h’ between 0 and 1, largely because the data you now have are extremely noisy. (I cannot guarantee that you will be able to fit your function reliably even then. It depends on how noisy your additional data are.)

[Mani Ahmadian]

Hi

I used above code on another data points and plot the result as bellow:

It's not a good fit based on my personal experiences. I think I have to change my view point to solve the problem.

Can I use regression functions to find the curve and use gamma2 as a limit for it?

gamma2 = @(b,h) b(1).*[(1.5*(h/b(2)) - 0.5*(h./b(2)).^3).*(h<=b(2)) +  (h>b(2))];

If it is possible, please help me to do.

Thanks so much

Mani

[Star Strider]

Your equation does not appear to accurately describe your data. I tried several options, including setting lower bounds on ‘b(1)’=C0=max(gam2data), and could not get a good fit.

I suggest you consider a different model. Your current model does not appear to describe the process that is generating the data you are fitting to it.

Answer 2

[Matt J]

Beware differentiability issues. Since your curve is not differentiable with respect to a0, the least squares cost function might not be either. It may be more trustworthy to use a derivative-free method like fminsearch instead,

gamma2 = @(b,h) b(1).*[(1.5*(h/b(2)) - 0.5*(h./b(2)).^3).*(h<=b(2)) +  (h>b(2))];

A0=fminsearch( @(a0) norm(gamma2(a0,X)-Y) , initial_a0  )

[Star Strider]

I agree, but we’ve been there before in set parameters of nlinfit function. I decided not to fight it this time.