MATLAB lsqcurvefit gives a line instead of a sinusoid

Question

I'm trying to use lsqcurvefit to fit a sinusoid with some disturbs added:

A=3.75;
omega=2;
phi=2;
t=1:10000;
y=A*sin(omega*t/1000+phi);

noise1=(rand(1,10000)-0.5)*0.2;

noise2=0.1*sin(2*t);

sig_out=A*sin(omega*t/1000+(phi-0.5))+noise1+noise2;
figure;
plot(t,y);
hold on;
plot(t,sig_out);
grid on;

The code above create a random noise, a small oscillation around the value of the theoretical curve, and a phase I've used this code to perform the fitting:

close all;
f=@(x,xdata) x(1)*sin(x(2)*xdata+x(3));
x0=[3.75 2 2];
%definiamo i limiti inferiore e superiore della regressione
lob=ones(1,10000)*x0(1)*1.05*-1;
upb=ones(1,10000)*x0(1)*1.05;

options=optimoptions('lsqcurvefit','Diagnostic','on','MaxIteration',1000000000,'Display','iter-detailed','FunctionTolerance',1e-100,'FiniteDifferenceType','central','StepTolerance',1e-100,'FiniteDifferenceStepSize',100);
[reg,EXITFLAG]=lsqcurvefit(f,x0,t,sig_out,lob,upb,options);

but the output gives me only a straight line centered around 0 with small oscillations.. I've tried different options. What am I missing? Many thanks

user_3181104 · Answer

It appears that lsqcurvefit cannot be used to solve the problem in general.

You can start with a very close initial condition:

x0=[4 1/500 2];
x = lsqcurvefit(f,x0,t,sig_out),
plot(t,f(x,t))

and see that you do get a solution of sort. But the solver still complains about local minimum.

And if you try other initial conditions, you would notice the resultant coefficient would barely change from whatever initial conditions you put in. That is a sign of the solver not being able to solve the problem at hand. To some extend, that is to be expected. The phase is periodical. So the objective function will be periodical relative to changes to phase. Frequency is not periodic but you can imagine varying frequency yields oscillatory changes in the objective function as well. That seems to be beyond what lsqcurvefit is designed for. Perhaps consider Fourier transform or other techniques for waveform analysis first.

user_3460361 · Answer

Fitting a sinusoid is always tricky. There are several methods but one of my favorite (which uses only built-in function) is descibed in this answer in MATLAB Central: Curve fitting to a sinusoidal function

With your data, it would go like this:

yu = max(y);
yl = min(y);
yr = (yu-yl);                               % Range of ‘y’
yz = y-yu+(yr/2);

zti = find(diff(sign(yz)));                 % find zero crossings
zt = t(zti);

per = 2*mean(diff(zt));                     % Estimate period
ym = mean(y);                               % Estimate offset
fit = @(b,t)  b(1).*(sin(2*pi*t./b(2) + 2*pi/b(3))) + b(4);    % Function to fit
fcn = @(b) sum((fit(b,t) - y).^2);                             % Least-Squares cost function
s = fminsearch(fcn, [yr;  per;  -1;  ym])                      % Minimise Least-Squares

The elements of output parameter vector s ( b in the function ) are:

s(1): sine wave amplitude (in units of y)
s(2): period (in units of x)
s(3): phase (phase is s(2)/(2*s(3)) in units of x)
s(4): offset (in units of y)

This will produces the following fit to your original data:

figure;
plot(t,y,'--k','linewidth',3,'DisplayName','Original sine');
hold on;
plot(t,sig_out,'b','DisplayName','Noisy sine');
tp = linspace(min(t),max(t));
plot(tp,fit(s,tp), 'r','DisplayName','fitted sine');
grid on ; legend show

MATLAB lsqcurvefit gives a line instead of a sinusoid

2 Answers