Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

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.

like image 504
Bethany Baxter Avatar asked Jan 31 '26 05:01

Bethany Baxter


1 Answers

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
like image 199
ViG Avatar answered Feb 03 '26 08:02

ViG