Solving closest point in the span of many vectors
Goal:
Special case: We can use the algorithm to determine whether b lies in Span {v1, . . . , vn}: If the vector in Span {v1, . . . , vn} closest to b is b itself then clearly b is in the span; if not, then b is not in the span.
Testing if b is in Span {v1, . . . , vn} is equivalent to testing if the equation Ax = b has a solution.
More generally: Even if Ax = b has no solution, we can use the algorithm to find the point in <math>\{Ax : x \in \mathbb{R}\}</math> closest to b.
Goal: find the coefficients x1, . . . , xn so that x1 v1 + · · · + xn vn is the vector in Span {v1, . . . , vn} closest to b.
Two methods:
Equivalent: Find the vector x that minimizes:
Given a matrix equation Ax = b that might have no solution, find the best solution available in the sense that the norm of the error b - Ax is as small as possible.
Two algo:
Find the vector x that minimizes:
This problem is called least squares (”methode des moindres carres”, due to Adrien-Marie Legendre but often attributed to Gauss).
High-dimensional projection onto/orthogonal to
For a vector b and a vector space V, we define the projection of b onto V (written <math>b^{||V}</math> ) and the projection of b orthogonal to V (written <math>b^{\perp V}</math> ) so that <MATH> b = b^{||V} + b^{\perp V} </MATH> where:
High-Dimensional Lemma: The point in vector space V closest to b is <math>b^{||V}</math> and the distance is the norm <math>\|b^{\perp V}\|</math> .
The below function will find the projection of b orthogonal to the space spanned by mutually orthogonal vectors.
# input: a vector b, and a list vlist of mutually orthogonal vectors
# output: the projection orthogonal of b (\perp) to the vectors in vlist
def project_orthogonal(b, vlist):
for v in vlist:
b = b - project_along(b, v)
return b
where:
The projection of a vector b orthogonal to a vector space is unique, so in principle the order of vectors in vlist doesn’t affect the output of project orthogonal(b, vlist).
The below function will return the <math>\sigma</math> of all generators.
def aug_project_orthogonal(b, vlist):
alphadict = {len(vlist):1}
for i in range(len(vlist)):
v = vlist[i]
sigma = (b*v)/(v*v) if v*v > 0 else 0
alphadict[i] = sigma
b = b - sigma*v
return (b, alphadict)
Exactly the same computations take place when orthogonalize is applied. Ie applying orthogonalize to the list of vector <math>[{v_1} , \dots , {v_n}, b]</math> will return <math>[v_1^* , \dots , v_n^*, b^{\perp}]</math>