矩陣的 QR 分解在電腦運算中可能造成誤差，本文探討一下一種改進版本的 Gram-Schmidt 正交化方法。

經典的 Gram-Schmidt 方法可能造成數值不穩定性。在電腦中，舍入誤差可能會累積，造成得到的正交基並不具有正交性。

經典的 Gram-Schmidt

我們先來回顧一下經典的 Gram-Schmidt 方法：

1
2
3
4
5
6

for j = 1 : n
    v_j = x_j
    for k = 1 : j - 1
        v_j = v_j - ( (v_k^T x_j) / (v_k^T v_k) ) * v_k
    endfor
endfor

最終我們將 vjv_j 歸一化得到標準正交基 qj=(vj)/(∣vj∣)q_j = (v_j)/(| v_j |) 。

注意：CGS 始終使用原始向量 xjx_j 與之前的基向量計算投影

數學證明

簡單進行數學證明

定義第 jj 步生成的向量 vjv_j

vj=xj−∑k=1j−1((vkTxj)/(vkTvk))vkv_j = x_j - \sum_{k=1}^{j-1} ( (v_k^T x_j) / (v_k^T v_k) ) v_k

選取任意一個之前的基向量 vmv_m（其中 1≤m<j1 \le m < j），計算它與 vjv_j 的內積 vmTvjv_m^T v_j：

vmTvj=vmT(xj−∑k=1j−1((vkTxj)/(vkTvk))vk)v_m^T v_j = v_m^T ( x_j - \sum_{k=1}^{j-1} ((v_k^T x_j) / (v_k^T v_k)) v_k )

利用內積的線性性質，將 vmTv_m^T 乘進去：

vmTvj=vmTxj−∑k=1j−1((vkTxj)/(vkTvk))(vmTvk)v_m^T v_j = v_m^T x_j - \sum_{k=1}^{j-1} ((v_k^T x_j) / (v_k^T v_k)) (v_m^T v_k)

由於前提假設 v1,…,vj−1v_1, \dots, v_{j-1} 是兩兩正交的，所以在求和符號 ∑\sum 中：

• 當 k≠mk \neq m 時，vmTvk=0v_m^T v_k = 0 。
• 當且僅當 k=mk = m 時，vmTvk=vmTvm≠0v_m^T v_k = v_m^T v_m \neq 0 。

因此，求和項中只剩下 k=mk=m 這一項：

vmTvj=vmTxj−((vmTxj)/(vmTvm))(vmTvm)v_m^T v_j = v_m^T x_j - ((v_m^T x_j) / (v_m^T v_m)) (v_m^T v_m)

分子分母中的標量 vmTvmv_m^T v_m 相互抵消：

vmTvj=vmTxj−vmTxjv_m^T v_j = v_m^T x_j - v_m^T x_j

vmTvj=0v_m^T v_j = 0

證畢

誤差分析

如果在 CGS（經典格拉姆-施密特正交化）中計算 q2q_2 時發生誤差，導致 q1Tq2=δq_1^T q_2 = \delta 是一個很小但非零的數值，這個誤差將不會在隨後的任何計算中被修正：

v3=x3−(q1Tx3)q1−(q2Tx3)q2v_3 = x_3 - (q_1^T x_3)q_1 - (q_2^T x_3)q_2

q2Tv3=q2Tx3−(q1Tx3)δ−(q2Tx3)=−(q1Tx3)δq_2^T v_3 = q_2^T x_3 - (q_1^T x_3)\delta - (q_2^T x_3) = -(q_1^T x_3)\delta

q1Tv3=q1Tx3−(q1Tx3)−(q2Tx3)δ=−(q2Tx3)δq_1^T v_3 = q_1^T x_3 - (q_1^T x_3) - (q_2^T x_3)\delta = -(q_2^T x_3)\delta

我們可以看出 v3v_3 與 q1q_1 或 q2q_2 都不正交。

改進的 Gram-Schmidt (MGS)

1
2
3
4
5
6
7
8
9
10

for j = 1 : n
    v_j = x_j
endfor

for j = 1 : n
    q_j = v_j / ||v_j||_2
    for k = j + 1 : n
        v_k = v_k - (q_j^T v_k) q_j
    endfor
endfor

最終我們將 vjv_j 歸一化得到標準正交基 qj=(vj)/(∣vj∣)q_j = (v_j)/(| v_j |) 。

區別在於，一旦算出來一個向量 qq，立即用它去更新後面所有的向量（減去向量在 qq 上的投影），使得後面所有的向量與這個向量 qq 正交。

誤差分析

假設 MGS 出現：q1Tq2=δq_1^T q_2 = \delta 很小但非零。

對於第三個向量 v3v_3 ：

• 初始狀態：v3(0)=x3v_3^{(0)} = x_3
• j=1:v3(1)=v3(0)−(q1Tv3(0))q1j = 1: v_3^{(1)} = v_3^{(0)} - (q_1^T v_3^{(0)})q_1
• j=2:v3=v3(1)−(q2Tv3(1))q2j = 2: v_3 = v_3^{(1)} - (q_2^T v_3^{(1)})q_2

此時 v3v_3 是第三個向量在歸一化之前的最終形式。

讓我們檢查正交性，假設不再產生其他誤差：

q2Tv3=q2Tv3(1)−(q2Tv3(1))=0q_2^T v_3 = q_2^T v_3^{(1)} - (q_2^T v_3^{(1)}) = 0

所以，我們保留了對 q2q_2 的正交性。

接下來看 q1q_1：

q1Tv3=q1Tv3(1)−(q2Tv3(1))δq_1^T v_3 = q_1^T v_3^{(1)} - (q_2^T v_3^{(1)})\delta

q2Tv3(1)=q2Tv3(0)−(q1Tv3(0))δq_2^T v_3^{(1)} = q_2^T v_3^{(0)} - (q_1^T v_3^{(0)})\delta

q1Tv3(1)=q1Tv3(0)−(q1Tv3(0))=0q_1^T v_3^{(1)} = q_1^T v_3^{(0)} - (q_1^T v_3^{(0)}) = 0

由此可得：

q1Tv3=−(q2Tv3(0)−(q1Tv3(0))δ)δ=−q2Tv3(0)δ+q1Tv3(0)δ2q_1^T v_3 = -(q_2^T v_3^{(0)} - (q_1^T v_3^{(0)})\delta)\delta= -q_2^T v_3^{(0)}\delta + q_1^T v_3^{(0)}\delta^2

由於 δ\delta 非常小，δ2\delta^2 就更小了。因此 q1Tv3q_1^T v_3 中的誤差和 CGS 差不多，但我們消除了 q2Tv3q_2^T v_3 中的誤差，使得誤差更小。

總結

在電腦中，由於浮點運算的特殊性，有許多數學方法需要重新考慮，儘量提高數值穩定性。這使得我們有必要重新審視我們學習過的數學知識，針對電腦的特殊性進行重新設計與優化。