前言

笔者认为格拉姆-施密特正交化与QR分解本质上没有区别，它们只是表示形式的不同，一个通常意义上说是代数学范畴，而另一个是矩阵分解范畴，故在本篇文章笔者不会去显式地去区分格拉姆-施密特正交化和QR分解

笔者发现格拉姆-施密特正交化只是QR分解的一种方法（简化QR分解），此外还可通过HouseHoldolder变换（最稳定，文章主要讲解这种方法）和Givens旋转进行（完整QR分解）

请注意本篇文章的大部分都专注于格拉姆—施密特正交化和QR分解背后的思维逻辑，笔者从不认为它们的结论是重要的，而思考为什么要这样做这个问题，对于我们很重要，所以在阅读此文章时笔者希望读者先暂时忘记结论，而尝试思考问题本身

为什么？

笔者认为，在研究一个问题之前了解我们为什么要这样做，这样做有什么好处和用途对于我们理解这个问题是至关重要的，事实上，笔者一直在思考的一个问题就是：
::: align-center
为什么我们一定要正交化？直接使用标准基\{e_1,e_2,\cdots,e_n\}不好吗？

:::
然而，事实上，标准基虽然方便，但在某些情况下可能不适用，原因就是它们太“标准”了，以至于我们在处理一些问题需要灵活的选取有利于问题的基底时再使用它们显得太过教条，后果可能就会带来表示的不直观和计算的困难。在文章后面我们将提供一些施密特正交化后得到基的一些很好的性质，相信读者在看到后会明白它在某些方面存在的必要性。

需要解决的问题

在前一节的基础上，我们需要暂时抛开有标准基这个工具的前提，想象
::: align-center
给定一个已知为线性无关的向量组\{a_1,a_2,\cdots,a_n\}，我们如何从它们之中得出一组相互正交的基？

:::
读者请自行思考为什么我们需要基是相互正交的，是否是因为在基相互正交时我们也有一些比较好的性质，这是否与我们前一节提到的问题有异曲同工之妙？这将在接下来的问题中至关重要

我们的第一个基

万事开头难，但对于这个问题笔者认为它开头难，但又不难。
在找到一整组基之前，我们首先需要考虑的一个问题就是
::: align-center
是不是要先把第一个基给找出来？

:::
怎么找呢？实际上谜底已经在谜面上了，笔者由衷想要告诉读者一定不要思维定式地由标准基的概念认为我们的第一个基也必须是某个元素为1，其他为0的样子，事实上若读者这样想则不妨回头看第一节的内容，千万不要忘记我们为何而出发。

事实上，我们的第一个基完全可以选择向量组中的任何一个向量不是吗？没有人强迫我们不能这样做，那我们的第一个基就完全可以选择此向量组的任何一个向量，而我们可以选择第一个基这件事就意味着这个正交化方法注定比用标准基更加灵活，更能适应我们的需求

那我们的第一个基\widetilde{q_1}=a_1

其他基？

在拥有了第一个基以后，读者自然而然地就会想到如何获得其他的基，事实上由于向量组的其他向量彼此之间的关系（这里主要指正交性）都是未知的，想从这些“杂乱无章”的向量中获得一个组相互正交的基看似并不容易，但笔者在这里不希望读者去过多思考
::: align-center
我该怎么拥有这组基？

:::
而读者应该在发现思考此问题较为困难后率先思考
::: align-center
如果我有了这组基我可以怎么样？

:::
诚然，在现实生活中不考虑过程直接去考虑结果是非常危险的，但在数学的世界里没有人限制你不能这么做

在中学的知识下我们知道，有了一组基我们可以表示这组基张成空间内的任意向量：
::: align-center
那如果向量组中的其他向量就是这组基张成的呢？

:::
接着想
::: align-center
如果我现在已经知道了结果（这些向量），我能知道它们从何而来吗（基）？

:::
事实上：
::: align-center
答案是可以的

:::

第二个基

我们不妨从最为简单的问题开始吧，如果现在我们要开始找第二个基，按照假设它已经已知了，而向量组中的某个向量（我们不妨令它就是第二个）就是它和第一个基张成的，我们可以得到第一个基吗？
::: align-center
当然可以！

:::
这里我们需要重新回归前面提到的正交性，这里我们需要提前知道在正交系下投影的计算非常容易，事实上在学习完内积的概念后读者应知道在正交系下a_2在\widetilde{q_1}下的投影向量为
::: align-center
\widetilde{q_1}'=\frac{\widetilde{q_1}^Ta_2}{\widetilde{q_1}^Tq_1}\widetilde{q_1}

:::
而又由向量的加减法我们知道，当\widetilde{q_1}'已知时我们可以由平行四边形法则方便的算出a_2在\widetilde{q_2}上的分量\widetilde{q_2}'

也就是
::: align-center
a_2-\widetilde{q_1}'=a_2-\frac{\widetilde{q_1}^Ta_2}{\widetilde{q_1}^Tq_1}\widetilde{q_1}

:::
这时我们就得到了第二个基“\widetilde{q_2}”，因为我们并不关心它是\widetilde{q_2}还是\widetilde{q_2}'，它们都是共线的，如果需要标准基我们只需要将它们单位化即可

更多的基

类似的，我们自然而然的可以将这个过程推广到寻找第三个,那么我们就可以假设向量组中的第三个向量a_3是被前两个基以及我们的第三个基表示的，如此我们就可以求得第三个基。类似地，我们可以推广至第四个...第n个基上，那么格拉姆，施密特正交化的基本方法也就被我们掌握了

格拉姆—施密特正交化

设a_1,a_2,\cdots,a_n是\mathbb{R}^n的一组基，格拉姆-施密特正交化的计算过程分为两步

正交化

\begin{aligned} &\widetilde{q_1}=a_1\\ &\widetilde{q_2}=a_2-\frac{\widetilde{q_1}^Ta_2}{\widetilde{q_1}^T\widetilde{q_1}}\widetilde{q_1}\\ &\widetilde{q_3}=a_3-\frac{\widetilde{q_1}^Ta_3}{\widetilde{q_1}^T\widetilde{q_1}}\widetilde{q_1}-\frac{\widetilde{q_2}^Ta_3}{\widetilde{q_2}^Tq_2}\widetilde{q_2}\\ &\vdots\\ &\widetilde{q_n}=a_n-\frac{\widetilde{q_1}^Ta_n}{\widetilde{q_1}^T\widetilde{q_1}}\widetilde{q_1}-\frac{\widetilde{q_2}^Ta_n}{\widetilde{q_2}^Tq_2}\widetilde{q_2}\cdots-\frac{\widetilde{q_{n-1}}^Ta_n}{\widetilde{q_{n-1}}\widetilde{q_{n-1}}}\widetilde{q_{n-1}}\\ \end{aligned}

单位化

第二步再单位化每个向量，得到标准正交基：
::: align-center
q_i=\frac{\widetilde{q_i}}{||\widetilde{q_i}||},i=1,2,\cdots,n

:::

QR分解

同LU分解一样，QR分解实际上是格拉姆-施密特正交化在矩阵运算上的体现，这一部分笔者将在近期以矩阵代数的形式给出

QR分解问题

::: align-center
对于列向量线性无关的A_{m\times n}，试求A=QR，其中Q_{m\times n}的列向量组成标准正交基，R_{m\times n}为上三角矩阵

:::

为什么？

这里读者可能又会询问
::: align-center
为什么我们要将A这样分解为QR？Q是写为这样是比较容易理解的，但为什么R要写为这样？

:::
事实上，由读者不妨思考如果我们需要在矩阵中执行格拉姆-施密特正交化时会发生什么，再次回顾格拉姆-施密特正交化的过程
第一步：
\begin{aligned} &\widetilde{q_1}=a_1\\ &\widetilde{q_2}=a_2-\frac{\widetilde{q_1}^Ta_2}{\widetilde{q_1}^T\widetilde{q_1}}\widetilde{q_1}\\ &\widetilde{q_3}=a_3-\frac{\widetilde{q_1}^Ta_3}{\widetilde{q_1}^T\widetilde{q_1}}\widetilde{q_1}-\frac{\widetilde{q_2}^Ta_3}{\widetilde{q_2}^Tq_2}\widetilde{q_2}\\ &\vdots\\ &\widetilde{q_n}=a_n-\frac{\widetilde{q_1}^Ta_n}{\widetilde{q_1}^T\widetilde{q_1}}\widetilde{q_1}-\frac{\widetilde{q_2}^Ta_n}{\widetilde{q_2}^Tq_2}\widetilde{q_2}\cdots-\frac{\widetilde{q_{n-1}}^Ta_n}{\widetilde{q_{n-1}}\widetilde{q_{n-1}}}\widetilde{q_{n-1}}\\ \end{aligned}
第二步：
::: align-center
q_i=\frac{\widetilde{q_i}}{||\widetilde{q_i}||},i=1,2,\cdots,n

:::

现在我们的A是一个这样的矩阵
::: align-center
A=\begin{bmatrix}a_1 & a_2 &\cdots & a_n\end{bmatrix}

:::
我们想要获得空间的一组基，这组基组成的矩阵实际上就是Q
::: align-center
\widetilde{Q}=\begin{bmatrix}\widetilde{q_1} & \widetilde{q_2} & \cdots & \widetilde{q_n}\end{bmatrix}

:::
如何将A化为\widetilde{Q}？，很显然我们需要对A进行列变换，而列变换的过程就是
\begin{aligned} &c_2-\frac{c_1^Tc_2}{c_1^Tc_1}c_1\\ &c_3-\frac{c_1^Tc_3}{c_1^Tc_1}c_1-\frac{c_2^Tc_3}{c_2^Tc_2}c_2\\ &\vdots\\ &c_n-\frac{c_1^Tc_n}{c_1^Tc_1}c_1-\frac{c_2^Tc_n}{c_2^Tc_2}c_2-\frac{c_{n-1}^Tc_n}{c_{n-1}^Tc_{n-1}}c_{n-1}\\ \end{aligned}

我们又知道，列变换可以使用初等矩阵来表示，对于上面属于的列倍加变换它应总是一个上三角矩阵，并且由于初等矩阵可互乘，最后我们得到的这个倍加变换矩阵应为\widetilde{R}^{-1}，它应该就是一个单位上三角矩阵，而变换得到的结果应该就是由未单位化基向量组成的矩阵\widetilde{Q_1}，由初等矩阵可逆，故有
::: align-center
A\widetilde{R}^{-1}=\widetilde{Q}\Rightarrow A=\widetilde{Q}\widetilde{R}(*)

:::
其中
::: align-center
\widetilde{R}=\begin{bmatrix} 1 & \frac{\widetilde{q_1}^Ta_2}{\widetilde{q_1}^Tq_1} &\cdots & \frac{\widetilde{q_1}^Ta_n}{\widetilde{q_1}^Tq_1} \\ & 1 & \ddots & \vdots \\ & & \ddots & \frac{\widetilde{q_{n-1}}^Ta_n}{\widetilde{q_{n-1}}^T\widetilde{q_{n-1}}}\\ & & & 1 \end{bmatrix}

:::
之后我们再对\widetilde{Q}内的矩阵进行单位化，即对\widetilde{Q}内的每个列向量的每个元素都除以这行列向量的范数，对于矩阵来说这相当于对每一列进行倍乘变换，同样地，我们有初等列倍乘变换矩阵diag(\frac{1}{||\widetilde{q_i}||})，并有
::: align-center
\widetilde{Q}diag(\frac{1}{||\widetilde{q_i}||})=Q \Rightarrow \widetilde{Q}=Qdiag(||\widetilde{q_i}||)

:::
代入(*)得
::: align-center
A=Qdiag(||\widetilde{q_i}||)\widetilde{R}=Q(diag(||\widetilde{q_i}||)\widetilde{R})=QR

:::
其中Q为由标准正交基为列组成的矩阵，R为上三角矩阵，可知R“记录了”我们正交化的过程

应试环节

又到了我们喜闻乐见但又恨之入骨的环节，笔者会提供一种不易出错的计算QR分解的方法，也可用于编写程序
我们以4x4QR分解为例

画表

画出一个表

	a_1	a_2	a_3	a_4
q_1^T
q_2^T
q_3^T
\widetilde{q}
取范数
做比

• 表最顶部的行：需要正交化的向量，即A的列向量
• 表最左边的列
  1. 正交化单位化后的基的转置，但要比顶部行数量少1
  2. \widetilde{q}，即未单位化的基
  3. 取范数和做比不做解释

填数

• 将R中的元素按照下面的规则填入表中
  • 对于前几行，r的下标是最左一列q的下标和顶部行a的行标，舍去对角线元素
  • 对于取范数这一行，r的下标是顶部行的行标X2，即对角线元素

	a_1	a_2	a_3	a_4
q_1^T		r_{12}	r_{13}	r_{14}
q_2^T			r_{23}	r_{24}
q_3^T				r_{34}
\widetilde{q}
取范数	r_{11}	r_{22}	r_{33}	r_{44}
做比

运算

• 从每一列开始，由上到下
• 将R中的元素按照下面的规则填入表中
  • 对于前几行，r的下标是最左一列q的下标和顶部行a的行标，舍去对角线元素
  • 对于\widetilde{q}这一行，r的下标是顶部行的行标X2，即对角线元素

从第一列开始

	a_1	a_2	a_3	a_4
q_1^T	\downarrow	r_{12}	r_{13}	r_{14}
q_2^T	\downarrow		r_{23}	r_{24}
q_3^T	\downarrow			r_{34}
\widetilde{q}	\widetilde{q_1}
取范数	r_{11}	r_{22}	r_{33}	r_{44}
做比	q_1

则\widetilde{q_1}=a_1，取范数得r_{11}=||\widetilde{q_1}||，二者做比得q_1=\frac{\widetilde{q_1}}{r_{11}}

第二列

	a_1	a_2	a_3	a_4
q_1^T	\downarrow	r_{12}	r_{13}	r_{14}
q_2^T	\downarrow	\downarrow	r_{23}	r_{24}
q_3^T	\downarrow	\downarrow		r_{34}
\widetilde{q}	\widetilde{q_1}	\widetilde{q_2}
取范数	r_{11}	r_{22}	r_{33}	r_{44}
做比	q_1	q_2

此时\widetilde{q_2}=a_2-q_1^Ta_2，取范数得r_{22}=||\widetilde{q_2}||，二者做比得q_1=\frac{\widetilde{q_2}}{r_{22}}

后几列类似

结果

	a_1	a_2	a_3	a_4
q_1^T	\downarrow	r_{12}	r_{13}	r_{14}
q_2^T	\downarrow	\downarrow	r_{23}	r_{24}
q_3^T	\downarrow	\downarrow	\downarrow	r_{34}
\widetilde{q}	\widetilde{q_1}	\widetilde{q_2}	\widetilde{q_3}	\widetilde{q_4}
取范数	r_{11}	r_{22}	r_{33}	r_{44}
做比	q_1	q_2	q_3	q_4

最后将r和q填入矩阵，就可得A的QR分解

CGS与MGS

我们上面讲的方法都属于简化QR分解中的经典格拉姆-施密特算法(Classic Gram-Schmidt, CGS)
在计算机中使用CGS时，由于每个\widetilde{q_i}都要直接减去所有其他以往方向上的投影，在浮点运算中可能导致误差累计从而导致Q最终的非正交性，由此改进格拉姆-施密特算法(Modified Gram-Schmidt, MGS)应运而生，他与CGS的不同之处在于改变了正交化中的计算次序
\begin{aligned} &\widetilde{q_1}=a_1\\ &\widetilde{q_2}'=a_2-\frac{\widetilde{q_1}^Ta_2}{\widetilde{q_1}^T\widetilde{q_1}}\widetilde{q_1}\\ &\widetilde{q_3}'=a_3-\frac{\widetilde{q_1}^Ta_3}{\widetilde{q_1}^T\widetilde{q_1}}\widetilde{q_1}\\ &\vdots\\ &\widetilde{q_n}'=a_n-\frac{\widetilde{q_1}^Ta_n}{\widetilde{q_1}^T\widetilde{q_1}}\widetilde{q_1}\\ \end{aligned}
此时\widetilde{q_2}到\widetilde{q_n}已经更新，然后再计算
\begin{aligned} &\widetilde{q_2}=\widetilde{q_2}'\\ &\widetilde{q_3}=\widetilde{q_3}'-\frac{\widetilde{q_2}^Ta_3}{\widetilde{q_2}^T\widetilde{q_2}}\widetilde{q_2}\\ &\vdots\\ &\widetilde{q_n}'=\widetilde{q_n}'-\frac{\widetilde{q_2}^Ta_n}{\widetilde{q_2}^T\widetilde{q_2}}\widetilde{q_2}\\ \end{aligned}
如此往复，进行n+(n-1)+\cdots+1=\frac{n(n+1)}{2}次循环后最终得到所有的基向量
因为每次都是先更新以往的向量后再计算新的，可以一定程度避免计算机内部的浮点数舍入误差累计，在计算机内使用简化QR分解时优先选择MGS

完整QR分解（HouseHolder)

反射变换

对于向量x和它的反射H_\theta(x)，存在关系
::: align-center
x-H_\theta(x)=2y

:::
其中y是x在直线的法方向v上的投影，由投影和内积定义知其为
::: align-center
y=(v^Tx)v

:::
整理得
::: align-center
H_\theta(x)=(I-2vv^T)x

:::
故得反射变换得表示矩阵为
::: align-center
I-2vv^T

:::

HouseHolder变换

对于任意\mathbb{R^n}中的向量v，都唯一决定以其为法向量的超平面\mathcal{N}(v^T)（\mathcal{N}指的是零向量的原像集，这里指的是线性映射Av^T=0,若考虑v^T与v正交，则Av^T=0表示的就是一个超平面)

令H_v:I-2vv^T，则H_v是\mathbb{R}^n上的反射变换，反射面是\mathcal{N}(v^T)
事实上，若我们取任意向量\omega，均可有恒等变形
::: align-center
\omega=I_n\omega=(I_n-vv^T)\omega+vv^T\omega

:::
其中vv^T\omega=(v^T\omega)v是\omega在v上的投影向量，它与v共线，而(I_n-vv^T)\omega=\omega-vv^T\omega是\omega减去它与v共线的部分，则这部分与v正交

此时我们考虑
::: align-center
(I_n-2vv^T)\omega=(I_n-vv^T)\omega-vv^T\omega

:::
即\omega在线性变换(I_n-2vv^T)下的像，它等于\omega与v正交的部分加上原来共线部分的反向的向量，可见这是一个反射变换，我们把这种变换叫做HouseHolder变换

::: align-center
给定\mathbb{R}^n中向量x,y，满足||x||=||y||，则存在HouseHolder反射变换H_v，其中v=\frac{y-x}{||y-x||}，使得H_v(x)=y

:::

QR分解与HouseHolder变换

（完整QR分解）对m\times n矩阵A，其中m\ge n，则存在m阶正交矩阵Q和m\times n矩阵R，使得A=QR，其中R=\begin{bmatrix}R_1\\O\end{bmatrix}，Q=\begin{bmatrix}Q_1 & Q_2\end{bmatrix}，即Q的列向量组由Q_1的列向量组扩充而成

对n用数学归纳法，
（i）当n=1时，若a_1\ne 0，则令Q=q_1=\frac{a_1}{||a_1||},R=||a_1||即可，若a_1=0，则任取一个向量作为q_1，令R=0即可

（ii）考虑矩阵A=\begin{bmatrix}a_1&a_2&\cdots&a_n\end{bmatrix}，其中的每个向量都是m维，
（a）若a_1\ne 0
现在考虑有一个HouseHolder变换Q_1，使得Q_1a_1=||a_1||e_1（我们把a_1和它等范数的正交基||a_1||e_1建立反射关系），因此Q_1=(I_n-2vv^T)，其中v=\frac{a_1-||a_1||e_1}{||a_1-||a_1e_1||||}，此时Q_1A=\begin{bmatrix}||a_1|| & b^T \\O& A_1\end{bmatrix}，其中A_1为(m-1)\times(n-1)矩阵，而m-1\ge n-1满足归纳假设，

即存在正交矩阵Q_2和上三角矩阵R_2，使得A_1=Q_2\begin{bmatrix}R_2\\O\end{bmatrix}，此时可再使用HouseHolder变换迭代，同时因为Q_2是正交矩阵，则Q_2^TA_1=\begin{bmatrix}R_2\\O\end{bmatrix}，

那么\begin{bmatrix}1 & \\ & Q_2^T\end{bmatrix}Q_1A=\begin{bmatrix}||a_1|| & b^T\\0 &R_2\\0 & O\end{bmatrix}，令Q=(\begin{bmatrix}1 & \\ & Q_2^T\end{bmatrix}Q_1)^{-1}=Q_1^T\begin{bmatrix}1 & \\ & Q_2^T\end{bmatrix},R=\begin{bmatrix}||a_1|| & b^T\\0 &R_2\end{bmatrix}即满足A=QR

（b）若a_1=a_2=\cdots=a_{j-1}=0但a_{j}\ne 0，则A余下的部分\begin{bmatrix}a_j&\cdots&a_n\end{bmatrix}=QR_j，因此A=\begin{bmatrix}O & QR_j\end{bmatrix}=Q\begin{bmatrix}O & R_j\end{bmatrix}，令R=\begin{bmatrix}O & R_j\end{bmatrix}即得A=QR

此方法相对于简化QR分解计算更为复杂，同时因为保留了矩阵中的零块，在计算机中占用的内存空间也会变大，但是它的好处就是计算的精度比简化QR分解高，因此在计算机中需要高精度QR分解时会采用这种方法