fitting scatter data into multiple cosine functions

hello everyone,

i am trying to fit a function:

y = 2*a^2*cos(w*t)+2*b^2*cos(w*t+2*phi)+2*a*b*cos(delta)*cos(w*t+phi);

with y and t known, and try to get the unknown parameters of a, b,w,phi and delta. the most important parameter is phi, the rest are less important.

the scatter data looks like as follows and is also attached in the matlab.mat file.

My question is that is it really possible to a get an unique/accurate result of the fitting?

thanks.

Answer 1

[Alex Sha]

Hi， how about the results below, the fitting is perfect, the reason for you cann't get good result is more possible the bad initial start valus. One more thing is the results are not unique.:

1:

Root of Mean Square Error (RMSE): 5.86503133415183E-10

Sum of Squared Residual: 6.91411710266713E-17

Correlation Coef. (R): 1

R-Square: 1

Adjusted R-Square: 1

Determination Coef. (DC): 1

Chi-Square: 1.19308902605578E-18

Parameter Best Estimate

---------- -------------

a 1.24243011796252

w -2000.00000005105

b -0.306494082231627

phi 0.448966482806276

delta 0.369303537459062

2:

Root of Mean Square Error (RMSE): 3.17092533626903E-9

Sum of Squared Residual: 2.02100826512676E-15

Correlation Coef. (R): 1

R-Square: 1

Adjusted R-Square: 1

Determination Coef. (DC): 1

Chi-Square: -1.07854914106232E-15

F-Statistic: 1.58977775365207E19

Parameter Best Estimate

---------- -------------

a 1.08573520679119

w 2000.0000002083

b -0.131006405724609

phi 3.73119497513834

delta 30.8055670286459

3:

Root of Mean Square Error (RMSE): 8.08540777797663E-10

Sum of Squared Residual: 1.31401376061692E-16

Correlation Coef. (R): 1

R-Square: 1

Adjusted R-Square: 1

Determination Coef. (DC): 1

Chi-Square: 2.45286536720573E-17

F-Statistic: 2.44514484842302E20

Parameter Best Estimate

---------- -------------

a 0.0596400967434023

w 2000.00000007458

b -1.14950810470384

phi -12.5353535070843

delta -17.9045320865089

[an gao]

thanks Alex, unfortunately the results are not what i am looking for. Is there a way to get a unique result?

Answer 2

[David Wilson]

Yes, but it is always tricky to fit such sinusoids without good starting estimates. I'm going to use lsqcurvefit from the optimisation toolbox.

%% fit a sinusoid 
load matlab 
plot(t,y)

It may well be prudent to start with a simplified version first. The trick is to get the dominant amplitude and frequency starting estimates to be reasonable. One way is simply to measure the time between peaks and use that for the period.

[PKS,LOCS]= findpeaks(y) 
kPeriod = mean(diff(LOCS))
Ts = mean(diff(t)); 
Period = kPeriod*Ts; 
w = 2*pi/Period

That gives us a frequency of about 2000. (I'm assuming t is regularly sampled, and you have findpeaks from the signal processing toolbox. If not, you just need to measure it yourself.)

The amplitude is about 3, so 3 = 2*a^2.

% Simplified version 
f = @(p,t) 2*p(1)^2*cos(p(2)*t); 
a = sqrt(3/2); w = 2000; 
p0 = [a, w]; 

p = lsqcurvefit(f, p0, t,y) 

ti = linspace(0,0.02)';
yi = f(p,ti); 
y0 = f(p0, ti);v% initial guess for comparison 
plot(t,y,'o', ti, yi,'-', ti, y0, '--')

Now we will try your (horribly overdone and over-parameterised version)

%% function of interest
f = @(p,t) 2*p(1)^2*cos(p(3)*t)+2*p(2)^2*cos(p(3)*t+2*p(4))+2*p(1)*p(2)*cos(p(5))*cos(p(3)*t+p(4));

a = sqrt(3/2); b = 0; w = 2000; phi = 0; delta = 0; 
p0 = [a, b,w, phi, delta]; 

p = lsqcurvefit(f, p0, t,y) 

ti = linspace(0,0.02)';
yi = f(p,ti); 
y0 = f(p0, ti)
plot(t,y,'o', ti, yi,'-', ti, y0, '--')

results look OK.

Answer 3

[Alex Sha]

So, an gao, what is the result you want? Unless you add range constraint for all parameters (except w), because there are trigonometric functions in the fit formula.