least_squares.m

function [a,fx,r2] = least_squares(Z,y)
%This function accepts a Z-matrix and dependent variable y
%as an input. This function returns the vector of best-fit 
%constants, a, the best-fit function evaluated at each
%point, fx, and the coefficient of determination, r2. 
n=length(Z);
Z=Z(:); 
y=y(:); %convert to column vectors
sZ=sum(Z); sy=sum(y);
sZ2=sum(Z.*Z); sZy=sum(Z.*y); sy2=sum(y.*y);
a(1)=(n*sZy-sZ*sy)/(n*sZ2-sZ^2);
a(2)=sy/n-a(1)*sZ/n;
r2=((n*sZy-sZ*sy)/sqrt(n*sZ2-sZ^2)/sqrt(n*sy2-sy^2))^2;
fx=Z*a;
end