Linear Algebra II: Project 1

Megan Wheeler

March 20, 2017

1 Abstract
The goal of this project was to find the best approximation to the least squares problem, min ||Ax b||2 , by
xR3
computing the QR-decomposition of the matrix A using Householder reflectors and Givens rotations.

2 Mathematical Background
2.1 Least Squares Approximation
Consider the equation Ax = b where A is an m n matrix. Then x must be a vector in Rn and b must be
a vector in Rm . If b is not in the column space (range) of A (b 6 Rn ), then there is no solution to this equation
then . A visual representation of this can be seen in Figure 1 below. When there is no solution to the equation
Ax = b, we can still find a best approximate solution, x where A x is in the range of A and is as close to b as
possible. Therefore, we want to find a vector x
that minimizes the length between the vector b and Ax. This will
be the solution to our least squares problem, min ||Ax b||2 . The closest vector to b that is in the range of A is
xR3
the projection of b onto the column space of A. Therefore, A
x = projRn b (seen in Figure 2).

Figure 1: In this image, b is not in the range of A Figure 2: A

x = projRn b.A x b is in the orthogonal
compliment of the column space of A

x b R(A) or A
We can see that A x b N (AT ). Thus
AT (A
x b) = 0
AT A
x AT b = 0
AT A
x = AT b.
If A has full column rank (the columns of A are linearly independent rankA = n), then AT A is invertible.
Therefore,
x = (AT A)1 AT b = A b. (1)
(AT A)1 AT is called the Moore-Penrose Pseudoinverse and is represented by the symbol A . This can be claculated
analytically but in practice, the inverse of a matrix can become very tedious to calculate for larger matrices (see
Figures 3 below which is only the inverse of a 4 4 matrix).

1
 Figure 3: The inverse of a 4 4 matrix. [2]

2.2 QR Decomposition
We can however, simplify our least squares problem by reducing our matrix A to upper trapezoidal form. If
A has full column rank, then we can factor A to,

Amn = Qmm Rmn

where Q is an orthogonal matrix (QT Q = I) and R is our upper trapezoidal matrix. Our least squares approxi-
mation then becomes:
QRx b (2)
Q can be written as:  
Q = q1 q2 ... qn qn+1 ... qm
where q1 , q2 , ...qn is an orthonormal basis for R(A) and qn+1 , ...qm is an orthonormal basis for R(A) . Q can
 

therefore be decomposed as,

h i
m(mn) mn m(mn)
Qmn
   
k Q where Q k = q 1 q 2 ... q n and Q = q n+1 ... q m

Since R is upper trapezoidal, we can write it as:

 nn 
R
R = (mn)n (see Figure 4 below for a clarifying example.)
0

Figure 4: A 4 3 trapezoidal matrix represented as a 3 3 triangular matrix and a 1 3 zero matrix.

2
 Our vector b can be expressed as:

b = bk + b , bk R(A) and b R(A)

Qk is a projection onto R(A) and Q is a projection onto R(A) , so b can be decompsed as:

bk = Qk QTk b

b = Q QT b
Plugging all of this into Eq:(2), we get:
h i  Rnn  h i QT b
mn m(mn) mn m(mn) k
Qk Q x Qk Q (3)
0(mn)n QT b

Because Q is an orthonomal matrix, if we multiply both sides of Eq:(3) by QT , we are left with:
 nn 
QTk b
 
R
x (4)
0(mn)n QT b

We have now successfully decoupled our system into two subsystems of linear equations:
nn
(1) R x = QTk b

(2) 0(mn)n x = QT b
The second equation cant be solved, but the first equation gives us the solution to our least squares problem:
1 T
x
= R Qk b (5)

Pluging Eq:(5) and A = QR into A

x:
1 T
A
x = Qk R R Qk b

x = Qk QTk b
A

x = bk
A (6)

2.3 Householder Reflections

There are a couple ways to compute Q and R to decouple a system into two subsystems, one of which uses
linear transformations called Householder reflectors. Applying these Householder transformations to both sides
of our equation Ax = b will reduce A to an upper trapezoidal matrix, R and decompose b into its orthogonal
components. A Householder reflector, H, has the form:
2
H=I vv T
||v||22

H is a not a projection operator because it is not idempotent (H 2 6= H). It is however and orthogonal projection
because H T = H (H T H = I). If we allow H to operate on any arbitrary vector, x Rn ,

hv, xi
Hv x = v
||v||22

then Hv can be seen as the orthogonal projection of x onto v (Figure 5).

3
 Figure 5: The projection of x onto v is the blue vector.

If Hv is an orthogonal projection operator that projects x onto the subspace of V = span{v}, then the operator
I Hv is a projection operator that projects x onto V (Figure 6).

Figure 6: The projection of x onto the subspace V is the blue vector.

Now if we write H = I Hv Hv , then to compute Hx we subtract the projection of x onto V (the blue vector
in Figure 5) from the projection of x onto V (blue vector in Figure 6). The result of this is the purple vector in
Figure 7. H is called a reflector because Hx reflects any vector x about span{v} .

Figure 7: The purple vector is the projection of x onto the subspace V minus the projection of x onto the
subspace V .

2.4 Givens Rotations

Another way to compute the QR factorization is by Givens rotations. Givens matrices, like householder
reflectors, are orthogonal matrices. The general form of a Givens rotation matrix is given by:

I 0 0 0 0
0 cos 0 sin 0

Gi,j,
0 0 I 0 0

0 sin 0 cos 0
0 0 0 0 I

[3] where i is the row which contains cos and sin and j is the row which contains sin and cos . The Givens
operator Gi,j, rotates a matrix in the (i, j) plan clockwise by radians, hence named Givens rotation. The Givens

4
operator Gi,j, will only effect rows i and j of the matrix it acts on. The Givens operators that we will use for our
4 3 matrix are listed below:

cos sin 0 0 1 0 0 0 1 0 0 0
sin cos 0 0 0 cos sin 0 0 1 0 0
G12, = , G23, = , G34, =
0 0 1 0 0 sin cos 0 0 0 cos sin
0 0 0 1 0 0 0 1 0 0 sin cos

cos and sin are chosen such that:

    
cos sin a r p
= , where r = a2 + b2
sin cos b 0

Therefore,
a b
cos = , sin =
r r
From this you can see that the Givens operator pushes the magnitude of [ a b ]T onto the first element in the new
vector. Using Givens matrices you always want to work from the bottom to top and then left to right. (See Figure
8.)

Figure 8: Work from the bottom to top and then left to right when working with Givens matrices. [4]

3 Methods

3 2 12 1
3 1 8 2
3 2 4 ,
A= b=
1

3 1 0 6
Lets first verify that the rank of A is 3. The rank is the maximum number of linearly independent vectors in a
matrix. This is equal to the number of non-zero rows in its row echelon matrix. Putting matrix A into a calculator
that reduces the matrix into row echelon form (see Figure 9), we see that each of the three columns in matrix A
are linearly independent. Therefore, the rank of A is 3 and the matrix has full column rank.

5
 Figure 9: Matrix A in echelon form. [5]

3.1 Using Household Reflectors:

The first household reflector is given by
2
H1 = I v 1 v T1 , where v 1 = a1 ||a1 ||2 e1
||v 1 ||22

Thus, v 1 is given by:

3 1 9
3 p 0 3
v1 = 2 2 2 2
3 + 3 + 3 + 3 + 3 0 = 3

3 0 3
The matrix H1 is then given by:

1 0 0 0 9
0 1 0 0
2 3 
9 3 3 3

H1 =
0 0 1 0 92 + 32 + 32 + 32 3
0 0 0 1 3

1 0 0 0 81 27 27 27
0 1 0 0 1 27
9 9 9
H1 =
0
0 1 0 54 27 9 9 9
0 0 0 1 27 9 9 9

1/2 1/2 1/2 1/2
1/2 5/6 1/6 1/6
H1 =
1/2 1/6 5/6

1/6
1/2 1/6 1/6 5/6
Multiplying A by H1 will put zeros in all rows of the first column following the first row:

1/2 1/2 1/2 1/2 3 2 12
1/2 5/6 1/6 1/6 3 1 8

H1 A = 1/2 1/6 5/6 1/6 3 2 4
1/2 1/6 1/6 5/6 3 1 0

6

6 1 12
0 2 0
H1 A =
0 1 4

0 2 8
And multiplying b by H1 gives us:

1/2 1/2 1/2 1/2 1
1/2 5/6 1/6 1/6 2
H1 b =
1/2 1/6 5/6 1/6 1
1/2 1/6 1/6 5/6 6

4
1/3
H1 b =
8/3

13/3
Next, we want to compute H2 . To do this, we compute the reflector R2 :
2
R2 = I v 2 v T2 , 2 ||
where v 2 = a a2 ||2 e1 .
||v 2 ||22

2
2 = 1
R2 is a 3 x 3 reflector, and a
2
Thus, v 2 is given by:
2 p 1 1
v 2 = 1 + (2)2 + 12 + (2)2 0 = 1
2 0 2
R2 is then given by:
1 0 0 1 
2 1 1 2

R2 = 0 1 0 2
1
1 + 12 + (2)2
0 0 1 2

1 0 0 1 1 2
1
R2 = 0 1 0 1 1 2
3
0 0 1 2 2 4

2/3 1/3 2/3
R2 = 1/3 2/3 2/3
2/3 2/3 1/3
The matrix H2 is then given by:
1 0 0 0
0 2/3 1/3 2/3
H2 =
0 1/3 2/3 2/3
0 2/3 2/3 1/3
Multiplying H1 A by H2 :
1 0 0 0 6 1 12
0 2/3 1/3 2/3 0 2 0
H2 H1 A =
0 1/3 2/3 2/3 0 1 4

0 2/3 2/3 1/3 0 2 8

7

6 1 12
0 3 4
H2 H1 A =
0 0 8

0 0 0
And multiplying H1 b by H2 gives us:

1 0 0 0 4
0 2/3 1/3 2/3 1/3
H2 H1 b =
0 1/3 2/3 2/3 8/3

0 2/3 2/3 1/3 13/3

4
4
H2 H1 b =
1

3
This gives us two equations:
6 1 12 x 1 4
0 3 4 x
2 = 4 (7)
0 0 8 x
3 1
0x
1 + 0 x 3 = 3
2 + 0 x (8)
The second equation has no solution, but the solution to the first equation gives us our least squares solution,

10/9
= 7/6
x
1/8

We can now calculate QT by multiplying the Householder reflectors together:

1/2 1/2 1/2 1/2 1 0 0 0 1/2 1/2 1/2 1/2
1/2 5/6 1/6 1/6 0 2/3 1/3 2/3 = 1/2 1/2 1/2 1/2
QT =

1/2 1/6 5/6 1/6 0 1/3 2/3 2/3 1/2 1/2 1/2 1/2
1/2 1/6 1/6 5/6 0 2/3 2/3 1/3 1/2 1/2 1/2 1/2

We can verify this by multiplying QT by A which should give us R.

We then take the transpose of QT to get Q:

1/2 1/2 1/2 1/2
1/2 1/2 1/2 1/2
(QT )T = Q =


1/2 1/2 1/2 1/2
1/2 1/2 1/2 1/2

3.2 Using Givens Rotations:

Using the G34, Givens rotation matrix to eliminate the 3 in the fourth row of the first column, we must first
determine our values for cos and sin .
p a b
r = a2 + b2 , cos = , sin = , where a = 3 and b = 3.
r r
Therefore,
p 3 2 3 2
r = 32 + 32 = 18 = 3 2, cos 1 = = , sin 1 = =
18 2 18 2

8
Plugging this in to our G34,1 Givens matrix and multiplying by our matrix A:

1 0 0 0 3 2 12 3 2 12
0 1 0 0
3 1 8
3
1 8
A2 = G34,1 A =
=

2 2
2
0 0 3 2 4 3 2 2 2

2
2
2

0 0 2 2 2 3 1 0 0 3 2
2 2 2
2

Now we want to eliminate
the 3 2 in the third row, first column of our new matrix, A2 . Our values for a and b
are now a = 3 and b = 3 2. Therefore,

q 3 3 3 2 6
2 2
r = 3 + (3 2) = 27 = 3 3, cos 2 = = , sin 2 = =
27 3 27 3

Plugging this in to our G23,2 Givens matrix and multiplying by our matrix A2 :

1 0 0 0 3 2 12 3 2 12
3

6 3 1 8 33 0

4 3
0 0

A3 = G23,2 A2 = 3
3 =

2 6
0 36 33 0 3 2 2 2 0 2 6

2
2

0 0 0 1 0 3 2 2 2 2 0 3 2 2 2 2

Next, well use the G12, Givens
matrix to eliminate the 3 3 in the second row, first column of A3 . Our values
for a and b are a = 3 and b = 3 3. Therefore,

q 1 3
2 2
r = 3 + (3 3) = 6, cos 3 = , sin 3 =
2 2
Plugging this in to our G12,3 Givens matrix and multiplying by our matrix A3 :

1 3
0 0 3 2 12 6 1 12
2 2
0 3 43
23 12 0 0 3 3
0 4 3
=
A4 = G12,3 A3 =

3 2 6
0 2 6 0 2 6

0 0 1 0 2 2


0 0 0 1 0 3 2 2 2 2 0 3 2 2 2 2

Now well use again the G34, Givens matrix, this time to eliminate the 3 2 2 in the fourth row, second column of

3 2
A4 Our values for a and b are a = 2 and b = 3 2 2 . Therefore,
s
3 2
3 2 2 3 2 2 1 322 3
r= ( ) + ( ) = 6, cos 4 = 2 = , sin 4 = =
2 2 6 2 6 2

Plugging this in to our G34,4 Givens matrix and multiplying by our matrix A4 :

1 0 0 0 6 1 12

6 1 12
0 1 0 0 0 3 4 3 0 3 4 3
A5 = G34,4 A4 =

1

3
=
0 0 2 2 0 3 22 2 6 0 6 0

0 0 23 1
2 0 322 2 2 0 0 4 2

To eliminate
the 6 in
the third row, second column of A5 we use the G23, Givens matrix. Our values for a and
b are a = 3 and b = 6. Therefore,
q
3 6
r = ( 3)2 + ( 6)2 = 3, cos 5 = , sin 5 =
3 3

9
Plugging this in to our G23,5 Givens matrix and multiplying by our matrix A5 :

1 0 0 0

6 1 12
6 1 12
0 33 6
0 0 3 4 3 = 0 3 4



A6 = G23,5 A5 = 3
0 36 33 0 0 6 0 0 0 4 2

0 0 0 1 0 0 4 2 0 0 4 2

the 4 2 in the
Finally, to get rid of fourth row, third column of A6 we use the G34, Givens matrix. Our values
for a and b are a = 4 2 and b = 4 2. Therefore,
q
2 2
r = (4 2)2 + (4 2)2 = 8, cos 6 = , sin 6 =
2 2
Plugging this in to our G34,6 Givens matrix and multiplying by our matrix A6 :

1 0 0 0

6 1 12 6 1 12
0 1 0 0 0 3 4
= 0 3 4

A7 = G34,6 A6 =

0 22 22 0 0 4 2 0

0 0 8
0 0 2 2 0 0 4 2 0 0 0
2 2

QT can now be calculated by multiplying all of the Givens rotations together: QT = (G34,1 G23,2 G12,3 G34,4 G23,5 G34,6 ).
Therefore,

1 0 0 0 1 0 0 0 1 3 1 0 0 0 1 0 0 0 1 0 0 0
2 2 0 0
0 1 0 0 3 6 0 1 0 0 3 6 0 1 0 0
0 3
3 1
0 0 0

QT =
3
3 2 2 0 0 3

2 2
0 0 6 3
0 1 0 0 0 2
1
23 0 6 3 0 0 0 2
22

2
2
0 3 3 0 0 3 3 2
0 0 22 22 0 0 0 1 0 0 0 1 0 0 23 1
2 0 0 0 1 0 0 2
2
2
2

1/2 1/2 1/2 1/2
1/2 1/2 1/2 1/2
QT =
1/2 1/2 1/2 1/2
1/2 1/2 1/2 1/2
We can verify this by multiplying QT by A which should give us R.

We take the transpose of QT to get Q:

1/2 1/2 1/2 1/2
1/2 1/2 1/2 1/2
(QT )T = Q =
1/2 1/2 1/2 1/2

1/2 1/2 1/2 1/2

Multiplying b by QT :
1/2 1/2 1/2 1/2 1 4
1/2 1/2 1/2 1/2 2 4
1/2 1/2 1/2 1/2 1 = 1

1/2 1/2 1/2 1/2 6 6

We again get two equations:
6 1 12 x
1 4
2 = 4
0 3 4 x (9)
0 0 8 x
3 1
0x
1 + 0 x 3 = 3
2 + 0 x (10)

10
The second equation has no solution, but the solution to the first equation gives us our least squares solution,

10/9
= 7/6
x
1/8

I know this result is not correct because I should get the same least squares solutions for both methods. If my
vector b for the Givens method were < 4, 4, 1 > rather than < 4, 4, 1 > I would get the same solutions.

4 Results
Compute x
using both Householder reflectors and Givens rotations:

10/9
= 7/6
x
1/8

Compute Qk and Q , as well as bk and b : Since the rank of A is 3, Qk and bk are calculated from the first 3
rows of Q, and Q and bk are calculated from the last row of Q. Calculating Q first since it is easiest:

Q = Qe4 eT4 QT

1/2 1/2 1/2 1/2 0 1/2 1/2 1/2 1/2
1/2 1/2 1/2 1/2  1/2 1/2 1/2 1/2
0 0 0 0

Q =
1/2 1
1/2 1/2 1/2 0 1/2 1/2 1/2 1/2
1/2 1/2 1/2 1/2 1 1/2 1/2 1/2 1/2

1/4 1/4 1/4 1/4
1/4 1/4 1/4 1/4
Q =
1/4 1/4 1/4

1/4
1/4 1/4 1/4 1/4
Thus, b is given by,
1/4 1/4 1/4 1/4 1 3/2
1/4 1/4 1/4 1/4 2 3/2
b =
1/4 1/4 1/4 1/4 1 = 3/2

1/4 1/4 1/4 1/4 6 3/2

Since b = bk + b :

1 3/2 1/2
2 3/2 7/2
1 3/2 = 1/2
bk =

6 3/2 9/2
Finally, Qk is given by:

1/2 1/2 1/2
1/2 1/2 1/2
Qk =
1/2 1/2 1/2

1/2 1/2 1/2

Define R4 , and compute A from Qk and R4 :

6 1 12
R4 = 0 3 4
0 0 8

11
 1 T
A = R4 Qk

1/6 1/18 2/9 1/2 1/2 1/2 1/2
1
R4 = 0 1/3 1/6 , and QTk = 1/2 1/2 1/2 1/2
0 0 1/8 1/2 1/2 1/2 1/2
Therefore,

1/6 1/18 2/9 1/2 1/2 1/2 1/2 1/18 0 1/6 2/9
A = 0 1/3 1/6 1/2 1/2 1/2 1/2 = 1/12 1/4 1/4 1/12
0 0 1/8 1/2 1/2 1/2 1/2 1/16 1/16 1/16 1/16

12
References
[1] Justin Eilertsen. Course Notes. [Householder Reflectors (Revisited), Appendix A: Givens Rotations]. Florida
State University, 2017.

[2] Daisuke Miyazaki. Inverse Matrix of 2 2 Matrix, 3 3 Matrix, 4 4 Matrix. Hiroshima City University,
n.d. Web. 20 Mar. 2017
http://www.cg.info.hiroshima-cu.ac.jp/ miyazaki/knowledge/teche23.html

[3] Che-Rung Lee, Sheng-Kai Wu. Lecture Notes 6: Givens Rotation. Hsinchu City, China, 2011.
http://www.cs.nthu.edu.tw/ cherung/teaching/2011anm/note07.pdf

[4] Tirupathi R. Chandruptla, Eric Constans. The Use of Householder Reflection in Place of Givens Rotation in
Matrix Decomposition. Advances in Engineering Software 32.8, 2001.
https://ocw.mit.edu/courses/mathematics/18-335j-introduction-to-numerical-methods-fall-2010/lectur

[5] P. Bogacki. Transforming a matrix to reduced row echelon form. Linear Algebra Toolkit, 2000-2013.
http://www.math.odu.edu/ bogacki/cgi-bin/lat.cgi

13