Code covered by the BSD License  

Highlights from
Data Driven Fitting with MATLAB

image thumbnail

Data Driven Fitting with MATLAB

by

 

Demostration code and data sets for the "Data Driven Fitting with MATLAB" webinar.

FititDemo.m
%% Demo for "Fitit"
% Copyright (c) 2011, The MathWorks, Inc.

%% Generate a set of random data

clear all
clc

s = RandStream('mt19937ar','seed',1971);
RandStream.setDefaultStream(s);

X = linspace(1,10,100);
X = X';

% Specify the parameters for a second order Fourier series

w = .6067;
a0 = 1.6345;
a1 = -.6235;
b1 = -1.3501;
a2 = -1.1622;
b2 = -.9443;

% Fourier2 is the true (unknown) relationship between X and Y
Y = a0 + a1*cos(X*w) + b1*sin(X*w) + a2*cos(2*X*w) + b2*sin(2*X*w);

% Add in a noise vector

K = max(Y) - min(Y);
noisy = Y +  .2*K*randn(100,1);

%% Generate a set of random data

s = RandStream('mt19937ar','seed',1973);
RandStream.setDefaultStream(s);

hold off

X1 = linspace(1,10,100);
X1 = X1';

% Specify the parameters for a second order Fourier series

w = rand;
a0 = randn;
a1 = randn;
b1 = randn;
a2 = randn;
b2 = randn;

% Fourier2 is the true (unknown) relationship between X and Y
Y1 = a0 + a1*cos(X*w) + b1*sin(X*w) + a2*cos(2*X*w) + b2*sin(2*X*w);

% Add in a noise vector

K = max(Y1) - min(Y1);
noisy1 = Y1 +  .2*K*randn(100,1);

foo = fitit(X1,noisy1)

% Compare the output of "fitit" with the original model

hold on
plot(X1,Y1, 'r', 'linewidth', 2)

%%  Working with the fit object

figure
plot(My_Fit)
My_Fit(noisy(1:10))

figure;
subplot(2,1,1); plot(X,My_Fit(X)); title('My Fit');
subplot(2,1,2); plot(X, differentiate(My_Fit, X)); refline(0,0),title('First derivative of My Fit');

methods(My_Fit)

%%  Add confidence intervals into the mix

[avgfit,fitstd1,fitstd2] = fitit(X,noisy,1000);

%%  Explain Localized Regression

figure
scatter(X, noisy)

fit1 = smooth(X,noisy,0.02,'loess')
hold on
plot(X, fit1, 'r', 'linewidth', 2)

fit2 = smooth(X,noisy,0.95,'loess')
plot(X, fit2, 'k', 'linewidth', 2)

fit3 = smooth(X,noisy,0.2,'loess')
plot(X, fit3, 'b', 'linewidth', 2)

%% Explain Cross Validation

% Divide data set into a test set and a training set

cp = cvpartition(size(noisy,1),'k',10);

trainingset1 = X(training(cp,1));
trainingset2 = noisy(training(cp,1));

testset1 = X(test(cp,1));
testset2 = noisy(test(cp,1));

figure
trainingscatter = scatter(trainingset1,trainingset2,'+', 'r');
hold on
testscatter = scatter(testset1, testset2,'FaceColor','b','Filled');
legend('Training Set','Test Set', 'location', 'NorthWest');

%% Use Training Set to create a LOESS smooth

delete(testscatter);
Z = smooth(trainingset1,trainingset2,.5,'loess');
ZPlot = plot(trainingset1, Z, 'color','r','linestyle','-','linewidth',2);

%% Use the test set to evaluate goodness of fit

foobar = fit(trainingset1, Z, 'cubicinterp')
scatter(testset1, testset2,'FaceColor','b','Filled')

for i = 1: length(testset1)
    
    plot([testset1(i) testset1(i)], [foobar(testset1(i)) testset2(i) ])
    
end

%%  Explain Bootstrap

figure
scatter(X, noisy, 'b')

% Create a new data set
index = randsample(length(X),length(X), 'true');
hold on
scatter(X(index), noisy(index), 'filled', 'r')
bar = fit(X(index), noisy(index), 'fourier2');
plot(bar)
legend('Original Data Set', 'Random Sample', 'Fit', 'location', 'NorthWest')


%%  Repeat

index = randsample(length(X),length(X), 'true');
scatter(X(index), noisy(index), 'filled', 'k')
bar = fit(X(index), noisy(index), 'fourier2');
plot(bar, 'k')
legend('off')

%% Repeat 20 more times

for i = 1:20
    
    index = randsample(length(X),length(X), 'true');
    bar = fit(X(index), noisy(index), 'fourier2');
    plot(bar, 'b')

end

legend off

Contact us