Gram-Schmidt Process - Matrix Computation - Lecture Slides, Slides of Advanced Computer Architecture

Download Gram-Schmidt Process - Matrix Computation - Lecture Slides and more Slides Advanced Computer Architecture in PDF only on Docsity!

• Lecture

Overview

Gram-Schmidt process

QR decomposition

Gram-Schmidt triangular orthogonalization

Examples

Rotation by θ in IR

2 is given by

y = Uθx, Uθ =

[

cos θ − sin θ

sin θ cos θ

]

since e 1 = [cos θ, sin θ]

> , e 2 = [− sin θ, cos θ]

>

Reflection across line x 2 = x 1 tan(θ/2) is given by

y = Uθx, Uθ =

[

cos θ sin θ

sin θ − cos θ

]

since e 1 = [cos θ, sin θ]

> , e 2 = [sin θ, − cos θ]

>

θ θ x 1 x 1

x 2 x 2

e 1 e 1

e 2 e 2

rotation reflection

can check that Uθ and Rθ are orthogonal

Gram-Schmidt process

Given independent vectors x 1 ,... , xn ∈ IR

m , Gram-Schmidt process

finds orthonormal vectors, q 1 ,... , qn such that

span(x 1 ,... , xr ) = span(q 1 ,... , qr )

Thus, q 1 ,... , qr are orthonormal basis for span(x 1 ,... , xr )

Idea: first orthogonalized each vector w.r.t. previous ones and then

normalize result to have unit norm

I (^) v 1 = x 1

I (^) q 1 = v^1 ‖v 1 ‖ (normalize) I (^) v 2 = x 2 − (q> 1 x^2 )q^1 (remove^ q^1 component from^ x^2 ) I (^) q 2 = v^2 ‖v 2 ‖

(normalize)

I (^) v 3 = x 3 − (q> 1 x^3 )q^1 −^ (q

> 2 x^3 )q^2 (remove,^ q^1 ,^ q^2 components) I (^) q 3 = v^3 ‖v 3 ‖ I (^) etc.

Find orthonormal basis incrementally

QR decomposition

In matrix form, A = QR where A ∈ IR

m×n , Q ∈ IR

m×n , R ∈ IR

n×n :

[

x 1 x 2 · · · xn

]

[

q 1 q 2 · · · qn

]

r 11 r 12 · · · r 1 n

0 r 22 · · · r 2 n

. . .

0 0 · · · rnn

A Q R

Q

> Q = I , and R is upper triangular and invertible

Usually computed using a variation of Gram-Schmidt process which is

less sensitive to numerical errors

Can also be computed by Householder transformation or Givens

rotations

Columns of Q are orthonormal basis in ran(A)

Q

> A = R

General Gram-Schmidt process

If x 1 ,... , xn are dependent, we find vj = 0 for some j, which means

xj is linearly dependent on x 1 ,... , xj− 1

Modified algorithm: when we have vj = 0, skip to the next vector

xj+1 and continue

k = 0

for i = 1,... , n

{

v = xi −

∑k

j=

qj q> j

xi ;

if v 6 = 0 {k = k + 1; qk =

v ‖v‖

Example (cont’d)

Then we have

R = Q

> A =

A =

√^1 2

√^1 3 0 1 0

0 0

√^1 2

√^1 3

Properties of QR decomposition

Find orthonormal basis for ran(A) directly

Let A = BC with B ∈ IR

m×p , C ∈ IR

p×n , p = rank(A)

To check whether y ∈ span(x 1 ,... , xn): apply Gram-Schmidt

procedure to [x 1 ,... , xn, y]

If A = Q 1 R 1 with rank(A) = p, the full QR factorization is

A =

[

Q 1 Q 2

]

[

R 1

]

where

[

Q 1 Q 2

]

is orthogonal, i.e., columns of Q 2 ∈ IR

m×(m−p) are

orthonormal, orthogonal to Q 1

To find Q 2 , one can use any matrix A˜ s.t. [A ˜A] is full rank (e.g.,

˜A = I ), and then apply general Gram-Schmidt process

Q 1 are orthonormal vectors obtained form columns of A

Q 2 are orthonormal vectors obtained from extra columns ˜A

Least squares via QR decomposition

A ∈ IR

m×n , A = QR with Q

> Q = I , R ∈ IR

n×n

The pseudo inverse is

(A

> A)

− 1 A

> = (R

> Q

> QR)

− 1 R

> Q

> = R

− 1 Q

>

The projected point, xls = R

− 1 Q

> y

The projection matrix

A(A

> A)

− 1 A

> = AR

− 1 Q

> = QQ

>

Least squares via full QR factorization

Full QR factorization

A =

[

Q 1 Q 2

]

[

R 1

]

where

[

Q 1 Q 2

]

∈ IR

m×m is orthogonal, R 1 ∈ IR

n×n is upper

triangular and invertible

Multiplication by orthogonal matrix

‖Ax − y‖

2

[

Q 1 Q 2

]

[

R 1

]

x − y

2

[

Q 1 Q 2

]> [

Q 1 Q 2

]

[

R 1

]

x −

[

Q 1 Q 2

]>

y

2

[

R 1 x − Q

> 1 y

−Q

> 2 y

] ∥

2

∥ R

1 x^ −^ Q

> 1

y

2

∥ Q>

2

y

2

Can be minimized by choosing xls = R

− 1 1

Q

> 1 y

Gram-Schmidt as triangular orthogonalization

Steps:

I (^) v 1 = x 1

I (^) q 1 = v^1 ‖v 1 ‖

(normalize)

I (^) v 2 = x 2 − (q> 1 x^2 )q^1 = (I^ −^ q^1 q

> 1 )x^2 (remove^ q^1 component from^ x^2 ) I (^) q 2 = v^2 ‖v 2 ‖ (normalize) I (^) v 3 = x 3 − (q> 1 x^3 )q^1 −^ (q

> 2 x^3 )q^2 = (I^ −^ q^1 q

> 1 −^ q^2 q

> 2 )x^3 = (I − q 1 q> 1 )(I − q 2 q> 2 )x 3 (remove, q 1 , q 2 components)

I (^) q 3 = v^3 ‖v 3 ‖ I (^) etc.

[

x 1 x 2 · · · xn

]

[

q 1 q 2 · · · qn

]

r 11 r 12 · · · r 1 n

r 22 · · ·

.. .

rnn

Gram-Schmidt as triangular orthogonalization (cont’d)

[

x 1 x 2 · · · xn

]

[

q 1 q 2 · · · qn

]

r 11 r 12 · · · r 1 n

r 22 · · ·

.. .

rnn

x 1 = r 11 q 1

x 2 = r 12 q 1 + r 22 q 2 =

r 12 r 11

x 1 + r 22 q 2 = (q

> 1 x^2 )q^1 +^ r^22 q^2

x 3 = r 13 q 1 + r 23 q 2 + r 33 q 3 = (q

> 1 x^3 )q^1 + (q

> 2 x^3 )q^2 +^ r^33 q^3 . . .

xn = r 1 nq 1 + r 2 nq 2 + · · · + rnnqn

xi = (q

> 1 xi^ )q^1 + (q

> 2 xi^ )q^2 +^ · · ·^ + (q

> i− 1

xi )qi− 1 + ‖ q˜i ‖qi

= r 1 i q 1 + r 2 i q 2 + · · · + rii qi

rij = q

> i xj^ (i^6 =^ j)

Gram-Schmidt as triangular orthogonalization (cont’d)

Let A ∈ IR

m×n , m ≥ n be a matrix of full rank with columns xi ,

consider the sequence of formulas

q 1 =

P 1 x 1

‖P 1 x 1 ‖

, q 2 =

P 2 x 2

‖P 2 x 2 ‖

,... , qn =

Pnxn

‖Pnxn‖

where Pj ∈ IR

m×m of rank m − (j − 1) projects x onto the space

orthogonal to q 1 ,... , qj− 1 (P 1 = I when j = 1)

Projection Pj can be represented explicitly

Let Qj− 1 denote the m × (j − 1) matrix containing the first j − 1

columns of Q

Qj− 1 =

[

q 1 q 2 · · · qj− 1

]

, then

Pj = I − Qj− 1 Q

> j− 1

Modified Gram-Schmidt process

Recall for rank-one orthogonal projection with q ∈ IR

m

Pq =

qq

>

q>q

The complements are rank m − 1 orthogonal projections

P⊥q = I −

qq

>

q>q

By definition of Pj

Pj = P⊥q j− 1

· · · P⊥q 2 P⊥q 1

where P 1 = I , and thus

vj = P⊥q j− 1

· · · P⊥q 2 P⊥q 1 xj