Approximation by polynomials using the least squares method

Using Matlab function polyfit

>> X =  [-1 0 0 1 1 2 4];  % x-coordinates of the points
>> Y =  [ 5 6 5 7 6 8 11]; % y-coordinates of the points
>> step = 0.2;             % step-size of the graph of the polynomial
>> xx = X(1):step:X(end);  % division on x-axis 

>> deg = 3;                % degree of the polynomial
>> p = polyfit(X, Y, deg)  % computation of coeff. of the polynomial

>> yy = polyval(p, xx);    % corresponding values of the polynomial
>> plot(X, Y, 'ro', xx, yy);

>> Yp = polyval(p, X);  % approximate values at given points
>> D = abs(Y-Yp);       % distances from given values
>> delta = sqrt( sum( D.^2) ) % quadratic deviation
>> max_err = max(D)     % maximal error 

Experiment with changing the degree of the polynomial from constant (i.e. average value, degree 0), a line (degree 1), up to interpolation polynomial (it can be obtained for distinct x-coordinates only). If you change degree of the polynomial only, you can skip first 4 commands.


Using a projection

You can use basics of Matlab instead of the function polyfit (the first five and the last six commands are the same as before):

>> X =  [-1 0 0 1 1 2 4];  % x-coordinates of the points
>> Y =  [ 5 6 5 7 6 8 11]; % y-coordinates of the points
>> step = 0.2;             % step-size of the graph of the polynomial
>> xx = X(1):step:X(end);  % division on x-axis 

>> deg = 3;                % degree of the polynomial

>>   %%% computation of coeff. of the polynomial:
>> n = length(X); 
>> A = ones(n,1);        % the first column - ones
>> for k = 1:deg
>>   A = [A , (X.^k)' ]; % for every power k, add one column
>> end
>>   %%%  A' realizes projection to lin. space
>>   %%%  generated by columns of A:
>> C = A' * A;          
>> b = A' * Y'; 
>> q = C\b;              % computation of coeff. of the polynomial
>> p = fliprl(q)         % reversing the elements in q        

>> yy = polyval(p, xx);  % corresponding values of the polynomial
>> plot(X, Y, 'ro', xx, yy);

>> Yp = polyval(p, X);  % approximate values at given points
>> D = abs(Y-Yp);       % distances from given values
>> delta = sqrt( sum( D.^2) ) % quadratic deviation
>> max_err = max(D)     % maximal error