Thread Subject:

lsqcurvefit gives different results with different formulations of the same function

Hi,

I call lsqcurvefit to fit some time data to a sum of cosine functions, with an exponential decay and a constant offset:

     input=[1 0.061 1 0.75 1 0.78 1 0.69 1 0.72 1 0.03 1 0.81 0 0];
     [coeffs, chisq, residual, exitflag]...
         = lsqcurvefit(@cosinesum7, input, time, A(peak,:), [], [], options);

where 'input' is the list of parameters for the function (cosine amplitudes etc., which are a good initial guess), 'time' is the x variable and 'A(peak,:)' is the y variable.

The user supplied function 'cosinesum7', which works fine, looks like this:

     function [y]=cosinesum7(input,x)

     y=(input(1).*cos((input(2)).*(x))+input(3).*cos((input(4)).*(x))...
         +input(5).*cos((input(6)).*(x))+input(7).*cos((input(8)).*(x))...
         +input(9).*cos((input(10)).*(x))+input(11).*cos((input(12)).*(x))...
         +input(13).*cos((input(14)).*(x)))...
         .*exp(-input(15).*(x))+input(16);

I wrote an alternative formulation, which is more flexible as you can use any number of cosine functions:

     function [y]=cosinesum_any(input,x)

     fre=zeros(size(x));
     for n=1:2:length(input)-3
          fre=fre+(input(n).*cos(input(n+1).*x));
     end
     pre=exp(-input(length(input)-1).*(x))+input(length(input));
     y=fre.*pre;

Now I know that this is an equivalent formulation, as if you put the same set of input parameters and time data into both functions, you get the same result! But lsqcurvefit behaves differently using the two functions, with the alternative formulation not fitting well to the data at all (many more iterations occur, depending on the tolerances, and the calculated adjusted R^2 is about 0.5 compared to 0.95 when fitting with the first function).

I do not understand how lsqcurve fit can behave differently in the two situations - as I understand it, it compares the output of the function to the raw data, and adjusts parameters to see if improvements can be made. If the output of the functions is the same, it should be treated the same, regardless of the inner workings of the function. lsqcurvefit does not see these, only the output.

Is this not the case? Does lsqcurvefit analyse the user supplied function in some way that is not compatible with the formulation of the second function? If so, any ideas on how I can formulate the function to work with lsqcurvefit and still be able to use any number of cosine functions?

Thanks very much with this! I can do it the first way, but it seems so much more elegant to do it the second and I'm curious as to why it doesn't work.

My guess is that it has something to do with input(16) being added to
everything else you're computing in cosinesum7, while in cosinesum_any
input(length(input)) is added to the exponential then multiplied by
everything else.

This doesn't show up for your current input vector since input(16) is 0.
But if you make it something else (which is what lsqcurvefit is doing to try
and find a better fit), then the results won't match up.




"Jonathan Midgley" <jonmidgley@hotmail.com> wrote in message
news:k2d17o$p0c$1@newscl01ah.mathworks.com...
> Hi,
>
> I call lsqcurvefit to fit some time data to a sum of cosine functions,
> with an exponential decay and a constant offset:
>
> input=[1 0.061 1 0.75 1 0.78 1 0.69 1 0.72 1 0.03 1 0.81 0 0];
> [coeffs, chisq, residual, exitflag]...
> = lsqcurvefit(@cosinesum7, input, time, A(peak,:), [], [],
> options);
>
> where 'input' is the list of parameters for the function (cosine
> amplitudes etc., which are a good initial guess), 'time' is the x variable
> and 'A(peak,:)' is the y variable.
>
> The user supplied function 'cosinesum7', which works fine, looks like
> this:
>
> function [y]=cosinesum7(input,x)
>
> y=(input(1).*cos((input(2)).*(x))+input(3).*cos((input(4)).*(x))...
> +input(5).*cos((input(6)).*(x))+input(7).*cos((input(8)).*(x))...
>
> +input(9).*cos((input(10)).*(x))+input(11).*cos((input(12)).*(x))...
> +input(13).*cos((input(14)).*(x)))...
> .*exp(-input(15).*(x))+input(16);
>
> I wrote an alternative formulation, which is more flexible as you can use
> any number of cosine functions:
>
> function [y]=cosinesum_any(input,x)
>
> fre=zeros(size(x));
> for n=1:2:length(input)-3 fre=fre+(input(n).*cos(input(n+1).*x));
> end
> pre=exp(-input(length(input)-1).*(x))+input(length(input));
> y=fre.*pre;
>
> Now I know that this is an equivalent formulation, as if you put the same
> set of input parameters and time data into both functions, you get the
> same result! But lsqcurvefit behaves differently using the two functions,
> with the alternative formulation not fitting well to the data at all (many
> more iterations occur, depending on the tolerances, and the calculated
> adjusted R^2 is about 0.5 compared to 0.95 when fitting with the first
> function).
>
> I do not understand how lsqcurve fit can behave differently in the two
> situations - as I understand it, it compares the output of the function to
> the raw data, and adjusts parameters to see if improvements can be made.
> If the output of the functions is the same, it should be treated the same,
> regardless of the inner workings of the function. lsqcurvefit does not see
> these, only the output.
>
> Is this not the case? Does lsqcurvefit analyse the user supplied function
> in some way that is not compatible with the formulation of the second
> function? If so, any ideas on how I can formulate the function to work
> with lsqcurvefit and still be able to use any number of cosine functions?
>
> Thanks very much with this! I can do it the first way, but it seems so
> much more elegant to do it the second and I'm curious as to why it doesn't
> work.

"Jonathan Midgley" <jonmidgley@hotmail.com> wrote in message
> Now I know that this is an equivalent formulation, as if you put the same set of input parameters and time data into both functions, you get the same result!

They are not. You should debug more carefully.

Bruno

Thanks for taking the time to point out my silly mistake!

I'll be more careful in future

"Seth Deland" <sdeland@mathworks.com> wrote in message <k2d2i2$14b$1@newscl01ah.mathworks.com>...
> My guess is that it has something to do with input(16) being added to
> everything else you're computing in cosinesum7, while in cosinesum_any
> input(length(input)) is added to the exponential then multiplied by
> everything else.
>
> This doesn't show up for your current input vector since input(16) is 0.
> But if you make it something else (which is what lsqcurvefit is doing to try
> and find a better fit), then the results won't match up.
>
>
>
>
> "Jonathan Midgley" <jonmidgley@hotmail.com> wrote in message
> news:k2d17o$p0c$1@newscl01ah.mathworks.com...
> > Hi,
> >
> > I call lsqcurvefit to fit some time data to a sum of cosine functions,
> > with an exponential decay and a constant offset:
> >
> > input=[1 0.061 1 0.75 1 0.78 1 0.69 1 0.72 1 0.03 1 0.81 0 0];
> > [coeffs, chisq, residual, exitflag]...
> > = lsqcurvefit(@cosinesum7, input, time, A(peak,:), [], [],
> > options);
> >
> > where 'input' is the list of parameters for the function (cosine
> > amplitudes etc., which are a good initial guess), 'time' is the x variable
> > and 'A(peak,:)' is the y variable.
> >
> > The user supplied function 'cosinesum7', which works fine, looks like
> > this:
> >
> > function [y]=cosinesum7(input,x)
> >
> > y=(input(1).*cos((input(2)).*(x))+input(3).*cos((input(4)).*(x))...
> > +input(5).*cos((input(6)).*(x))+input(7).*cos((input(8)).*(x))...
> >
> > +input(9).*cos((input(10)).*(x))+input(11).*cos((input(12)).*(x))...
> > +input(13).*cos((input(14)).*(x)))...
> > .*exp(-input(15).*(x))+input(16);
> >
> > I wrote an alternative formulation, which is more flexible as you can use
> > any number of cosine functions:
> >
> > function [y]=cosinesum_any(input,x)
> >
> > fre=zeros(size(x));
> > for n=1:2:length(input)-3 fre=fre+(input(n).*cos(input(n+1).*x));
> > end
> > pre=exp(-input(length(input)-1).*(x))+input(length(input));
> > y=fre.*pre;
> >
> > Now I know that this is an equivalent formulation, as if you put the same
> > set of input parameters and time data into both functions, you get the
> > same result! But lsqcurvefit behaves differently using the two functions,
> > with the alternative formulation not fitting well to the data at all (many
> > more iterations occur, depending on the tolerances, and the calculated
> > adjusted R^2 is about 0.5 compared to 0.95 when fitting with the first
> > function).
> >
> > I do not understand how lsqcurve fit can behave differently in the two
> > situations - as I understand it, it compares the output of the function to
> > the raw data, and adjusts parameters to see if improvements can be made.
> > If the output of the functions is the same, it should be treated the same,
> > regardless of the inner workings of the function. lsqcurvefit does not see
> > these, only the output.
> >
> > Is this not the case? Does lsqcurvefit analyse the user supplied function
> > in some way that is not compatible with the formulation of the second
> > function? If so, any ideas on how I can formulate the function to work
> > with lsqcurvefit and still be able to use any number of cosine functions?
> >
> > Thanks very much with this! I can do it the first way, but it seems so
> > much more elegant to do it the second and I'm curious as to why it doesn't
> > work.