約束最小平方問題

本文的閱讀等級:中級

A 為一個 m\times n 階實矩陣,\mathbf{b}\in\mathbb{R}^m。如果線性方程 A\mathbf{x}=\mathbf{b} 是不一致的 (即不存在解),實務的作法是將線性方程問題改為最小平方近似問題:

\displaystyle  \min_{\mathbf{x}}\Vert A\mathbf{x}-\mathbf{b}\Vert^2

其中 \Vert A\mathbf{x}-\mathbf{b}\Vert 是2-範數 (見“向量範數”),即 A\mathbf{x}\mathbf{b} 的歐幾里得距離。根據正交原則,最小平方解 \hat{\mathbf{x}} 滿足正規方程 (normal equation) A^TA\hat{\mathbf{x}}=A^T\mathbf{b} (見“從線性變換解釋最小平方近似”)。若 \hbox{rank}A=n,也就是說 A 的行向量 (column vector) 構成一個線性獨立集合,則存在唯一的最小平方解

\displaystyle  \hat{\mathbf{x}}=(A^TA)^{-1}A^T\mathbf{b}

如果最小平方解必須滿足某些束縛條件,則稱為約束最小平方問題 (constrained least-squares problem)。本文討論兩種常出現在多種應用場合的約束形式。線性約束最小平方問題是指束縛條件為線性方程[1]

\displaystyle  \min_{C\mathbf{x}=\mathbf{d}}\Vert A\mathbf{x}-\mathbf{b}\Vert^2

其中 C 是一個 p\times n 階實矩陣,\mathbf{d}\in\mathbb{R}^p。正則 (regularized) 最小平方問題限制未知向量的長度必須固定:

\displaystyle  \min_{\Vert\mathbf{x}\Vert=d}\Vert A\mathbf{x}-\mathbf{b}\Vert^2

 
線性約束最小平方問題

考慮下面的問題:

\displaystyle  \min_{C\mathbf{x}=\mathbf{d}}\Vert A\mathbf{x}-\mathbf{b}\Vert^2

其中 A 是一個 m\times n 階實矩陣,C 是一個 p\times n 階實矩陣,\mathbf{b}\in\mathbb{R}^m\mathbf{d}\in\mathbb{R}^p。我們要求限制條件數不大於未知向量的維數,即 p\le n。極小範數問題 (minimum norm problem) 是一種特殊的線性約束最小平方問題:

\displaystyle  \min_{C\mathbf{x}=\mathbf{d}}\Vert \mathbf{x}\Vert^2

\hbox{rank}C=p,則存在唯一的極小範數解 (見“極小範數解”)

\displaystyle  \hat{\mathbf{x}}=C^T(CC^T)^{-1}\mathbf{d}

底下介紹一般的線性約束最小平方問題的解法。

 
Lagrange 乘數法是最常被使用的約束優化問題解法 (見“Lagrange 乘數法”)。定義 Lagrangian 函數

L(\mathbf{x},\boldsymbol{\alpha})=\Vert A\mathbf{x}-\mathbf{b}\Vert^2+2\boldsymbol{\alpha}^T(C\mathbf{x}-\mathbf{d})

其中 \boldsymbol{\alpha}=(\alpha_1,\ldots,\alpha_p)^T 是 Lagrange 乘數組成的向量。展開 Lagrangian 函數,

\displaystyle\begin{aligned}  L(\mathbf{x},\boldsymbol{\alpha})&=  (A\mathbf{x}-\mathbf{b})^T(A\mathbf{x}-\mathbf{b})+2\boldsymbol{\alpha}^T(C\mathbf{x}-\mathbf{d})\\  &=\mathbf{x}^TA^TA\mathbf{x}-2\mathbf{x}^TA^T\mathbf{b}+\mathbf{b}^T\mathbf{b}+2\boldsymbol{\alpha}^T(C\mathbf{x}-\mathbf{d}),\end{aligned}

接著計算偏導數 (見“矩陣導數”),

\displaystyle  \begin{aligned}  \frac{\partial L}{\partial \mathbf{x}}&=2A^TA\mathbf{x}-2A^T\mathbf{b}+2C^T\boldsymbol{\alpha}\\  \frac{\partial L}{\partial\boldsymbol{\alpha}}&=2(C\mathbf{x}-\mathbf{d}).  \end{aligned}

令上面兩式等於零可得優化條件式:

\begin{bmatrix}  A^TA&C^T\\  C&0  \end{bmatrix}\begin{bmatrix}  \mathbf{x}\\  \boldsymbol{\alpha}  \end{bmatrix}=\begin{bmatrix}  A^T\mathbf{b}\\  \mathbf{d}  \end{bmatrix}

我們可以證明 (n+p)\times (n+p) 階係數矩陣 \begin{bmatrix}  A^TA&C^T\\  C&0  \end{bmatrix} 是可逆的一個充分條件為 \hbox{rank}A=n\hbox{rank}C=p。考慮齊次方程

\displaystyle  \begin{bmatrix}  A^TA&C^T\\  C&0  \end{bmatrix}\begin{bmatrix}  \mathbf{x}\\  \mathbf{y}  \end{bmatrix}=\mathbf{0}

A^TA\mathbf{x}+C^T\mathbf{y}=\mathbf{0}C\mathbf{x}=\mathbf{0}。第一式左乘 \mathbf{x}^T 並使用第二式,可得

0=\mathbf{x}^TA^TA\mathbf{x}+\mathbf{x}^TC^T\mathbf{y}=(A\mathbf{x})^T(A\mathbf{x})+(C\mathbf{x})^T\mathbf{y}=\Vert A\mathbf{x}\Vert^2

A\mathbf{x}=\mathbf{0} 並有 C^T\mathbf{y}=-A^TA\mathbf{x}=\mathbf{0}。若 \hbox{rank}A=n\hbox{rank}C=\hbox{rank}C^T=p,推論 \mathbf{x}=\mathbf{0}\mathbf{y}=\mathbf{0},即證明 \begin{bmatrix}  A^TA&C^T\\  C&0  \end{bmatrix} 是一個可逆矩陣。在此情況下,解開優化條件式可得唯一的約束最小平方解 \hat{\mathbf{x}},不過這個方法的主要缺點在於 A^TA 容易引入誤差。

 
除了以 Lagrangian 函數推導,另一個作法使用 QR 分解 (見“線代膠囊──QR 分解”)。為清楚說明,我們假設 \hbox{rank}A=n\hbox{rank}C=pn\ge p。設 n\times p 階矩陣 C^T 的 QR 分解為

C^T=Q\begin{bmatrix}  R\\  0  \end{bmatrix}

其中 Q 是一個 n\times n 階實正交矩陣,Q^T=Q^{-1}Rp\times p 階上三角可逆矩陣。令

AQ=\begin{bmatrix}  A_1&A_2  \end{bmatrix},~~Q^T\mathbf{x}=\begin{bmatrix}  \mathbf{x}_1\\  \mathbf{x}_2  \end{bmatrix}

其中 A_1m\times p 階,A_2m\times (n-p) 階,\mathbf{x}_1\in\mathbb{R}^p\mathbf{x}_2\in\mathbb{R}^{n-p}。寫出

A\mathbf{x}-\mathbf{b}=AQQ^T\mathbf{x}-\mathbf{b}=\begin{bmatrix}  A_1&A_2  \end{bmatrix}\begin{bmatrix}  \mathbf{x}_1\\  \mathbf{x}_2  \end{bmatrix}-\mathbf{b}=A_1\mathbf{x}_1+A_2\mathbf{x}_2-\mathbf{b}

以及

C\mathbf{x}=\begin{bmatrix}  R^T&0  \end{bmatrix}Q^T\mathbf{x}=\begin{bmatrix}  R^T&0  \end{bmatrix}\begin{bmatrix}  \mathbf{x}_1\\  \mathbf{x}_2  \end{bmatrix}=R^T\mathbf{x}_1

原線性約束最小平方問題等價於

\displaystyle  \min_{R^T\mathbf{x}_1=\mathbf{d}}\left\Vert A_1\mathbf{x}_1+A_2\mathbf{x}_2-\mathbf{b}\right\Vert^2

因此,\mathbf{x}_1 由一致的線性方程 R^T\mathbf{x}_1=\mathbf{d} 決定,而 \mathbf{x}_2 則為下列無約束最小平方解:

\displaystyle  \min_{\mathbf{x}_2}\left\Vert A_2\mathbf{x}_2-(\mathbf{b}-A_1\mathbf{x}_1)\right\Vert^2

合併以上結果可得 \hat{\mathbf{x}}=Q\begin{bmatrix}  \mathbf{x}_1\\  \mathbf{x}_2  \end{bmatrix}

 
附帶一提,約束優化問題還可用無約束優化問題近似。考慮

\displaystyle  \min_{\mathbf{x}}\Vert A\mathbf{x}-\mathbf{b}\Vert^2+\beta^2\Vert C\mathbf{x}-\mathbf{d}\Vert^2

其中 \Vert C\mathbf{x}-\mathbf{d}\Vert^2 稱為懲罰項,\beta 是一個很大的數,稱為懲罰係數。上式可改寫為標準的最小平方問題

\displaystyle  \min_{\mathbf{x}}\left\Vert \begin{bmatrix}  A\\  \beta C  \end{bmatrix}\mathbf{x}-\begin{bmatrix}  \mathbf{b}\\  \beta\mathbf{d}  \end{bmatrix}\right\Vert^2

\beta\to\infty\Vert C\mathbf{x}-\mathbf{b}\Vert\to 0,故最小平方解趨於受限的最小平方解。看一個例子[2],線性約束最小平方問題

\displaystyle  \min_{x_1=x_2}\left\Vert\begin{bmatrix}  1&2\\  3&4\\  5&6  \end{bmatrix}\begin{bmatrix}  x_1\\  x_2  \end{bmatrix}-\begin{bmatrix}  7\\  1\\  3  \end{bmatrix}\right\Vert^2

的最佳解為 \hat{\mathbf{x}}=(0.3407821,0.3407821)^T。如果以下列最小平方問題近似

\displaystyle  \min_{\mathbf{x}}\left\Vert\left[\!\!\begin{array}{rr}  1&2\\  3&4\\  5&6\\  1000&-1000  \end{array}\!\!\right]\begin{bmatrix}  x_1\\  x_2  \end{bmatrix}-\begin{bmatrix}  7\\  1\\  3\\  0  \end{bmatrix}\right\Vert^2

最佳解則為 \hat{\mathbf{x}}'=(0.3407810,0.3407829)^T

 
正則最小平方問題

給定一個 m\times n 階實矩陣 A\hbox{rank}A=n,以及 \mathbf{b}\in\mathbb{R}^m,考慮這個問題:

\displaystyle  \min_{\Vert\mathbf{x}\Vert=d}\Vert A\mathbf{x}-\mathbf{b}\Vert^2

下面分開兩種情況討論:\mathbf{b}=\mathbf{0}\mathbf{b}\neq\mathbf{0}。若 \mathbf{b}=\mathbf{0},在不失一般性的原則下,正則最小平方問題可表示為

\displaystyle  \min_{\Vert\mathbf{x}\Vert=1}\Vert A\mathbf{x}\Vert^2

A 的奇異值分解為 A=U\Sigma V^T,其中 Um\times m 階實正交矩陣,U^T=U^{-1}Vn\times n 階實正交矩陣,V^T=V^{-1}\Sigma=\begin{bmatrix}  D\\  0  \end{bmatrix}m\times n 階矩陣,D=\hbox{diag}(\sigma_1,\ldots,\sigma_n)\sigma_1\ge\cdots\ge\sigma_n>0 (見“奇異值分解的幾何意義”)。令 \mathbf{y}=V^T\mathbf{x}。因此,\Vert\mathbf{y}\Vert^2=\mathbf{x}^TVV^T\mathbf{x}=\Vert\mathbf{x}\Vert^2=1\Vert A\mathbf{x}\Vert^2=\Vert U\Sigma V^T\mathbf{x}\Vert^2=\Vert \Sigma \mathbf{y}\Vert^2。我們得到約化的等價問題:

\displaystyle  \min_{\Vert\mathbf{y}\Vert=1}\Vert \Sigma\mathbf{y}\Vert^2=\min_{y_1^2+\cdots+y_n^2=1}(\sigma_1y_1)^2+\cdots+(\sigma_ny_n)^2

明顯地,\mathbf{y}=\mathbf{e}_n=(0,\ldots,0,1)^T 最小化目標函數。寫出正交矩陣 V 的行向量表示式 V=\begin{bmatrix}  \mathbf{v}_1&\cdots&\mathbf{v}_n  \end{bmatrix},正則最小平方解為 V 的最末一個行向量,\hat{\mathbf{x}}=V\mathbf{e}_n=\mathbf{v}_n。在此情形下,最小的目標函數值為 \Vert A\hat{\mathbf{x}}\Vert^2=\Vert\Sigma \mathbf{e}_n\Vert^2=\sigma_n^2

 
假設 \mathbf{b}\neq\mathbf{0}。定義 Lagrangian 函數

L(\mathbf{x},\lambda)=\Vert A\mathbf{x}-\mathbf{b}\Vert^2+\lambda(\Vert\mathbf{x}\Vert^2-d^2)

其中 \lambda 是 Lagrange 乘數。將上式寫為

\displaystyle\begin{aligned}  L(\mathbf{x},\lambda)&=(A\mathbf{x}-\mathbf{b})^T(A\mathbf{x}-\mathbf{b})+\lambda(\mathbf{x}^T\mathbf{x}-d^2)\\  &=\mathbf{x}^TA^TA\mathbf{x}-2\mathbf{x}^TA^T\mathbf{b}+\mathbf{b}^T\mathbf{b}^2+\lambda(\mathbf{x}^T\mathbf{x}-d^2).  \end{aligned}

求出偏導數並設為零,

\displaystyle  \frac{\partial L}{\partial\mathbf{x}}=2A^TA\mathbf{x}-2A^T\mathbf{b}+2\lambda\mathbf{x}=\mathbf{0}

可得優化條件式[3]

(A^TA+\lambda I)\mathbf{x}=A^T\mathbf{b}

如前,代入 A 的奇異值分解 A=U\Sigma V^T,則

A^TA+\lambda I=V\Sigma^T U^TU\Sigma V^T+\lambda I=V(\Sigma^T\Sigma +\lambda I)V^T=V(D^2+\lambda I)V^T

A^T\mathbf{b}=V\begin{bmatrix}  D&0  \end{bmatrix}U^T\mathbf{b}。設 \mathbf{c}=(c_1,\ldots,c_m)^T=U^T\mathbf{b},優化條件式左乘 V^T,就得到簡化條件式

\displaystyle  (D^2+\lambda I)V^T\mathbf{x}=\begin{bmatrix}  D&0  \end{bmatrix}\mathbf{c}=\begin{bmatrix}  \sigma_1c_1\\  \vdots\\  \sigma_nc_n  \end{bmatrix}

若係數矩陣 D^2+\lambda I 是可逆的,則存在唯一的正則最小平方解

\displaystyle  \hat{\mathbf{x}}=V(D^2+\lambda I)^{-1}\begin{bmatrix}  \sigma_1c_1\\  \vdots\\  \sigma_nc_n  \end{bmatrix}=V\begin{bmatrix}  \frac{\sigma_1c_1}{\sigma_1^2+\hat\lambda}\\  \vdots\\  \frac{\sigma_nc_n}{\sigma_n^2+\hat\lambda}  \end{bmatrix}=\sum_{i=1}^n\frac{\sigma_ic_i}{\sigma_i^2+\lambda}\mathbf{v}_i

接下來要找出 \hat\lambda 使滿足 \Vert\hat{\mathbf{x}}\Vert^2=d^2。使用上面結果,寫出

\displaystyle  d^2=\Vert\hat{\mathbf{x}}\Vert^2=\Vert V^T\hat{\mathbf{x}}\Vert^2=\left\Vert\begin{bmatrix}  \frac{\sigma_1c_1}{\sigma_1^2+\lambda}\\  \vdots\\  \frac{\sigma_nc_n}{\sigma_n^2+\lambda}  \end{bmatrix}\right\Vert^2=\sum_{i=1}^n\left(\frac{\sigma_ic_i}{\sigma_i^2+\lambda}\right)^2

稱為特徵方程 (secular equation)。定義

\displaystyle  \phi(\lambda)=\sum_{i=1}^n\left(\frac{\sigma_ic_i}{\sigma_i^2+\lambda}\right)^2

\lambda>0\phi(\lambda) 是一個單調遞減函數。若 \phi(0)>d^2,則存在唯一的 \hat\lambda>0 滿足 \phi(\hat\lambda)=d^2,因此保證 D^2+\hat\lambda I 是可逆的。特徵方程 \phi(\lambda)=d^2 是一個非線性方程,使用數值方法可求出 \hat{\lambda} (例如牛頓法,見“牛頓法──非線性方程的求根方法”),代回即得到正則最小平方解 \hat{\mathbf{x}}

 
舉一個例子[4]

\displaystyle  \min_{\Vert\mathbf{x}\Vert=1}\left\Vert\begin{bmatrix}  2&0\\  0&1\\  0&0  \end{bmatrix}\mathbf{x}-\begin{bmatrix}  4\\  2\\  3  \end{bmatrix}\right\Vert^2

的特徵方程為

\displaystyle  \left(\frac{8}{4+\lambda}\right)^2+\left(\frac{2}{1+\lambda}\right)^2=1

這個問題的解是 \hat{\lambda}=4.57132\hat{\mathbf{x}}=(0.93334,0.35898)^T

 
註解
[1] 這裡使用較簡潔的寫法,約束優化問題的傳統表達式為

\displaystyle  \begin{array}{ll}  \hbox{minimize}&\Vert A\mathbf{x}-\mathbf{b}\Vert^2\\  \hbox{subject to}&C\mathbf{x}=\mathbf{d}.\end{array}

[2] 取自 G.H. Golub and C.F. Van Loan, Matrix Computations, 2nd edition, 1989, pp 568.
[3] 對於 \lambda>0,Tikhonov 正則化問題

\displaystyle  \min_{\mathbf{x}}\Vert A\mathbf{x}-\mathbf{b}\Vert^2+\lambda\Vert\mathbf{x}\Vert^2

的最佳解滿足 (A^TA+\lambda I)\mathbf{x}=A^T\mathbf{b} (見“每週問題 September 5, 2016”)。
[4] 同 [1], pp 565.

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

4 則回應給 約束最小平方問題

  1. Meiyue Shao 說道:

    那个n+p阶矩阵可逆的充要条件应当改成充分条件,必要性没有验证(事实上必要性不成立,比如n=p, A=0, C=I)。
    或者把rank(A)=n作为大前提,这样那个n+p阶矩阵可逆的充要条件就是rank(C)=p。

  2. Lin 說道:

    rank(A)=n这个条件太关键了。如果这个条件不成立,很多推导都不成立了。

發表迴響

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

WordPress.com Logo

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

Twitter picture

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

Facebook照片

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

Google+ photo

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

連結到 %s