MATLAB Examples

Robust Fitting

This example shows how to compare the effects of excluding outliers and robust fitting. The example shows how to exclude outliers at an arbitrary distance greater than 1.5 standard deviations from the model. The steps then compare removing outliers with specifying a robust fit which gives lower weight to outliers.

Create a baseline sinusoidal signal:

xdata = (0:0.1:2*pi)';
y0 = sin(xdata);

Add noise to the signal with nonconstant variance.

% Response-dependent Gaussian noise
gnoise = y0.*randn(size(y0));

% Salt-and-pepper noise
spnoise = zeros(size(y0));
p = randperm(length(y0));
sppoints = p(1:round(length(p)/5));
spnoise(sppoints) = 5*sign(y0(sppoints));

ydata = y0 + gnoise + spnoise;

Fit the noisy data with a baseline sinusoidal model, and specify 3 output arguments to get fitting information including residuals.

f = fittype('a*sin(b*x)');
[fit1,gof,fitinfo] = fit(xdata,ydata,f,'StartPoint',[1 1]);

Examine the information in the fitinfo structure.

fitinfo
fitinfo = 

  struct with fields:

           numobs: 63
         numparam: 2
        residuals: [63x1 double]
         Jacobian: [63x2 double]
         exitflag: 3
    firstorderopt: 0.0883
       iterations: 5
        funcCount: 18
     cgiterations: 0
        algorithm: 'trust-region-reflective'
         stepsize: 4.1539e-04
          message: 'Success, but fitting stopped because change in residuals less than tolerance (TolFun).'

Get the residuals from the fitinfo structure.

residuals = fitinfo.residuals;

Identify "outliers" as points at an arbitrary distance greater than 1.5 standard deviations from the baseline model, and refit the data with the outliers excluded.

I = abs( residuals) > 1.5 * std( residuals );
outliers = excludedata(xdata,ydata,'indices',I);

fit2 = fit(xdata,ydata,f,'StartPoint',[1 1],...
           'Exclude',outliers);

Compare the effect of excluding the outliers with the effect of giving them lower bisquare weight in a robust fit.

fit3 = fit(xdata,ydata,f,'StartPoint',[1 1],'Robust','on');

Plot the data, the outliers, and the results of the fits. Specify an informative legend.

plot(fit1,'r-',xdata,ydata,'k.',outliers,'m*')
hold on
plot(fit2,'c--')
plot(fit3,'b:')
xlim([0 2*pi])
legend( 'Data', 'Data excluded from second fit', 'Original fit',...
    'Fit with points excluded', 'Robust fit' )
hold off

Plot the residuals for the two fits considering outliers:

figure
plot(fit2,xdata,ydata,'co','residuals')
hold on
plot(fit3,xdata,ydata,'bx','residuals')
hold off