## Let's see if my soln is the same
function least_squares_transform(A)
T = Matrix(1.0I,size(A,1),size(A,1))
Z = Matrix(1.0I,size(A,1),size(A,1))
R = Matrix(1.0I,size(A,1),size(A,1))
Q = Matrix(0.0I,size(A,1),size(A,2))
for i=1:size(A,2)
a = A[:,i]
na = norm(a)
a = a/na
Q[:,i] = a
R = R*(Matrix(1.0I,size(A,1),size(A,1)) - a*a')
T = (Matrix(1.0I,size(A,1),size(A,1)) - a*a')*T
A = (Matrix(1.0I,size(A,1),size(A,1)) - a*a')*A
Z = (Matrix(1.0I,size(A,1),size(A,1)) - 2*a*a')*Z
end
return T,R,Q,Z
end
T,R,Q,Z = least_squares_transform(A)
display(T*b)