Bethany Baxter
Bethany Baxter

Reputation: 31

Matlab Curve Fitting via Optimization

I have tried to follow this tutorial to fit a curve to my dataset. The equation for the curve should be

f(t) = log10((wpmcoeff./(t.^2)) + 
       ((1.038+3.*log(2.*pi.*1e6.*t)).*fpmcoeff./(t.^2))+(wfmcoeff./t) + 
       (ffmcoeff)+(rwfmcoeff.*t)).

I have created the following code:

clock='atomicclockgpsworld.txt';
data=importdata(clock);

carrier=10e6;

sig=data(:,2);
t=data(:,1);
sigsq=log10(sig.^2);
fun = @(coeff)sseval(coeff,t,sigsq);
x0 = rand(5,1);
bestx = fminsearch(fun,x0);
wpmcoeff = bestx(1);
fpmcoeff = bestx(2);
wfmcoeff = bestx(3);
ffmcoeff = bestx(4);
rwfmcoeff = bestx(5);
yfit=log10((wpmcoeff./(t.^2))+((1.038+3.*log(2.*pi.*1e6.*t)).*fpmcoeff./(t.^2))+(wfmcoeff./t))+(ffmcoeff)+(rwfmcoeff.*t);
semilogx(t,sigsq,'x');
hold on
semilogx(t,yfit);
saveas(gcf,'fit','png');

and the corresponding function

function sse = sseval(coeff,t,sigsq)
    wpmcoeff = coeff(1);
    fpmcoeff = coeff(2);
    wfmcoeff = coeff(3);
    ffmcoeff = coeff(4);
    rwfmcoeff = coeff(5);
    sse = sum(sigsq - (log10((wpmcoeff./(t.^2))+((1.038+3.*log(2.*pi.*1e6.*t)).*fpmcoeff./(t.^2))+(wfmcoeff./t)+(ffmcoeff)+(rwfmcoeff.*t))));
end

But the fit produced is horrible (my y data should vary between approximately -20 to -22 but the fit produces a curve that reaches 1e59!). Can anyone suggest where I may be going wrong?

Current output vs data:

Link to the plot that is created.

Upvotes: 3

Views: 683

Answers (1)

ViG
ViG

Reputation: 1868

In your function sseval the function is different then the one in your first script. In the function you take the log10 of the whole equation, while in the script log10 stops at (wfmcoeff./t):

first script:

log10((wpmcoeff./(t.^2))+((1.038+3.*log(2.*pi.*1e6.*t)).*fpmcoeff./(t.^2))+(wfmcoeff./t))+(ffmcoeff)+(rwfmcoeff.*t)

function:

log10((wpmcoeff./(t.^2))+((1.038+3.*log(2.*pi.*1e6.*t)).*fpmcoeff./(t.^2))+(wfmcoeff./t)+(ffmcoeff)+(rwfmcoeff.*t))

A second thing is that you didn't take the square of the difference. So in your function change the last line to

sse = sum((sigsq - (log10((wpmcoeff./(t.^2))+((1.038+3.*log(2.*pi.*1e6.*t)).*fpmcoeff./(t.^2))+(wfmcoeff./t))+(ffmcoeff)+(rwfmcoeff.*t))).^2);

Note: the fitting is sometimes good, sometimes it's really rubbish, so run the script a couple of times to get a good fitting. After a few times I got the following:

enter image description here

The parameters found are:

wpmcoeff: 10.898483535691309
wfmcoeff: 22.933722400252414
rwfmcoeff: 3.601059651996531e-05
fpmcoeff: -0.324473127267299
ffmcoeff: -21.497862719234053

Upvotes: 3

Related Questions