x=0:.05:1;
y= [.486 .866 .944 1.144 1.103 1.202 1.166 1.191 1.124 1.095 1.122 1.102 1.099 1.017 1.111 1.117 1.152 1.265 1.38 1.575 1.857]';

plot(x,y,'x')
title('Data for cubic Least Squares example');

%input(''); % pause for 

z= [ones(size(x)); x; x.^2; x.^3]';

%coefficient matrix- symmetric
A= z'*z
b= z'*y

%find the solution
a= A\b % = (z'*z)\(z'*y)
% a =
%     0.5747
%     4.7259
%   -11.1282
%     7.6687

hold on
plot (x, a(4)*x.^3 + a(3)*x.^2 + a(2)*x + a(1))

% check the error for least squares- find E_r
er= sum((y-z*a).^2)

% the coefficient of determination
r2= 1-er/sum((y-mean(y)).^2)
% therefore 97.2% of the original uncertainty has been explained by the
% nonlinear model

% but is our function well-conditioned?
cond= norm(z'*z, inf) * norm(inv(z'*z),inf)
% no- 2.1981e+04 == 21,981
log10(cond)
% therefore, expect 4 digits to be off 
% means that a small change in b will make a large change to the solution

% try it & see
b2= b+[0.01 -.01 .01 -.01]'
a2= A\b2
% a2 =
%    0.74078985507247
%    2.68251884101982
%   -6.15375923339264
%    4.45500003852375
% small change in b had a large change in our solution