The QR Factorization - USM

1 Jim LambersMAT 610 Summer Session 2009-10 Lecture 9 NotesThese notes correspond to Section in the FactorizationLet be an matrix with full column rank. The factorizationof is a decomposition = , where is an orthogonal matrix and is an upper triangular are three ways to compute this decomposition:1. Using Householder matrices, developed by Alston S. Householder2. Using Givens rotations, also known as Jacobi rotations, used by W. Givens and originallyinvented by Jacobi for use with in solving the symmetric eigenvalue problem in A third, less frequently used approach is theGram-Schmidt RotationsWe illustrate the process in the case where is a 2 2 matrix. In Gaussian elimination, we compute 1 = where 1is unit lower triangular and is upper triangular.

2 Specifically, 10 211 11 12 21 22 =" (2)11 (2)120 (2)22#, 21= 21 contrast, the decomposition computes = , or 11 12 21 22 = 11 120 22 ,where 2+ 2= the relationship 11+ 21= 0 we obtain 21= 11 2 221= 2 211= (1 2) 211which yields = 11p 221+ is conventional to choose the + sign. Then, we obtain 2= 1 2= 1 211 221+ 211= 221 221+ 211,or = 21p 221+ , we choose the + sign. As a result, we have 11= 11 11p 221+ 211+ 21 21p 221+ 211=q 221+ matrix = is called aGivens rotation. It is called a rotation because it is orthogonal, and therefore length-preserving, and also because there is an angle such that sin = and cos = , and its effect isto rotate a vector clockwise through the angle . In particular, = 0 where =p 2+ 2, = cos and = sin.

3 It is easy to verify that the product of tworotations is itself a rotation. Now, in the case where is an matrix, suppose that we havethe vector ..2 Then .. = .. 0 ..So, in order to transform into an upper triangular matrix , we can find a product of rotations such that = . It is easy to see that ( 2) rotations are required. Each rotation takes ( )operations, so the entire process of computing the Factorization requires ( 3) is important to note that the straightforward approach to computing the entries and ofthe Givens rotation, = p 2+ 2, = p + 2 + 2,is not always advisable, because in floating-point arithmetic, the computation ofp 2+ 2couldoverflow.

4 To get around this problem, suppose that . Then, we can instead compute = , =1 1 + 2, = ,which is guaranteed not to overflow since the only number that is squared is less than one inmagnitude. On the other hand, if , then we compute = , =1 1 + 2, = .Now, we describe the entire algorithm for computing the Factorization using Givens rota-tions. Let [ , ] =givens( , ) be aMatlab-style function that computes and such that = 0 , =p 2+ , let ( , , , ) be the Givens rotation matrix that rotates the th and th elements of avectorvclockwise by the angle such that cos = and sin = , so that if = and = ,and [ , ] =givens( , ), then in the updated vectoru= ( , , , ) v, = = 2+ 2and = 0. The Factorization of an matrix is then computed as = = for = 1 : dofor = : 1 : + 1do[ , ] =givens( 1, , ) = ( , , , ) = ( , , , )endendNote that the matrix is accumulated bycolumnrotations of the identity matrix, because thematrix by which is multiplied to reduce to upper-triangular form, a product of row rotations,is.

5 ExampleWe use Givens rotations to compute the Factorization of = .First, we compute a Givens rotation that, when applied to 41and 51, zeros 51: = .Applying this rotation to rows 4 and 5 yields 1 0 0000 1 0000 0 1000 0 0 0 0 = .Next, we compute a Givens rotation that, when applied to 31and 41, zeros 41: = .Applying this rotation to rows 3 and 4 yields 1 000 00 100 00 0 00 0 00 000 1 = .4 Next, we compute a Givens rotation that, when applied to 21and 31, zeros 31: =.

6 Applying this rotation to rows 2 and 3 yields 100 0 00 0 00 0 0000 1 0000 0 1 = .To complete the first column, we compute a Givens rotation that, when applied to 11and 21,zeros 21: = .Applying this rotation to rows 1 and 2 yields 0 0 0 0 000 1 0 000 0 1 000 0 0 1 = .Moving to the second column, we compute a Givens rotation that, when applied to 42and 52,zeros 52: = .Applying this rotation to rows 4 and 5 yields 1 0 0000 1 0000 0 1000 0 0 0 = .Next, we compute a Givens rotation that, when applied to 32and 42, zeros 42: =.

7 5 Applying this rotation to rows 3 and 4 yields 1 000 00 100 00 00 0 00 000 1 = .Next, we compute a Givens rotation that, when applied to 22and 32, zeros 32: = .Applying this rotation to rows 3 and 4 yields 100 0 00 0 00 0 0000 1 0000 0 1 = .Moving to the third column, we compute a Givens rotation that, when applied to 43and 53, zeros 53: = .Applying this rotation to rows 4 and 5 yields 1 0 0000 1 0000 0 1000 0 0 0 0 = .Finally, we compute a Givens rotation that, when applied to 33and 43, zeros 43: =.

8 Applying this rotation to rows 3 and 4 yields 1 000 00 100 00 0 00 0 00 000 1 = = .6 Applying the transpose of each Givens rotation, in order, to thecolumnsof the identity matrixyields the matrix = such that = is upper triangular. Householder ReflectionsIt is natural to ask whether we can introduce more zeros with each orthogonal rotation. To thatend, we examineHouseholder reflections. Consider a matrix of the form = uu , whereu =0and is a nonzero constant. It is clear that is a symmetric rank-1 change of . Can wechoose so that is also orthogonal? From the desired relation = we obtain = ( uu ) ( uu )= 2 uu + 2uu uu = 2 uu + 2(u u)uu = + ( 2u u 2 )uu = + ( u u 2)uu.

9 It follows that if = 2/u u, then = for any nonzerou. Without loss of generality, we canstipulate thatu u= 1, and therefore takes the form = 2vv , wherev v= is the matrix called a reflection? This is because for any nonzero vectorx, xis thereflection ofxacross the hyperplane that is normal tov. To see this, we consider the 2 2 caseand setv= 1 0 andx= 1 2 . Then = 2vv = 2 10 1 0 = 1 00 1 2 1 00 0 = 1 001 Therefore x= 1 001 12 = 12 .7 Now, letxbe any vector. We wish to construct so that x= e1for some , wheree1= 1 0 0 . From the relations x 2= x 2, e1 2= e1 2= ,we obtain = x 2. To determine , we begin with the equation x= ( 2vv )x=x 2vv x= , we obtain12(x e1) = (v x) follows that the vectorv, which is a unit vector, must be a scalar multiple ofx e1.

10 Therefore,vis defined by the equations 1= 1 x e1 2= 1 p x 22 2 1+ 2= 1 2 2 2 1= 1p2 ( 1)= sgn( )r 12 , 2= 2p2 ( 1)= 22 1,.. = 2 avoid catastrophic cancellation, it is best to choose the sign of so that it has the opposite signof 1. It can be seen that the computation ofvrequires about 3 that the matrix is never formed explicitly. For any vectorb, the product bcan becomputed as follows: b= ( 2vv )b=b 2(v b) process requires only 4 operations. It is easy to see that we can represent simply by , suppose that thatx=a1is the first column of a matrix . Then we construct aHouseholder reflection 1= 2v1v 1such that x= e1, and we have (2)= 1 = 11 12 1 (2)2: ,2 a(2)2: , 0 .where we denote the constant by 11, as it is the (1,1) element of the updated matrix (2).