條件數

本文的閱讀等級:高級

當一線性系統受到極微小的擾動即可引發方程解劇烈變化時,我們將無從信任計算結果,便稱它是病態系統 (見“病態系統”)。條件數 (condition number) 是矩陣運算誤差分析的基本工具,它可以度量矩陣對於數值計算的敏感性和穩定性,也可以用來檢定病態系統。本文通過一個簡單的線性方程擾動問題介紹條件數的推導過程,基本推演工具是矩陣範數 \Vert A\Vert 的定義所含的兩個不等式 (見“矩陣範數”):

\Vert A+B\Vert\le\Vert A\Vert+\Vert B\Vert

\Vert AB\Vert\le\Vert A\Vert\cdot\Vert B\Vert

 
問題一:考慮線性方程 A\mathbf{x}=\mathbf{b},其中 An\times n 階可逆矩陣。假設 \mathbf{b} 受到微小擾動成為 \mathbf{b}+\delta\mathbf{b},方程解隨之由 \mathbf{x} 變為 \mathbf{x}+\delta\mathbf{x},試計算相對誤差 \frac{\Vert\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert} 的可能範圍。

 
考慮包含擾動及誤差的線性方程:

A(\mathbf{x}+\delta\mathbf{x})=\mathbf{b}+\delta\mathbf{b}

與原方程式 A\mathbf{x}=\mathbf{b} 相減可得 A(\delta\mathbf{x})=\delta\mathbf{b},方程解的誤差為 \delta\mathbf{x}=A^{-1}(\delta\mathbf{b})。利用不等式 \Vert\delta\mathbf{x}\Vert=\Vert A^{-1}(\delta\mathbf{b})\Vert\le\Vert A^{-1}\Vert\cdot\Vert\delta\mathbf{b}\Vert\Vert\mathbf{b}\Vert=\Vert A\mathbf{x}\Vert\le\Vert A\Vert\cdot\Vert\mathbf{x}\Vert,可得

\displaystyle\frac{\Vert\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert}=\frac{\Vert A^{-1}(\delta\mathbf{b})\Vert}{\Vert\mathbf{x}\Vert}\le\frac{\Vert A^{-1}\Vert\cdot\Vert\delta\mathbf{b}\Vert}{\Vert\mathbf{x}\Vert}\le\left(\Vert A\Vert\cdot\Vert A^{-1}\Vert\right)\frac{ \Vert\delta\mathbf{b}\Vert}{\Vert\mathbf{b}\Vert}

另一方面,我們考慮不等式 \Vert\delta\mathbf{b}\Vert=\Vert A(\delta\mathbf{x})\Vert\le\Vert A\Vert\cdot\Vert\delta\mathbf{x}\Vert\Vert\mathbf{x}\Vert=\Vert A^{-1}\mathbf{b}\Vert\le\Vert A^{-1}\Vert\cdot\Vert\mathbf{b}\Vert,就有

\displaystyle \frac{\Vert\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert}\ge\frac{\Vert\delta\mathbf{b}\Vert}{\Vert A\Vert\cdot\Vert\mathbf{x}\Vert}\ge\left(\Vert A\Vert\cdot\Vert A^{-1}\Vert\right)^{-1}\frac{\Vert\delta\mathbf{b}\Vert}{\Vert\mathbf{b}\Vert}

上面二式指出相對誤差 \frac{\Vert\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert} 的上下界由矩陣範數乘積 \Vert A\Vert \cdot\Vert A^{-1}\Vert 決定,根據這個結果我們定義方陣 A 的條件數為

\kappa(A)=\Vert A\Vert\cdot\Vert A^{-1}\Vert

也就是說,\delta\mathbf{b} 引起的相對誤差大小由內部因素 \kappa(A) 和外部因素 \frac{\Vert\delta\mathbf{b}\Vert}{\Vert\mathbf{b}\Vert} 共同決定:

\displaystyle \kappa(A)^{-1}\frac{\Vert\delta\mathbf{b}\Vert}{\Vert\mathbf{b}\Vert}\le\frac{\Vert\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert}\le\kappa(A)\frac{\Vert\delta\mathbf{b}\Vert}{\Vert\mathbf{b}\Vert}

 
接著討論條件數的計算方法。設 n\times n 階矩陣 A 的奇異值分解為 A=U\Sigma V^T,其中 \Sigma=\mathrm{diag}(\sigma_1,\ldots,\sigma_n)\sigma_1\ge\cdots\ge\sigma_n\ge 0UV 為正交矩陣 (見“奇異值分解 (SVD)”)。若 A 為可逆矩陣,\mathrm{rank}A=n,所以 \sigma_i>0i=1,\ldots,n。因為正交矩陣 UV 不改變向量長度,最大奇異值和最小奇異值分別滿足

\displaystyle \sigma_1=\max_{\Vert\mathbf{x}\Vert=1}\Vert\Sigma\mathbf{x}\Vert=\max_{\Vert\mathbf{x}\Vert=1}\Vert A\mathbf{x}\Vert=\Vert A\Vert_2

\displaystyle \sigma_n=\min_{\Vert\mathbf{x}\Vert=1}\Vert\Sigma\mathbf{x}\Vert=\min_{\Vert\mathbf{x}\Vert=1}\Vert A\mathbf{x}\Vert=\frac{1}{\Vert A^{-1}\Vert_2}

以2-範數 (2-norm) 計算的條件數也就為

\displaystyle \kappa(A)=\Vert A\Vert_2\cdot\Vert A^{-1}\Vert_2=\frac{\sigma_1}{\sigma_n}\ge 1

A 為對稱可逆矩陣,計算

A^2=AA^T=(U\Sigma V^T)(V\Sigma U^T)=U\Sigma^2U^T

因此 \sigma_i=\vert\lambda_i\verti=1,\ldots,n\lambda_iA 的特徵值,條件數可表示為

\displaystyle \kappa(A)=\left|\frac{\lambda_1}{\lambda_n}\right|

又若 A 為正交矩陣,A^T=A^{-1},則 \sigma_1=\cdots=\sigma_n=1,故 \kappa(A)=1,這說明正交矩陣具有最佳的數值穩定性。

 
下面我們進一步討論奇異值分解與條件數上下界的關係。奇異值分解的 U 矩陣其行向量 \mathbf{u}_1,\ldots,\mathbf{u}_n 稱為左奇異向量,V 矩陣的行向量 \mathbf{v}_1,\ldots,\mathbf{v}_n 稱為右奇異向量,左右奇異向量的關係如下:

A\mathbf{v}_j=\sigma_j\mathbf{u}_j,~~~j=1,\ldots,n

考慮 \mathbf{b}=\alpha\mathbf{u}_1\delta\mathbf{b}=\beta\mathbf{u}_n,方程解和誤差可分別表示為

\displaystyle \mathbf{x}=A^{-1}\mathbf{b}=A^{-1}(\alpha\mathbf{u}_1)=\frac{\alpha}{\sigma_1}\mathbf{v}_1

\displaystyle \delta\mathbf{x}=A^{-1}(\delta\mathbf{b})=A^{-1}(\beta\mathbf{u}_n)=\frac{\beta}{\sigma_n}\mathbf{v}_n

計算向量長度並相除,

\displaystyle\frac{\Vert\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert}=\left(\frac{\sigma_1}{\sigma_n}\right)\frac{\vert\beta\vert}{\vert\alpha\vert}=\kappa(A)\frac{\Vert\delta\mathbf{b}\Vert}{\Vert\mathbf{b}\Vert}

得知 \mathbf{b}=\alpha\mathbf{u}_1\delta\mathbf{b}=\beta\mathbf{u}_n 滿足相對誤差的上界。同樣方式也可以證明 \mathbf{b}=\alpha\mathbf{u}_n\delta\mathbf{b}=\beta\mathbf{u}_1 滿足相對誤差的下界,亦即

\displaystyle\frac{\Vert\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert}=\left(\frac{\sigma_n}{\sigma_1}\right)\frac{\vert\beta\vert}{\vert\alpha\vert}=\kappa(A)^{-1}\frac{\Vert\delta\mathbf{b}\Vert}{\Vert\mathbf{b}\Vert}

 
考慮下面兩個取自“病態系統”一文的例子:

A=\begin{bmatrix}  .540&.387\\  .647&.323  \end{bmatrix},~B=\begin{bmatrix}  .540&.323\\  .647&.387  \end{bmatrix}

計算奇異值後代入條件數公式:

\begin{aligned} \kappa(A)&=.978920/.077605=12.614\\  \kappa(B)&=.981991/.000001=981,991\approx 10^{6}\end{aligned}

此例中,\kappa(A) 小 (指 10^{\gamma} 的冪 \gamma),A 是良置的,小擾動 \delta\mathbf{b} 不至於對方程解造成大誤差。但 \kappa(B) 很大,B 是病態的,小擾動 \delta\mathbf{b} 可能引起大 \delta\mathbf{x},但大擾動 \delta\mathbf{b} 也可能產生微小的 \delta\mathbf{x}。由於難以掌握 \delta\mathbf{b} 的變化方向,在面對病態系統時我們總是要預防最壞的情形發生。

 
問題二:如果線性方程 A\mathbf{x}=\mathbf{b} 的係數矩陣 A 和常數向量 \mathbf{b} 都受到雜音干擾,這對方程解 \mathbf{x} 的誤差有何影響?

 
考慮 A\mathbf{b} 分別受到 \delta A\delta\mathbf{b} 干擾,則 \delta\mathbf{x} 滿足

(A+\delta A)(\mathbf{x}+\delta\mathbf{x})=\mathbf{b}+\delta\mathbf{b}

為簡化問題,以下假設 \frac{\Vert\delta A\Vert}{\Vert A\Vert}\le\epsilon\frac{\Vert\delta\mathbf{b}\Vert}{\Vert\mathbf{b}\Vert}\le\epsilon。令 k=\epsilon\kappa(A),當 \epsilon<\frac{1}{\kappa(A)} 時,下面三個性質成立。

 
性質一A+\delta A 為可逆矩陣。

已知 A 是可逆的,由假設條件可得

\Vert A^{-1}(\delta A)\Vert\le\Vert A^{-1}\Vert\cdot\Vert\delta A\Vert\le\epsilon\Vert A^{-1}\Vert\cdot\Vert A\Vert=k<1

根據 Neumann 無窮級數性質可知 I+A^{-1}(\delta A) 是可逆的,且 \Vert(I+A^{-1}(\delta A))^{-1}\Vert\le\frac{1}{1-k} (見“Neumann 無窮級數”),故 A(I+A^{-1}(\delta A))= A+\delta A 亦為可逆矩陣。

 
性質二\displaystyle\frac{\Vert\mathbf{x}+\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert}\le\frac{1+k}{1-k}

(A+\delta A)(\mathbf{x}+\delta\mathbf{x})=\mathbf{b}+\delta\mathbf{b} 等號兩邊左乘 A^{-1}

(I+A^{-1}(\delta A))(\mathbf{x}+\delta\mathbf{x})=\mathbf{x}+A^{-1}(\delta\mathbf{b})

性質一指出 I+A^{-1}(\delta A) 是可逆的,就有

\mathbf{x}+\delta\mathbf{x}=( I+A^{-1}(\delta A))^{-1}(\mathbf{x}+A^{-1}(\delta\mathbf{b}))

利用性質一的不等式,可得

\displaystyle \Vert\mathbf{x}+\delta\mathbf{x}\Vert\le\Vert(I+A^{-1}(\delta A))^{-1}\Vert\cdot(\Vert\mathbf{x}\Vert+\epsilon\Vert A^{-1}\Vert\cdot\Vert\mathbf{b}\Vert)\le\frac{1}{1-k}\left(\Vert\mathbf{x}\Vert+k\frac{\Vert\mathbf{b}\Vert}{\Vert A\Vert}\right)

\Vert\mathbf{b}\Vert=\Vert A\mathbf{x}\Vert\le\Vert A\Vert\cdot\Vert\mathbf{x}\Vert,所以

\displaystyle \Vert\mathbf{x}+\delta\mathbf{x}\Vert\le\frac{1}{1-k}(\Vert\mathbf{x}\Vert+k\Vert\mathbf{x}\Vert)=\frac{1+k}{1-k}\Vert\mathbf{x}\Vert

 
性質三\displaystyle\frac{\Vert\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert}\le\frac{2\epsilon}{1-k}\kappa(A)

性質二的關係式可寫為 \delta\mathbf{x}=A^{-1}(\delta\mathbf{b})-A^{-1}(\delta A)(\mathbf{x}+\delta\mathbf{x}),計算 \delta\mathbf{x} 長度:

\Vert\delta\mathbf{x}\Vert\le\Vert A^{-1}\Vert\cdot\Vert\delta\mathbf{b}-(\delta A)(\mathbf{x}+\delta\mathbf{x})\Vert

利用已知條件與三角不等式,可得

\Vert\delta\mathbf{x}\Vert\le\epsilon\Vert A^{-1}\Vert\cdot\Vert\mathbf{b}\Vert+\epsilon\Vert A^{-1}\Vert\cdot\Vert A\Vert\cdot\Vert\mathbf{x}+\delta\mathbf{x}\Vert

將上式除以 \Vert\mathbf{x}\Vert 並利用性質二,

\displaystyle\frac{\Vert\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert}\le\epsilon\kappa(A)\frac{\Vert\mathbf{b}\Vert}{\Vert A\Vert\cdot\Vert\mathbf{x}\Vert}+\epsilon\kappa(A)\frac{\Vert\mathbf{x}+\delta\mathbf{x}\Vert}{\Vert\mathbf{x}\Vert}\le\epsilon\kappa(A)\left(1+\frac{1+k}{1-k}\right)=\frac{2\epsilon}{1-k}\kappa(A)

 
方程解的相對誤差來自 A\mathbf{b} 的相對擾動量和 2\epsilon,條件數 \kappa(A) 再將此誤差放大。換句話說,執行數值計算時,因捨入誤差而造成的方程解誤差有兩個來源:一是線性系統自身的敏感性,以 \kappa(A) 表示,另一則是真實誤差 \delta A\delta\mathbf{b}。當我們用 LU 分解計算時,受限於電腦的精準度,實際得到的是近似矩陣 \tilde{L}\tilde{U},因此無可避免地使用了錯誤矩陣 A+\delta A=\tilde{L}\tilde{U}。英國數學家威爾金森 (James H. Wilkinson) 證明部分軸元法 (partial pivoting,見“PA=LU 分解”) 確實能夠有效控制 \delta A,故捨入誤差的主要決定因素為係數矩陣的條件數。感謝威爾金森提供的定心丸,除非碰上病態系統,否則我們大可放心地使用高斯消去法。

 
本文參考:
Gene H. Golub, Charles F. Van Loan,Matrix Computations,2nd ed.,1989。

Advertisements
本篇發表於 線性代數專欄, 數值線性代數 並標籤為 , , , , , , , , 。將永久鏈結加入書籤。

2 則回應給 條件數

  1. edge 說道:

    可以請老師提供partial pivoting的證明嗎??
    需要這顆定心丸

發表迴響

在下方填入你的資料或按右方圖示以社群網站登入:

WordPress.com Logo

您的留言將使用 WordPress.com 帳號。 登出 / 變更 )

Twitter picture

您的留言將使用 Twitter 帳號。 登出 / 變更 )

Facebook照片

您的留言將使用 Facebook 帳號。 登出 / 變更 )

Google+ photo

您的留言將使用 Google+ 帳號。 登出 / 變更 )

連結到 %s