Strang's Linear Algebra: homework

In Strang, section 4.2 we're doing projections. Example A (p. 213) asks:

Project the vector b = (3,4,4) onto the line through a = (2,2,1), and then onto the plane that also contains a* = (1,0,0).
We're supposed to find p, the projection onto a, which is equal to x̂ (x with a little hat on it) times a. It is some fraction of a. The "error" is the part of b that is perpendicular to the projection:

x̂ a + e = b p + e = b e = b - p

The basic equations are:

x̂ = a • b / a • a (for vectors) x̂ = AT b / AT A (for matrices) AT A x̂ = AT b

This comes from

e = b - x̂ a a • e = 0 a • (b - p) = 0 a • b - a • x̂ a = 0

Furthermore

p = x̂ a P = A (AT A)-1 AT

For the first part, with a = (2,2,1), we just have

x̂ = (2,2,1)•(3,4,4) / (2,2,1)•(2,2,1) = 18/9 = 2 p = x̂ a = (4,4,2) e = b - p = (3,4,4) - (4,4,2) = (-1,0,2)

We can check that

e • a (or e • x̂a) = 0

---

Part 2

Now we also consider a* = (1,0,0) to give a plane formed from the 2 vectors.
We construct the matrix:

A = [a a*] = [ 2 1 ] [ 2 0 ] [ 1 0 ]

The fundamental equation is:

AT A x̂ = AT b

We could solve for x̂, then for p = x̂ a, then for P b = p.
Instead of doing this, Strang uses the equations:

P = A (AT A)-1 AT p* = P b e* = b - p*

I don't want to do the arithmetic, so we'll use numpy:

>>> import numpy as np >>> a = np.array([2,2,1]) >>> b = np.array([3,4,4]) >>> a2 = np.array([1,0,0]) >>> a.shape = (3,1) >>> b.shape = a2.shape = a.shape >>> A = np.hstack([a,a2]) >>> print A [[2 1] [2 0] [1 0]] >>> >>> Atrans = A.transpose() >>> temp = np.linalg.inv(np.dot(Atrans,A)) >>> P = np.dot(A,np.dot(temp,Atrans)) >>> P array([[ 1. , 0. , 0. ], [ 0. , 0.8, 0.4], [ 0. , 0.4, 0.2]]) >>> >>> p2 = np.dot(P,b) >>> p2 array([[ 3. ], [ 4.8], [ 2.4]]) >>> >>> e2 = b - p2 >>> e2 array([[ 0. ], [-0.8], [ 1.6]]) >>>

We confirm that e* (e2 in the code) is perpendicular to both a and a*.
Also, for reasons that I don't understand yet, P2 = P.
[UPDATE: The reason is simple. Suppose we do P b to get the projection of b in the plane. What happens if we do it again? Answer: nothing, we're already in the plane! So P b must equal P P b.

>>> np.dot(e2.transpose(),a) array([[ 2.66453526e-15]]) >>> np.dot(e2.transpose(),a2) array([[ 0.]]) >>> >>> np.dot(P,P) array([[ 1. , 0. , 0. ], [ 0. , 0.8, 0.4], [ 0. , 0.4, 0.2]]) >>>