-92第三讲线性最小二乘问题3.3QR分解QR分解是将一个矩阵分解一个正交矩阵(酉矩阵)和一个三角矩阵的乘积.QR分解被广泛应用于线性最小二乘问题的求解和矩阵特征值的计算3.3.1QR分解的存在性与唯一性定理3.9(QR分解)[73] 设 A e Cmxn (m ≥ n),则存在一个单位列正交矩阵 Q e Cmxn (即Q*Q=Inxn)和一个上三角矩阵RECnxn,使得(3.8)A=QR.若A列满秩,则存在一个具有正对角线元素的上三角矩阵R使得(3.8)成立,且此时QR分解唯,即Q和R都唯一.证明.设A=[a1,a2,...,an)ECmxn.若A列满秩,即rank(A)=n.则QR分解(3.8)就是对A的列向量组进行Gram-Schmidt正交化过程的矩阵描述(见算法3.3)算法3.3.Gram-Schmidt正交化过程1: r11 = a1ll22:q1=α1/r113: for j = 2 to n do4:qj = aj5.fori= l to j-l do6:%计算内积,用于正交化rij=(aj,qi)1qj=qj-rqiend for8:9:r=lqil210:qi=qi/rjjI1: end for由算法3.3可知r1jr2jaj = r1jq1 +r2jq2 +.*++ rjqi = [q1,q2, .., qi]=2,3,...,n.al=r11q1,.rii记Q=[91,q2.., 9n], R=[ri]nxn,其中fori≤jqtaj,(3.9)rij0,fori>j
· 92 · 第三讲 线性最小二乘问题 3.3 QR 分解 QR 分解是将一个矩阵分解一个正交矩阵 (酉矩阵) 和一个三角矩阵的乘积. QR 分解被广泛 应用于线性最小二乘问题的求解和矩阵特征值的计算. 3.3.1 QR 分解的存在性与唯一性 定理 3.9 (QR 分解) [73] 设 A ∈ C m×n (m ≥ n), 则存在一个单位列正交矩阵 Q ∈ C m×n (即 Q∗Q = In×n) 和一个上三角矩阵 R ∈ C n×n , 使得 A = QR. (3.8) 若 A 列满秩, 则存在一个具有正对角线元素的上三角矩阵 R 使得 (3.8) 成立, 且此时 QR 分解唯 一, 即 Q 和 R 都唯一. 证明. 设 A = [a1, a2, . . . , an] ∈ C m×n . 若 A 列满秩, 即 rank(A) = n. 则 QR 分解 (3.8) 就是对 A 的列向量组进行 Gram-Schmidt 正交化过程的矩阵描述 (见算法 3.3). 算法 3.3. Gram-Schmidt 正交化过程 1: r11 = ∥a1∥2 2: q1 = a1/r11 3: for j = 2 to n do 4: qj = aj 5: for i = 1 to j − 1 do 6: rij = (aj , qi) % 计算内积, 用于正交化 7: qj = qj − rijqi 8: end for 9: rjj = ∥qj∥2 10: qj = qj/rjj 11: end for 由算法 3.3 可知 a1 = r11q1, aj = r1jq1 + r2jq2 + · · · + rjjqj = [q1, q2, . . . , qj ] r1j r2j . . . rjj , j = 2, 3, . . . , n. 记 Q = [q1, q2, . . . , qn], R = [rij ]n×n, 其中 rij = q ∗ i aj , for i ≤ j 0, for i > j (3.9)
3.3QR分解.93.于是Gram-Schmidt正交化过程可表示为r11r12.Tinr22r2n即A=QR[a1, a2,..,an] = [q1, q2,..,qn]rn-Inn如果A不是列满秩,我们可以通过下面的方式做类似的正交化过程:·如果a1=0,则令q1=0;否则令q1=a1/la1l2;·对于=2,3.,n,计算=aj--(gia,)i如果j=0,则表明aj可以由a1,a2....,aj-1线性表出,令qj=0.否则令qj=/g;ll2于是我们有A=QR,其中Q=[91,92,.….,9n]列正交(但不是单位列正交),其列向量要么是单位向量,要么就是零向量.R=[ri]nxn的定义同(3.9).需要注意的是,如果Q的某一列qk=0,那么R中对应的第k行就全部为0.设rank(A)=1<n,则Q有1个非零列,设为qi,qiz.,ir它们形成Cm中的一个单位正交向量组,所以我们可以将其扩展成Cm中的一组标准正交基,即lii,liz....,qinq1.,qm-然后我们用91替换Q中的第一个零列,用2替换Q中的第二个零列,依此类推,将Q中的所有零列都替换掉.将最后得到的矩阵记为Q,则QeCmxn单位列正交,且QR=QR.这是由于Q中的新添加的列向量正好与R中的零行相对应.所以我们有QR分解A=QR.下面证明满秩矩阵QR分解的存在唯一性存在性:由于A列满秩,由Gram-Schmidt正交化过程(算法3.3)可知,存在上三角矩阵R=[ri]nxn满足rji>0,使得A=QR,其中Q单位列正交唯一性:假设A存在QR分解A= QiR1 = Q2R2其中Q1,Q2ECmxn单位列正交,R1,R2ECnxn为具有正对角元素的上三角矩阵.则有Q1= Q2R2R-1(3.10)于是1=Q1ll2=IlQ2R2R-12=R2R-12又R1,R2均为上三角矩阵,所以R2R-1也是上三角矩阵,且其对角线元素为R2(i,i)/Ri(i,i)
3.3 QR 分解 · 93 · 于是 Gram-Schmidt 正交化过程可表示为 [a1, a2, . . . , an] = [q1, q2, . . . , qn] r11 r12 · · · r1n r22 · · · r2n . . . rn−1,n rnn , 即 A = QR. 如果 A 不是列满秩, 我们可以通过下面的方式做类似的正交化过程: • 如果 a1 = 0, 则令 q1 = 0; 否则令 q1 = a1/∥a1∥2 ; • 对于 j = 2, 3, . . . , n, 计算 q˜j = aj − Pj−1 i=1 (q ∗ i aj )qi . 如果 q˜j = 0, 则表明 aj 可以由 a1, a2, . . . , aj−1 线性表出, 令 qj = 0. 否则令 qj = ˜qj/∥q˜j∥2 . 于是我们有 A = QR, 其中 Q = [q1, q2, . . . , qn] 列正交 (但不是单位列正交), 其列向量要么是单位向量, 要么就是零向量. R = [rij ]n×n 的定义同 (3.9). 需要注意的是, 如果 Q 的某一列 qk = 0, 那么 R 中对应的第 k 行就 全部为 0. 设 rank(A) = l < n, 则 Q 有 l 个非零列, 设为 qi1 , qi2 , . . . , qil . 它们形成 C m 中的一个单位正 交向量组, 所以我们可以将其扩展成 C m 中的一组标准正交基, 即 qi1 , qi2 , . . . , qil , q˜1, . . . , q˜m−l . 然后我们用 q˜1 替换 Q 中的第一个零列, 用 q˜2 替换 Q 中的第二个零列, 依此类推, 将 Q 中的所有 零列都替换掉. 将最后得到的矩阵记为 Q˜, 则 Q˜ ∈ C m×n 单位列正交, 且 QR˜ = QR. 这是由于 Q˜ 中的新添加的列向量正好与 R 中的零行相对应. 所以我们有 QR 分解 A = QR. ˜ 下面证明满秩矩阵 QR 分解的存在唯一性. 存在性: 由于 A 列满秩, 由 Gram-Schmidt 正交化过程 (算法 3.3) 可知, 存在上三角矩阵 R = [rij ]n×n 满足 rjj > 0, 使得 A = QR, 其中 Q 单位列正交. 唯一性: 假设 A 存在 QR 分解 A = Q1R1 = Q2R2, 其中 Q1, Q2 ∈ C m×n 单位列正交, R1, R2 ∈ C n×n 为具有正对角元素的上三角矩阵. 则有 Q1 = Q2R2R −1 1 . (3.10) 于是 1 = ∥Q1∥2 = ∥Q2R2R −1 1 ∥2 = ∥R2R −1 1 ∥2. 又 R1, R2 均为上三角矩阵, 所以 R2R −1 1 也是上三角矩阵, 且其对角线元素为 R2(i, i)/R1(i, i)
-94.第三讲线性最小二乘问题i=1,2...,n.因此R2(i,i)≤p(R2R-)≤IR2R-ll2≤1,i=1,2,...,nRi(i,i)同理可证Ri(i,i)/R2(i,i)≤1.所以Ri(i,i)= R2(i,i), i= 1,2,...,n又Q1呢=tr(QQ1)=n,所以由(3.10)可知R2R-=Q2R2R-呢=Q1=n由于R2R-1的对角线元素都是1,所以R2R-1只能是单位矩阵,即R2=R1.因此Q2=AR-1=口AR-1=Q1,即A的QR分解是唯一的.有时也将QR分解定义为:存在西矩阵QeCmxm使得A= QR,R1其中RECmxn是上三角矩阵如果A是实矩阵,则上面证明中的运算都可以在实数下进行,因此Q和R都可以是实矩阵否如果A是非奇异的方阵,则QR分解也可以用来求解线性方程组Ar=bGram-Schmidt正交化与正交投影Gram-Schmidt正交化过程的第i步:qi=air1i9i-r2i92-..:-rj-1,jqj-1可以看作是a,减去其在span[q1,92,….,qj-1)内的正交投影.事实上,rij9i就是a,在qi方向的正交投影,即Pspanfqjaj = (qit)aj =(qTaj)qi =rijqi下面给出QR分解的具体实现方法,分别基于MGS正交化过程,Householder变换和Givens变换3.3.2基于MGS的QR分解在证明QR分解的存在性时,我们利用了Gram-Schmidt正交化过程.但由于数值稳定性方面的原因,在实际计算中,我们一般不采用Gram-Schmidt正交化过程,取而代之的是修正的Gram-Schmidt正交化过程(modifiedGram-Schmidtprocess,MGS),即对正交化过程做如下修改:。Gram-Schmidt正交化过程的第j步:(1)计算rj=(aj,q),i=1,2,...,j-1;(2)计算qj=aj—r1jq1-r2jq2—...-rj-1,sqj-1;(3)计算=gll,q=i/s
· 94 · 第三讲 线性最小二乘问题 i = 1, 2, . . . , n. 因此 R2(i, i) R1(i, i) ≤ ρ(R2R −1 1 ) ≤ ∥R2R −1 1 ∥2 ≤ 1, i = 1, 2, . . . , n. 同理可证 R1(i, i)/R2(i, i) ≤ 1. 所以 R1(i, i) = R2(i, i), i = 1, 2, . . . , n. 又 ∥Q1∥ 2 F = tr(Q∗ 1Q1) = n, 所以由 (3.10) 可知 ∥R2R −1 1 ∥ 2 F = ∥Q2R2R −1 1 ∥ 2 F = ∥Q1∥ 2 F = n. 由于 R2R −1 1 的对角线元素都是 1, 所以 R2R −1 1 只能是单位矩阵, 即 R2 = R1. 因此 Q2 = AR−1 2 = AR−1 1 = Q1, 即 A 的 QR 分解是唯一的. □ b 有时也将 QR 分解定义为: 存在酉矩阵 Q ∈ C m×m 使得 A = QR, 其中 R = " R11 0 # ∈ C m×n 是上三角矩阵. b 如果 A 是实矩阵, 则上面证明中的运算都可以在实数下进行, 因此 Q 和 R 都可以是实矩阵. b 如果 A 是非奇异的方阵, 则 QR 分解也可以用来求解线性方程组 Ax = b. Gram-Schmidt 正交化与正交投影 Gram-Schmidt 正交化过程的第 j 步: q˜j = aj − r1jq1 − r2jq2 − · · · − rj−1,jqj−1 可以看作是 aj 减去其在 span{q1, q2, . . . , qj−1} 内的正交投影. 事实上, rijqi 就是 aj 在 qi 方向 的正交投影, 即 Pspan{qi}aj = (qiq T i )aj = (q T i aj )qi = rijqi . 下面给出 QR 分解的具体实现方法, 分别基于 MGS 正交化过程, Householder 变换和 Givens 变换. 3.3.2 基于 MGS 的 QR 分解 在证明 QR 分解的存在性时, 我们利用了 Gram-Schmidt 正交化过程. 但由于数值稳定性方面 的原因, 在实际计算中, 我们一般不采用 Gram-Schmidt 正交化过程, 取而代之的是 修正的 GramSchmidt 正交化过程 (modified Gram-Schmidt process, MGS ), 即对正交化过程做如下修改: • Gram-Schmidt 正交化过程的第 j 步: (1) 计算 rij = (aj , qi), i = 1, 2, . . . , j − 1; (2) 计算 q˜j = aj − r1jq1 − r2jq2 − · · · − rj−1,jqj−1; (3) 计算 rjj = ∥q˜j∥, qj = ˜qj/rjj ;
3.3QR分解.95.·MGS正交化过程的第i步:(1)令qj=aj;(2)计算rij=(qj,qi),qj=qj-rjqi,=1,2,...,j-1;(3)计算r=llill,gj=g/rjs可以证明,数学上这两个算法完全等价,即ri和qi都一样.但在数值上,Gram-Schmidt正交化过程是不稳定的,而MGS是向后稳定的[98].算法3.4.基于MGS的QR分解% Given A e Cmxn, compute Q = [q1, .-*, 9n] e Rmxn and R e Rnxn such that A = QR1: Set R = [ri] = Onxn (the n x n zero matrix)2: if ai=0 then3:q1 = 04:else5:r11=aill2,q1=a1/llaill26: end if7: for j = 2 to n do8:qi=aj9:fori=ltoj-ldog%MGS,注意与GS的区别10:r=(q,q),q=q-i911:end for12:if qi0then13:rj=llgill2,qi=qi/r314:end if15: end for白本算法的运算量大约为2mn2凸算法中R是一列一列计算的,因此也称为按列MGS(column-orientedMGS).按列MGS不适合列主元,因此,如果需要计算列主元QR时需采用按行MGS(row-orientedMGS),可参见相关资料 由 MGS得到的 QR分解中,QRmxn,RERnxn3.3.3基于Householder变换的QR分解由定理3.5可知,通过Householder变换,我们可以将任何一个非零变量rERn转化成Ilrl2e1,即除第一个元素外,其它都为零.下面我们就考虑通过Householder变换来实现矩阵的QR分解
3.3 QR 分解 · 95 · • MGS 正交化过程的第 j 步: (1) 令 q˜j = aj ; (2) 计算 rij = (˜qj , qi), q˜j = ˜qj − rijqi , i = 1, 2, . . . , j − 1; (3) 计算 rjj = ∥q˜j∥, qj = ˜qj/rjj ; 可以证明, 数学上这两个算法完全等价, 即 rij 和 qj 都一样. 但在数值上, Gram-Schmidt 正交化过 程是不稳定的, 而 MGS 是向后稳定的 [98]. 算法 3.4. 基于 MGS 的 QR 分解 % Given A ∈ C m×n , compute Q = [q1, . . . , qn] ∈ R m×n and R ∈ R n×n such that A = QR 1: Set R = [rij ] = 0n×n (the n × n zero matrix) 2: if a1 = 0 then 3: q1 = 0 4: else 5: r11 = ∥a1∥2, q1 = a1/∥a1∥2 6: end if 7: for j = 2 to n do 8: qj = aj 9: for i = 1 to j − 1 do % MGS, 注意与 GS 的区别 10: rij = (qj , qi), qj = qj − rijqi 11: end for 12: if qj ̸= 0 then 13: rjj = ∥qj∥2, qj = qj/rjj 14: end if 15: end for b 本算法的运算量大约为 2mn2 . b 算法中 R 是一列一列计算的, 因此也称为 按列 MGS (column-oriented MGS). 按列 MGS 不 适合列主元, 因此, 如果需要计算列主元 QR 时需采用 按行 MGS (row-oriented MGS), 可参 见相关资料. b 由 MGS 得到的 QR 分解中, Q ∈ R m×n , R ∈ R n×n . 3.3.3 基于 Householder 变换的 QR 分解 由定理 3.5 可知, 通过 Householder 变换, 我们可以将任何一个非零变量 x ∈ R n 转化成 ∥x∥2e1, 即除第一个元素外, 其它都为零. 下面我们就考虑通过 Householder 变换来实现矩阵的 QR 分解
-96.第三讲线性最小二乘问题设AERmxn(m≥n),构造Householder变换HiERmxm,使得oa11a21Hi:.0aml于是a1a12r10HiA=A20其中A2E R(m-1)x(n-1).同样地,我们可以构造一个Householder变换H2E R(m-1)x(m-1),将 A2的第一列中除第一个元素外的所有元素都化为0,即a23a2nr20H,A2=A3..0,则H2EIRmxm也是Householder变换(见习题3.8),且令H2=0H2a12ria13ain0a23a2nr20H2HiA=0A3::00不断重复上述过程,我们就可以得到一系列的Householder变换[Ik-1 0]Hk E IR(n+1-k)x(n+1-k),Hk=k =1,2,...,n,Hk0使得[ra12ain...0r2.a2n.....:R.H.H,HiA-00In000.....+:<o00由于Householder变换都是正交矩阵,令Q≤(Hn..-H2Hi)-1=H-1H,1...H-1=HiH2...Hn
· 96 · 第三讲 线性最小二乘问题 设 A ∈ R m×n (m ≥ n), 构造 Householder 变换 H1 ∈ R m×m, 使得 H1 a11 a21 . . . am1 = r1 0 . . . 0 . 于是 H1A = r1 a˜12 · · · a˜1n 0 . . . A˜ 2 0 , 其中 A˜ 2 ∈ R (m−1)×(n−1) . 同样地, 我们可以构造一个 Householder 变换 H˜ 2 ∈ R (m−1)×(m−1) , 将 A˜ 2 的第一列中除第一个元素外的所有元素都化为 0, 即 H˜ 2A˜ 2 = r2 a˜23 · · · a˜2n 0 . . . A˜ 3 0 . 令 H2 = " 1 0 0 H˜ 2 # , 则 H2 ∈ R m×m 也是 Householder 变换 (见习题 3.8), 且 H2H1A = r1 a˜12 a˜13 · · · a˜1n 0 r2 a˜23 · · · a˜2n 0 0 . . . . . . A˜ 3 0 0 . 不断重复上述过程, 我们就可以得到一系列的 Householder 变换 Hk = " Ik−1 0 0 H˜ k # , H˜ k ∈ R (n+1−k)×(n+1−k) , k = 1, 2, . . . , n, 使得 Hn · · · H2H1A = r1 a˜12 · · · a˜1n 0 r2 · · · a˜2n . . . . . . . . . 0 0 · · · rn 0 0 · · · 0 . . . . . . . . . 0 0 · · · 0 ≜ R. 由于 Householder 变换都是正交矩阵, 令 Q ≜ (Hn · · · H2H1) −1 = H −1 1 H −1 2 · · · H−1 n = H1H2 · · · Hn