【矩阵计算】QR分解-基于Householder变换

一、QR分解

QR分解是将一个矩阵分解为正交矩阵和三角矩阵的乘积。

QR分解被广泛应用于线性最小二乘问题的求解和矩阵特征值的计算。

定义2.4.1 如果实矩阵A∈R^(m×n)能化成正交矩阵Q∈R^(m×m)与上三角矩阵R∈R^(m×n)的乘积，即A=QR，则称其为A的QR分解。

二、QR分解存在性证明：基于Householder变化实现

已知，通过Householder变换，我们可以将任何一个非零向量x∈R^n转化为‖x‖2*e1，即除第一个元素外，其它元素均为零。

下面通过Householder变化来实现矩阵的QR分解。

仅考虑m=n时的情形。设矩阵A∈R^(n×n)，令H1∈R^(n×n)为一个Householder变化，满足

于是，利用分块矩阵的运算，得到对矩阵左乘一次Householder变换的结果（使该矩阵第一列变成只有第一个元素非零的向量）

其中，A2阶数为n-1。

类似地，继续对A2左乘Householder矩阵，直到最后整个矩阵变成一个n阶的下三角阵。

不断重复上述过程，这样，我们得到一系列Householder矩阵

最终得到一个下三角阵

回看QR分解的主要目的：将矩阵分解成正交阵与三角阵相乘。

我们现在已经得到三角阵，即QR分解中的R，那么Q是什么呢？

已知A=QR, Hn-1...H2H1A=R，并且Hi正交且正交，得到A=H1H2...Hn-1R于是记Q=H1H2...Hn-1，并且可知Q是正交矩阵，即A=QR

于是得到利用Householder变换进行QR分解的思路：逐一左乘上Householder矩阵，Householder矩阵的阶数从一开始与A一致的n递减，矩阵规模逐步缩小至1阶，结束循环，得到三角阵R。

那么，Q要怎么返回呢？

已知Q是n个Householder矩阵的乘积，即Q=Q1Q2Q3...Qn，其中

Qj=Im-βj*v^(j)*[v^(j)]^T

其中，v^(j)是[0,...,0,vj+1^(j),...,vm(j)]T

假设有矩阵A，在A(j+1:m,j)中存储了第j个Householder向量的基本部分v^(j)(j+1:m)

为了算Q（Householder矩阵乘积），有两种做法：向前累积、向后累积。一般会使用向后累积。

三、矩阵Q的向后累积算法

结合Qj的计算公式

Qj=Im-βj*v^(j)*[v^(j)]^T

求解QR分解中矩阵Q的算法如下

给定矩阵A∈R^(m×n)，假设m≥n

Q=Im

for j=n:-1:1 %因为Q是若干个从1阶到n阶的Householder矩阵的乘积，所以Qj的下标j也从n开始，步长-1，到1

Q=QjQ(左乘Househloder矩阵，改变行)

Q(j:m,j:m)=Q(j:m,j:m)-(βj*v(j:m))(v(j:m)^T*Q(j:m,j:m))(β*v)(v^T*β)

%从上面推导过程的图可以看到，第j个Householder矩阵，是与Q的每一行从第j列开始相乘，直到m，所以是j:m；同样，是与Q的每一列从第j行开始相乘直到m，所以是j:m。其他变量的范围也与这个对应。

四、矩阵R的算法

注意，此处的v与β都是使用之前博文中的Householder变换得到的返回值。

五、算法不足

上述算法适用情况是m≥n，即行数≥列数的情况，对于m＜n的情况，matlab会报错，需要进一步完善算法讨论。