Curve fitting: Getting very poor results on seemingly simple data set ("fit" function)

28 views (last 30 days)
Hello everyone,
I will be the first to admit I am a complete newbie with the "fit" function and "cfun" data type, and find the documentation available on it much less than user friendly, so please refrain from directing me to it - I have been through it front to back and have not found a solution.
I am going through a dataset in which I am analyzing hundreds of discrete blocks, each of which has an essentially periodic behaviour. I need to find the amplitude of this periodic behaviour. I have included an example.
My problem is this. I used a very clear sample dataset containing 13 samples (one of the best I have) with a very obvious cosine behaviour. I set the fit function up and ran it, and got parameters that give me almost a straight line. I tried using repmat thinking maybe it was my sample size, and this does almost nothing.
I manually solved the periodic fit and have included the plot with the original data, the empirically derived fit, and the results from the 'fit' function below. I have also provided what I believe to be the relevant chunk of code, but I could be wrong. Please let me know what else if anything I need to include. As I have to perform this operation on literally hundreds or thousands of these blocks, this is not something I can feasibly do manually, so if someone can direct me as to what I'm doing wrong, I would be infinitely grateful.
Thank you kindly for your time.
My image should be available here: http://i42.tinypic.com/jf7oud.jpg
Here is what cfun spits out:
cfun =
General model Sin1:
cfun(x) = a1*sin(b1*x+c1)
Coefficients (with 95% confidence bounds):
a1 = 0.8323 (-22.99, 24.66)
b1 = 0.004459 (-0.9294, 0.9383)
c1 = 1.936 (-72.97, 76.84)
And here is the relevant snippet:
%%%%%%%%%%%%%%%%%
% CURVE FITTING %
%%%%%%%%%%%%%%%%%
x = 2*pi.*((0:num_phases-1)./(num_phases-1));
num_phases2 = num_phases*10;
x2 = 2*pi.*((0:num_phases2-1)./(num_phases2-1));
x = x';
x2 = x2';
y = results(working_block*num_phases+1:working_block*num_phases+num_phases,1);
y2 = repmat(y,10,1);
cfun = fit(x,y,'sin1');
test = 0.025*cos(2*pi.*((0:num_phases-1)./(num_phases-1)))+0.775;
figure(1);
plot(2*pi.*((0:num_phases-1)./(num_phases-1)),results(working_block*num_phases+1:working_block*num_phases+num_phases,1)); hold on;
plot(2*pi.*((0:num_phases-1)./(num_phases-1)),test,'-r');
plot(cfun,'-g');
xlabel('Block');
ylabel('Correlation');
legend('Original Data','0.025cos(2pi*x)+0.775','Fitted Data');

Accepted Answer

dpb
dpb on 18 Aug 2013
Edited: dpb on 18 Aug 2013
Your model without an intercept term doesn't fit the data pattern worth a hoot...
With no intercept there's no offset to shift an initial zero value of sin(0) and the amplitude of the sin/cos has to be about +/-0.025 to fit the range of the data.
Try something more like
g = fittype( @(a,b,c,d,x) a.*sin(b.*x+c)+d )
With that I get a reasonable fit. Of course, w/ four parameters and only 13 data points you're using up a sizable fraction of the available DOF...
ADDENDUM:
Came to see if had been any followup and noticed what hadn't actually seen before...
I see in your legend text you don't have a phase angle in the equation--altho you don't show it, in that case you mis-specified the model. As the tool showed you, the model you fitted was
cfun(x) = a1*sin(b1*x+c1)
whereas you wanted (apparently)
cfun(x) = a1*sin(b1*x)+c1
based on the text in
legend('Original Data','0.025cos(2pi*x)+0.775','Fitted Data');
IOW, you had a misplaced ')' in the function definition...
  2 Comments
Matthew
Matthew on 19 Aug 2013
Edited: Matthew on 19 Aug 2013
Thank you, that definitely makes much more sense. However, 0.025*cos(2*pi.*((0:num_phases-1)./(num_phases-1)))+0.775; will not always be the appropriate parameters. I assume that @(a,b,c,d,x) are function handlers to find the appropriate values?
Additionally, when I try that snippet, I get this error:
??? Error using ==> fittype.fittype>fittype.fittype at 201
First input argument must be a string or cell array.
Error in ==> Results_Analysis_working_copy at 95
g = fittype( @(a,b,c,d,x) 'a.*sin(b.*x+c)+d' )
Could you provide an extended snippet of how you utilized this code to get a reasonable fit?
Matthew
Matthew on 19 Aug 2013
Edited: Matthew on 19 Aug 2013
My ignorance, it works now.
Per your ADDENDUM, that I was using the 'sin1' model, which would be a.*sin(b.*x+c) model. I tried using the 'a.*sin(b.*x+c)+d' model as well, and it works. Thank you kindly.

Sign in to comment.

More Answers (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!