廣義特徵值問題

本文的閱讀等級:高級

對於 n\times n 階矩陣 A,一般特徵值問題欲解出 A\mathbf{x}=\lambda\mathbf{x},其中 \lambdaA 的特徵值,\mathbf{x} 是對應的特徵向量。在一些工程和統計問題中,譬如,自由振動系統,譜聚類分析 (spectral clustering)[1],我們面對的是廣義特徵值 (generalized eigenvalue) 問題:A\mathbf{x}=\lambda B\mathbf{x},其中 AB 是兩個 n\times n 階 Hermitian 矩陣 (或實對稱矩陣),\lambda 稱為 AB 的廣義特徵值,\mathbf{x} 是對應的廣義特徵向量[2]。在多數的應用場合,B 是一正定矩陣。本文將推導自由振動系統的動態方程 (譜聚類分析較為複雜,他日另文介紹),證明優化廣義 Rayleigh 商 (quotient) 等價於廣義特徵值問題,並討論廣義特徵值與廣義特徵向量的性質與算法。

 

Spring mass vibrations

Spring-mass vibrations

見上圖,考慮質量為 m_1m_2 的兩個緊密物體,彼此間由一彈簧連接 (彈力常數為 k_2),且兩物體各自連接至一端固定於牆面的彈簧 (彈力常數同為 k_1)。令 x_1x_2 表示兩物體的平衡穩定位置,並令向右位移為正。在彈簧形變不大的情況下,虎克定律 (Hooke’s law) 給出彈力與彈簧伸縮長度 x 之間的線性關係 f(x)=-kx,其中 k 為彈力常數,由彈簧材質決定,上式負號表示彈簧所產生的彈力 f 與其伸長 (或壓縮) 的方向相反。圖中顯示左中右三個彈簧的伸長量分別為 x_1x_2-x_1-x_2。因為物體各自受左右兩個彈簧作用,根據牛頓第二運動定律,物體受彈力產生加速度 ,於是得到描述兩物體振動的諧振子方程式:

\displaystyle\begin{aligned}  m_1{x}''_1&=-k_1x_1+k_2(x_2-x_1)=-(k_1+k_2)x_1+k_2x_2\\  m_2{x}''_2&=-k_2(x_2-x_1)+k_1(-x_2)=k_2x_1-(k_1+k_2)x_2.\end{aligned}

將上面兩式合併為常微分矩陣方程

\displaystyle  \begin{bmatrix}  m_1&0\\  0&m_2  \end{bmatrix}\begin{bmatrix}  x''_1\\  x''_2  \end{bmatrix}=\begin{bmatrix}  -(k_1+k_2)&k_2\\  k_2&-(k_1+k_2)  \end{bmatrix}\begin{bmatrix}  x_1\\  x_2  \end{bmatrix}

或簡明表示成

M\mathbf{x}''=K\mathbf{x}

其中 M=\begin{bmatrix}  m_1&0\\  0&m_2  \end{bmatrix}K=\begin{bmatrix}  -(k_1+k_2)&k_2\\  k_2&-(k_1+k_2)  \end{bmatrix}\mathbf{x}=\begin{bmatrix}  x_1\\  x_2  \end{bmatrix}。我們猜測二階微分方程的解為 \mathbf{x}(t)=e^{i\omega t}\mathbf{c},其中 \mathbf{c} 是初始條件,i=\sqrt{-1}。求二次導數 \mathbf{x}''=-\omega^2e^{i\omega t}\mathbf{c}=-\omega^2\mathbf{x},代回微分方程,即得廣義特徵值問題:

K\mathbf{x}=\lambda M\mathbf{x}

其中 \lambda=-\omega^2。因為 m_1m_2 都是正數,可知 M 是正定矩陣。如果 m_1=m_2=1,廣義特徵值問題退化為一般特徵值問題 K\mathbf{x}=\lambda\mathbf{x},詳細討論請見“自由振動系統的特徵值與特徵向量”。

 
考慮下列最佳化問題 (見“每週問題 July 18, 2011”):

\displaystyle  \max\frac{\mathbf{x}^\ast A\mathbf{x}}{\mathbf{x}^\ast B\mathbf{x}}

分式 \frac{\mathbf{x}^\ast A\mathbf{x}}{\mathbf{x}^\ast B\mathbf{x}} 稱為廣義 Rayleigh 商,其中 AB 是 Hermitian 矩陣,且 B 是正定矩陣。若 B=I,廣義 Rayleigh 商退化為 Rayleigh 商 \frac{\mathbf{x}^\ast A\mathbf{x}}{\mathbf{x}^\ast\mathbf{x}}。Hermitian 矩陣 A 的特徵值是實數 (“特殊矩陣 (9):Hermitian 矩陣”)。當 \mathbf{x} 是對應 A 的最大特徵值,記為 \lambda_{\max},的特徵向量時,Rayleigh 商有極大值 \lambda_{\max} (見“Hermitian 矩陣特徵值的變化界定”),所以最大化 Rayleigh 商 \frac{\mathbf{x}^\ast A\mathbf{x}}{\mathbf{x}^\ast\mathbf{x}} 等價於求解 A\mathbf{x}=\lambda\mathbf{x}。下面我們證明最大化廣義 Rayleigh 商 \frac{\mathbf{x}^\ast A\mathbf{x}}{\mathbf{x}^\ast B\mathbf{x}} 等價於求解廣義特徵值問題 A\mathbf{x}=\lambda B\mathbf{x}

 
Hermitian 正定矩陣 (以下簡稱正定矩陣) 的特徵值都是正數 (“正定矩陣的性質與判別方法”)。令 \mu_1,\ldots,\mu_n>0 是正定矩陣 B 的特徵值。將 B 正交對角化為 B=UDU^\ast,其中 D=\hbox{diag}(\mu_1,\ldots,\mu_n)U 是么正矩陣 (unitary matrix) 滿足 U^\ast=U^{-1}。將 B 分解如下:

\displaystyle  B=UDU^\ast=UD^{1/2}U^\ast UD^{1/2}U^\ast=C^\ast C

上式中,我們令 D^{1/2}=\hbox{diag}(\sqrt{\mu_1},\ldots,\sqrt{\mu_n})C=UD^{1/2}U^\ast 也是正定矩陣。令 \mathbf{y}=C\mathbf{x}。廣義 Rayleigh 商可改寫成

\displaystyle\begin{aligned}  \frac{\mathbf{x}^\ast A\mathbf{x}}{\mathbf{x}^\ast B\mathbf{x}}&=\frac{\mathbf{x}^\ast A\mathbf{x}}{\mathbf{x}^\ast C^\ast C\mathbf{x}}=\frac{(C^{-1}\mathbf{y})^\ast A(C^{-1}\mathbf{y})}{\mathbf{y}^\ast\mathbf{y}}\\  &=\frac{\mathbf{y}^\ast (C^\ast)^{-1}AC^{-1}\mathbf{y}}{\mathbf{y}^\ast\mathbf{y}}=\frac{\mathbf{y}^\ast H\mathbf{y}}{\mathbf{y}^\ast\mathbf{y}},\end{aligned}

上面我們令 H=(C^\ast)^{-1}AC^{-1},就有 A=C^\ast HC。不難確認 H 是 Hermitian 矩陣。經過變數變換,最大化廣義 Rayleigh 商可化簡為最大化 Rayleigh 商:

\displaystyle  \max\frac{\mathbf{y}^\ast H\mathbf{y}}{\mathbf{y}^\ast\mathbf{y}}

考慮特徵方程 H\mathbf{y}=\lambda\mathbf{y},即有 HC\mathbf{x}=\lambda C\mathbf{x}。等號兩邊同時左乘 C^\ast,可得

\displaystyle  A\mathbf{x}=C^\ast HC\mathbf{x}=\lambda C^\ast C\mathbf{x}=\lambda B\mathbf{x}

因此證明 H 的特徵值等於 AB 的廣義特徵值,換句話說,最大化廣義 Rayleigh 商等價於廣義特徵值問題。

 
正定矩陣 B 是可逆矩陣,廣義特徵方程 A\mathbf{x}=\lambda B\mathbf{x} 同義於 B^{-1}A\mathbf{x}=\lambda\mathbf{x},所以 AB 的廣義特徵值即為 B^{-1}A 的特徵值。然而 B^{-1}A 不再是 Hermitian 矩陣,除非 AB^{-1}=B^{-1}A。為了保留 Hermitian 矩陣的美好性質,我們比較喜愛 H\mathbf{y}=\lambda\mathbf{y},或表示為

\displaystyle  P^\ast AP\mathbf{y}=\lambda\mathbf{y}

其中 P=C^{-1}=UD^{-1/2}U^\ast。令 \lambda_1,\ldots,\lambda_nH=P^\ast AP 的特徵值,\mathbf{y}_1,\ldots,\mathbf{y}_n 是對應的單範正交 (orthonormal) 特徵向量。據此,

\displaystyle  A\mathbf{x}_i=\lambda_i B\mathbf{x}_i,~~i=1,\ldots,n

其中 \mathbf{x}_i=P\mathbf{y}_ii=1,\ldots,n,是對應的線性獨立廣義特徵向量 (因為 \{\mathbf{y}_i\} 是線性獨立集且 P 是可逆矩陣)。令 \Lambda=\hbox{diag}(\lambda_1,\ldots,\lambda_n)S=\begin{bmatrix}  \mathbf{x}_1&\cdots&\mathbf{x}_n  \end{bmatrix}。上面 n 個廣義特徵方程可合併成矩陣表達式:

\displaystyle\begin{aligned}  AS&=A\begin{bmatrix}  \mathbf{x}_1&\cdots&\mathbf{x}_n  \end{bmatrix}=\begin{bmatrix}  A\mathbf{x}_1&\cdots&A\mathbf{x}_n  \end{bmatrix}\\  &=\begin{bmatrix}  \lambda_1 B\mathbf{x}_1&\cdots&\lambda_n B\mathbf{x}_n  \end{bmatrix}\\  &=\begin{bmatrix}  B\mathbf{x}_1&\cdots&B\mathbf{x}_n  \end{bmatrix}\begin{bmatrix}  \lambda_1&&\\  &\ddots&\\  &&\lambda_n  \end{bmatrix}\\  &=BS\Lambda.\end{aligned}

廣義特徵值問題解法如下:(1) 分解 B=UDU^\ast,設 P=UD^{-1/2}U^\ast,算出 H=P^\ast AP;(2) 解特徵方程 H\mathbf{y}=\lambda\mathbf{y} 得到 (廣義) 特徵值 \lambda_1,\ldots,\lambda_n 和特徵向量 \mathbf{y}_1,\ldots,\mathbf{y}_n;(3) 計算廣義特徵向量 \mathbf{x}_i=P\mathbf{y}_ii=1,\ldots,n

 
最後我們整理廣義特徵值與廣義特徵向量的性質。

  1. 廣義特徵值 \lambda_1,\ldots,\lambda_n 都是實數。因為 H 是 Hermitian 矩陣,故特徵值皆為實數。
  2. 矩陣 A 相合於 H=P^\ast AP。根據 Sylvester 慣性定律 (law of inertia),AH 有相同的慣性,也就是說,AH 有相同的正特徵值個數,負特徵值個數與零特徵值個數 (見“相合變換”)。
  3. 因為 \{\mathbf{y}_i\} 是一個單範正交集,

    \displaystyle  \mathbf{x}_i^\ast B\mathbf{x}_j=\mathbf{x}_i^\ast C^\ast C\mathbf{x}_j=\mathbf{y}_i^\ast\mathbf{y}_j=\delta_{ij}

    其中 \delta_{ij} 是 Kronecker 記號 (\delta_{ij}=1i=j\delta_{ij}=0i\neq j)。上式即為 S^\ast BS=I

  4. 套用廣義特徵方程和性質 (3),可得

    S^\ast AS=S^\ast BS\Lambda=\Lambda

    合併以上結果,A 相合於廣義特徵值矩陣 \LambdaB 相合於單位矩陣 I。通過廣義特徵向量矩陣 SAB 可同時對角化,只不過這不是一般常見的相似變換。

 
註解:
[1] 維基百科:Spectral Clustering
[2] 本文所稱的廣義特徵向量與 Jordan 典型形式的廣義特徵向量不同,詳見“Jordan 形式大解讀之尋找廣義特徵向量“。

Advertisements
This entry was posted in 特徵分析, 線性代數專欄 and tagged , , , . Bookmark the permalink.

8 則回應給 廣義特徵值問題

  1. 張盛東 說:

    老師,不知可否順便介紹一下偽特征值,和偽譜(pseudospectrum)呢?

  2. Siang 說:

    老師,可以請教你一個小問題嗎~ 就是:兩物體振動的諧振子方程式第一式
    等號右邊 -k1x1+k2(x2-x1) 而不是 -k1x1-k2(x2-x1)呢?

    • Siang 說:

      因為這個問題應該是屬於物理的範疇~"~ 不知道可不可以問您><"

    • ccjou 說:

      假設k1=k2=k。當x2=0,不論x1是多少,你寫出的公式為-kx1-k(x2-x1)=0,這並不正確。

      • Siang 說:

        真的耶!?
        那老師針對剛剛諧振子方程式
        第一式等號右邊最後整理出
        -(k1+k2)x1+k2x2
        對於中間彈簧往右x2位移,彈力應該是向左吧?
        那對於中間與左邊彈簧擁有x1的位移,彈力也是向左
        為什麼整理出來的式子,是異號呢?(我一直卡在正負號的問題~"~)

  3. Meiyue Shao 說:

    实际计算的时候第一步通常对 B 做 Cholesky 分解, 这比计算平方根要便宜不少.

發表迴響

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

WordPress.com Logo

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

Twitter picture

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

Facebook照片

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

Google+ photo

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

連結到 %s