least square fitting and phase shift

Hi, I have few doubts in Least square curve fitting and phase shift in Matlab !!

Here is my problem !

I have 100 values data for time and current (which is not pure sinusoidal, slightly increasing in time). I have to make a least square fitted curve out of these data.

1. what is the significance of initial values ('x0')giving in the 'lsqcurvefit'?

2. How can I give the initial values x0=[x(1) x(2) x(3)], where x(1) =dc offset, x(2) =amplitude of sine, x(3)=phase shift ? Because it is the current signal which we are going to find out, so how can we predict this one step before !!

> Equation I selected is '  x(1)+x(2)*sin(2*pi*100*t+x(3)) '

3. If I am not wrong , we are optimizing the values given to the 'x0', pls correct me if I am wrong.

4. Does the 'output' values give an optimized value for x(1),x(2),x(3) ?

5. Does the x(3) from 'output' give the phase shift of current to voltage ?

6. Does the x(1)+x(2) from 'output' of last period give the maximum value of the current ?

My ultimate aim is to find the phase shift ! (Im not using FFT)

Here is my code .

'data' has 3 columns, Time, Current, Voltage

        T= data{1} ;   % time 100 values for 10 periods
        I= data{2} ;     % current 100 values
        V= data{3} ;     % voltage 100 values (sine voltage)

        L2_T=T(end-19:end); % time points in last 2 periods 
        L2_I=I(end-19:end); % current points  in last 2 periods (selected bcy which are more sinusoidal than previous periods)
        L2_V=V(end-19:end); % voltage points  in last 2 periods (sine voltage)

        options = optimset('Display','off');

        x0 = [5 1 pi/3 ] ; % initial  values for lease sqr fit

        F = @(x,xdata)  x(1)+x(2)*sin(2*pi*100*xdata+x(3)); % least sq SINE 

        output = lsqcurvefit(F,x0,L2_T,L2_I,[],[],options); % least sq  result

        time = L2_T(1)+1e-7 :1e-7:L2_T(end); % to make the least sq current more smooth and for close fitting

        Lst_sq_I=F(output,time );

You can suggest any other method for more precise result.!

Answer 1

[Matt J]

The significance of the initial guess x0 is to help lsqcurvefit decide where to start its iterative search for the minimum of the least squares cost function. The cost function is non-convex in your case, so choosing a good initial approximate x0 of the final solution can be essential to helping the iterations avoid local minima.

Your code looks fine, except that I don't know how you decided on your initial guess x0=[5 1 pi/3 ]. If you know a priori that these values are pretty close to what the final values should be, then fine. But, for the curve equation you've shown, it should be possible to derive an initial guess more systematically, e.g,

     avg=mean(ydata);
     peak=max(ydata-avg);
     phase=asin((ydata(i0)-avg)./peak); %where i0 corresponds to time t=0

     x0=[avg, peak, phase];

[Matt J]

It would also be worth trying fminspleas ( Download ). It can take advantage of the linear dependence of the curve on x(1) and x(2). Below, the variable "phase" is your initial guess for the phase, generated using my suggestion above or however else you want.

      flist={1, @(x3,t) sin(2*pi*100*t+x3)};

      [x3,x1x2]=fminspleas(flist,phase,L2_T,L2_I);

      x(1)=x1x2(1);
      x(2)=x1x2(2);
      x(3)=x3;

[Vipin Mohan]

Hi Matt,

Thanks!!

I chose x0=[5 1 pi/3 ] just as a guess. Not a reliable guess. I forget to mention it in the question. Ur explanation is good..!! But still I have doubts on my last 2 3 questions.

1. Does the 'phase' from ur code "x0=[avg, peak, phase];" give the phase shift of current from voltage ?

2. Does avg+peak give the max point ?

Thank you.!

[Matt J]

1. The "phase" from my code is an attempted initial estimate of the unknown x(3) in the curve equation that you provided x(1)+x(2)*sin(2*pi*100*t+x(3)). Since you provided the equation, you should know what it means physically!

2. It is an initial approximation of x(1)+x(2). The value x(1)+x(2) would definitely be the max value of the curve x(1)+x(2)*sin(2*pi*100*t+x(3)). Looking at the plot you recently added to your post, however, I'm not sure you are modeling the data with the correct equation. It looks like the data has a linearly growing term and would be better modeled by,

   x(1)+x(2)*t+x(3)*sin(2*pi*100*t+x(4))

[Vipin Mohan]

Thanx Matt, In ur previous comment

 phase=asin((ydata(i0)-avg)./peak);

Whats it physical meaning.

Inverse sine gives the arc length. but whats the use dividing the first mean value with peak !! whats it actually meant!!

[Matt J]

As you can see from

x0=[avg, peak, phase];

the variables avg, peak, and phase are estimates of x(1),x(2),x(3) respectively. The estimate of the phase comes from writing your original model equation at t=0

   ydata(t==0)=x(1)+x(2)*sin(x(3))

and solving for x(3),

   x(3)=asin((ydata(t==0)-x(1))/x(2))

Answer 2

[Sami Gernaz]

hi Mohan, I'm trying to run your script in Matlab but it is not working any help plz? still new to matlab. thank you