正交 Procrustes 問題

本文的閱讀等級:高級

普洛克路斯忒斯 (Procrustes) 是希臘神話中海神波塞頓 (Poseidon) 的兒子。他在雅典到埃萊夫西納 (Eleusis) 的神聖之路 (The Sacred Way) 上開設一間黑店,向路過的旅人謊稱店內設有一張適合所有人的鐵床。旅客投宿時,普洛克路斯忒斯將身高者截斷雙足,身矮者則強行拉長,使之與床的長短相同。從來沒有一個人的身長與鐵床的長度相同而免於凌遲,因為他暗地裡準備了兩張床[1]。後人於是以 Procrustean 表示「削足適履,殺頭便冠」,意思是將不同的長度、大小或屬性安裝到一個任意的標準。

正交 Procrustes 問題:給定兩個 m\times n 階實矩陣 AB,求一個 n\times n 階實正交矩陣 QQ^T=Q^{-1},使得 \Vert B-AQ\Vert_F 具有最小值[2],其中 \Vert\cdot\Vert_F 是 Frobenius 範數。

 
正交 Procrustes 問題是一個最小平方矩陣近似問題,可以這麼解讀:A 是旅人,B 是鐵床,正交 Procrustean 變換 (包含旋轉和鏡射) Q 即為施予旅人的肢體酷刑。我們的問題是求出一酷刑使旅人變形後的身長與鐵床的長度最為吻合。以下討論限定於實矩陣。

Procrustes 和他的客人(點圖看動畫) From http://www.mythweb.com/today/media/procrustes07.gif

 
對於一個 m\times n 階矩陣 A=[a_{ij}],Frobenius 範數定義為 (見“矩陣範數”)

\displaystyle \Vert A\Vert_{F}=\sqrt{\mathrm{trace}(A^TA)}=\sqrt{\sum_{i=1}^m\sum_{j=1}^n a_{ij}^2}

使用跡數循環不變性和轉置不變性 (見“跡數的性質與應用”),

\displaystyle\begin{aligned} \Vert B-AQ\Vert_F^2&=\hbox{trace}\left((B-AQ)^T(B-AQ)\right)\\ &=\hbox{trace}(B^TB-B^TAQ-Q^TA^TB+Q^TA^TAQ)\\ &=\hbox{trace}(B^TB)-\hbox{trace}(B^TAQ)-\hbox{trace}(Q^TA^TB)+\hbox{trace}(Q^TA^TAQ)\\ &=\hbox{trace}(B^TB)-\hbox{trace}(QB^TA)-\hbox{trace}(A^TBQ^T)+\hbox{trace}(A^TAQQ^T)\\ &=\hbox{trace}(B^TB)-\hbox{trace}(QB^TA)^T-\hbox{trace}(A^TBQ^T)+\hbox{trace}(A^TA)\\ &=\Vert B\Vert_F^2-2\hbox{trace}(A^TBQ^T)+\Vert A\Vert_F^2. \end{aligned}

因為 AB 是給定的矩陣,最小化 \Vert B-AQ\Vert_F 等價於最大化 \hbox{trace}(A^TBQ^T)。另外,設 B=0,上式表明 \Vert AQ\Vert_F=\Vert A\Vert_F,正交矩陣乘法不改變 Frobenius 範數。

 
正交 procrustes 問題可轉換為底下的約束最佳化問題:

\displaystyle \max_{QQ^T=I}\hbox{trace}(A^TBQ^T)

使用 Lagrange 乘數法 (見“Lagrange 乘數法”)。令 Lagrangian 函數為

\displaystyle L=\hbox{trace}(A^TBQ^T)-\hbox{trace}(\Lambda(QQ^T-I))

其中 \Lambdan\times n 階 Lagrange 乘數矩陣。利用跡數的導數公式 (見“跡數與行列式的導數”,tr-6,tr-13),計算 L 對於 Q 的導數:

\displaystyle\begin{aligned} \frac{\partial L}{\partial Q}&=\frac{\partial\hbox{trace}(A^TBQ^T)}{\partial Q}-\frac{\partial\hbox{trace}(Q^T\Lambda Q)}{\partial Q}+\frac{\partial\hbox{trace}(\Lambda)}{\partial Q}\\ &=A^TB-(\Lambda+\Lambda^T)Q. \end{aligned}

設上式為零,可得 A^TB=(\Lambda+\Lambda^T)Q,或 A^TBQ^T=\Lambda+\Lambda^T。設 n\times n 階矩陣 A^TB 的奇異值分解為 A^TB=U\Sigma V^T,其中 UVn\times n 階正交矩陣,\Sigma=\hbox{diag}(\sigma_1(A^TB),\ldots,\sigma_n(A^TB)),這裡 \sigma_i(A^TB)\ge 0A^TB 的奇異值 (見“奇異值分解 (SVD)”)。因為 \Lambda+\Lambda^T 是對稱矩陣,可知 A^TBQ^T=Q(A^TB)^T。代入 A^TB 的奇異值分解,

\displaystyle U\Sigma V^TQ^T=QV\Sigma^TU^T=QV\Sigma U^T

比較等號兩邊可得 QV=U,故 Q=UV^T 為一個最佳解,而且

\displaystyle\begin{aligned} \hbox{trace}(A^TBQ^T)&=\hbox{trace}(U\Sigma V^TVU^T)=\hbox{trace}(U\Sigma U^T)\\ &=\hbox{trace}(\Sigma UU^T)=\hbox{trace}\Sigma=\sum_{i=1}^n\sigma_i(A^TB).\end{aligned}

矩陣 B 的最小平方近似矩陣是 AQ=AUV^T,最小誤差為

\displaystyle \Vert B-AQ\Vert_F=\Vert B-AUV^T\Vert_F=\sqrt{\Vert A\Vert_F^2+\Vert B\Vert_F^2-2\sum_{i=1}^n\sigma_i(A^TB)}

上式提供了 B=AQ 的一個充要條件。若 B=AQ,則 \Vert B\Vert_F=\Vert AQ\Vert_F=\Vert A\Vert_F,就有 \Vert A\Vert_F^2=\Vert B\Vert_F^2=\sum_{i=1}^n\sigma_i(A^TB)。反向的推論十分明顯。

 
補充一個較快捷的推導方式。如果不採行 Lagrange 乘數法,直接將 A^TB=U\Sigma V^T 代入目標函數,如下:

\hbox{trace}(A^TBQ^T)=\hbox{trace}(U\Sigma V^TQ^T)=\hbox{trace}(\Sigma V^TQ^TU)=\hbox{trace}(\Sigma P)

上面令 P=[p_{ij}]=V^TQ^TU。因為 U, V, Q 都是正交矩陣,不難證明 P^T=P^{-1}P 也是一個正交矩陣。明確寫出

\displaystyle \hbox{trace}(\Sigma P)=\sum_{i=1}^n\sigma_i(A^TB)p_{ii}

\sum_{i=1}^np^2_{ij}=1\sum_{j=1}^np^2_{ij}=1,上式的最大值必發生於 p_{ii}=1i=1,\ldots,n,即 P=I。因此,V^TQ^TU=I,解得 Q=UV^T

 
舉一個例子[3],若

\displaystyle A=\begin{bmatrix} 1.2&2.1\\ 2,9&4.3\\ 5.2&6.1\\ 6.8&8.1 \end{bmatrix},~~ B=\begin{bmatrix} 1&2\\ 3&4\\ 5&6\\ 7&8 \end{bmatrix}

則正交矩陣

\displaystyle Q=\left[\!\!\begin{array}{cr} 0.9999 &-0.0126\\ 0.0126&0.9999 \end{array}\!\!\right]

最小化 \Vert B-AQ\Vert_F,最小誤差為 0.4661

 
接下來我們探討正交 Procrustes 問題的一個變形問題,稱為兩面正交 Procrustes 問題。

兩面正交 Procrustes 問題:給定兩個 m\times n 階實矩陣 AB,求一個 n\times n 階實正交矩陣 Q 和一個 m\times m 階實正交矩陣 P,使得 \Vert B-PAQ\Vert_F 具有最小值。

 
仿造前面的作法,計算目標函數

\displaystyle\begin{aligned} \Vert B-PAQ\Vert_F^2&=\hbox{trace}\left((B-PAQ)^T(B-PAQ)\right)\\ &=\hbox{trace}(B^TB-B^TPAQ-Q^TA^TP^TB+Q^TA^TP^TPAQ)\\ &=\hbox{trace}(B^TB)-\hbox{trace}(B^TPAQ)-\hbox{trace}(Q^TA^TP^TB)+\hbox{trace}(Q^TA^TAQ)\\ &=\hbox{trace}(B^TB)-\hbox{trace}(QB^TPA)-\hbox{trace}(A^TP^TBQ^T)+\hbox{trace}(A^TAQQ^T)\\ &=\hbox{trace}(B^TB)-\hbox{trace}(QB^TPA)^T-\hbox{trace}(A^TP^TBQ^T)+\hbox{trace}(A^TA)\\ &=\Vert B\Vert_F^2-2\hbox{trace}(A^TP^TBQ^T)+\Vert A\Vert_F^2, \end{aligned}

所以最小化 \Vert B-PAQ\Vert_F 等價於最大化 \hbox{trace}(A^TP^TBQ^T)

 
AB 的奇異值分解為 A=U_1\Sigma_1V_1^TB=U_2\Sigma_2V_2^T,其中 U_1U_2m\times m 階正交矩陣,V_1V_2n\times n 階正交矩陣,\Sigma_1=\begin{bmatrix} S_1&0\\ 0&0 \end{bmatrix}\Sigma_2=\begin{bmatrix} S_2&0\\ 0&0 \end{bmatrix}m\times n 階奇異值矩陣 (零列和零行必須隨 mn 的大小調整),S_1=\hbox{diag}(\sigma_1(A),\ldots,\sigma_s(A))S_2=\hbox{diag}(\sigma_1(B),\ldots,\sigma_s(B))s=\min\{m,n\}。代入奇異值分解,使用跡數循環不變性,

\displaystyle\begin{aligned} \hbox{trace}(A^TP^TBQ^T)&=\hbox{trace}(V_1\Sigma_1^TU_1^TP^TU_2\Sigma_2V_2^TQ^T)\\ &=\hbox{trace}(\Sigma_1^TU_1^TP^TU_2\Sigma_2V_2^TQ^TV_1)\\ &=\hbox{trace}(\Sigma_1^TX^T\Sigma_2Y^T), \end{aligned}

上面令 m\times m 階正交矩陣 X^T=U_1^TP^TU_2n\times n 階正交矩陣 Y^T=V_2^TQ^TV_1。兩式取轉置,解得 P=U_2XU_1^TQ=V_1YV_2^T。剩下的工作是找出正交矩陣 XY 使得 \hbox{trace}(\Sigma_1^TX^T\Sigma_2Y^T) 有最大值。

 
考慮約束最佳化問題

\displaystyle \max_{XX^T=I\atop{YY^T=I}}\hbox{trace}(\Sigma_1^TX^T\Sigma_2Y^T)

使用 Lagrange 乘數法。寫出 Lagrangian 函數

\displaystyle L=\hbox{trace}(\Sigma_1^TX^T\Sigma_2Y^T)-\hbox{trace}(\Lambda(XX^T-I))-\hbox{trace}(M(YY^T-I))

其中 \Lambdam\times m 階且 Mn\times n 階 Lagrange 乘數矩陣。利用跡數的導數公式,計算 L 對於 XY 的導數:

\displaystyle\begin{aligned} \frac{\partial L}{\partial X}&=\frac{\partial\hbox{trace}(X^T\Sigma_2Y^T\Sigma_1^T)}{\partial X}-\frac{\partial\hbox{trace}(X^T\Lambda X-\Lambda)}{\partial X}-\frac{\partial\hbox{trace}(Y^TMY-M)}{\partial X}\\ &=\Sigma_2Y^T\Sigma_1^T-(\Lambda+\Lambda^T)X, \end{aligned}

\displaystyle\begin{aligned} \frac{\partial L}{\partial Y}&=\frac{\partial\hbox{trace}(\Sigma_1^TX^T\Sigma_2Y^T)}{\partial Y}-\frac{\partial\hbox{trace}(X^T\Lambda X-\Lambda)}{\partial Y}-\frac{\partial\hbox{trace}(Y^TMY-M)}{\partial Y}\\ &=\Sigma_1^TX^T\Sigma_2-(M+M^T)Y. \end{aligned}

設上面兩式為零,可得最佳解的必要條件:

\displaystyle\begin{aligned} \Sigma_2Y^T\Sigma_1^TX^T&=\Lambda+\Lambda^T\\ \Sigma_1^TX^T\Sigma_2Y^T&=M+M^T. \end{aligned}

因為 \Lambda+\Lambda^TM+M^T 是對稱矩陣,

\displaystyle\begin{aligned} \Sigma_2Y^T\Sigma_1^TX^T&=X\Sigma_1Y\Sigma_2^T\\ \Sigma_1^TX^T\Sigma_2Y^T&=Y\Sigma_2^TX\Sigma_1 .\end{aligned}

第二式右乘 \Sigma_1^TX^T

\displaystyle \Sigma_1^TX^T\Sigma_2Y^T\Sigma_1^TX^T=Y\Sigma_2^TX\Sigma_1\Sigma_1^TX^T

將第一式代入等號左邊,

\displaystyle \Sigma_1^TX^T\Sigma_2Y^T\Sigma_1^TX^T=\Sigma_1^TX^TX\Sigma_1Y\Sigma_2^T=\Sigma_1^T\Sigma_1Y\Sigma_2^T

合併上面兩式,右乘 X,就有

\displaystyle (Y\Sigma_2^TX)(\Sigma_1\Sigma_1^T)=(\Sigma_1^T\Sigma_1)(Y\Sigma_2^TX)

其中 \Sigma_1\Sigma_1^Tm\times m 階對角矩陣,\Sigma_1^T\Sigma_1n\times n 階對角矩陣,兩者有相同的主對角元 (扣除零元) \sigma_i^2(A)1\le i\le s。比較等號兩邊的 (i,j) 元,可得

\displaystyle (Y\Sigma_2^TX)_{ij}\sigma_j^2(A)=\sigma_i^2(A)(Y\Sigma_2^TX)_{ij}

對於 i\neq j,若 \sigma_i(A)\neq\sigma_j(A),則 (Y\Sigma_2^TX)_{ij}=0。這個結果提示 X=I_mY=I_n 為一組最佳解 (事實上,所有的非零純量矩陣 X=cI_mY=(1/c)I_n 皆為最佳解),也就有 P=U_2U_1^TQ=V_1V_2^T,目標函數可化簡成

\displaystyle \hbox{trace}(A^TP^TBQ^T)=\hbox{trace}(\Sigma_1^T\Sigma_2)=\sum_{i=1}^s\sigma_i(A)\sigma_i(B)

根據排序不等式 (rearrangement inequality)[4],上式的最大值發生於 \sigma_1(A)\ge\cdots\ge\sigma_s(A)\sigma_1(B)\ge\cdots\ge\sigma_s(B)。因此,矩陣 B 的最小平方近似矩陣為

\displaystyle PAQ=U_2U_1^TU_1\Sigma_1V_1^TV_1V_2^T=U_2\Sigma_1V_2^T

最小誤差則為

\displaystyle\begin{aligned} \Vert B-PAQ\Vert_F&=\Vert B-U_2\Sigma_1V_2^T\Vert_F=\sqrt{\Vert A\Vert_F^2+\Vert B\Vert_F^2-2\sum_{i=1}^s\sigma_i(A)\sigma_i(B)}\\ &=\sqrt{\sum_{i=1}^s\sigma_i^2(A)+\sum_{i=1}^s\sigma_i^2(B)-2\sum_{i=1}^s\sigma_i(A)\sigma_i(B)}\\ &=\sqrt{\sum_{i=1}^s\left(\sigma_i(A)-\sigma_i(B)\right)^2}, \end{aligned}

上面我們使用了 Frobenius 範數的奇異值恆等式 (見“SVD 於矩陣近似的應用”)

\displaystyle \Vert A\Vert_F=\sqrt{\sum_{i=1}^s\sigma_i^2(A)}

在兩面正交 Procrustes 問題中,從最小誤差表達式可知 B=PAQ 的充分必要條件是 AB 有相同的奇異值。

 
直觀上,兩面正交 Procrustes 問題 \min_{PP^T=I,QQ^T=I}\Vert B-PAQ\Vert_F 的極小值應比正交 Procrustes 問題 \min_{QQ^T=I}\Vert B-AQ\Vert_F 的極小值來得小,這意味

\displaystyle \sum_{i=1}^n\sigma_i(A^TB)\le\sum_{i=1}^{\min\{m,n\}}\sigma_i(A)\sigma_i(B)

這個不等式的證明過程相當冗長,在此從略。

 
參考來源:
[1] 維基百科:Procrustes
[2] 正交 Procrustes 問題也可以表示為求一個 m\times m 階正交矩陣 Q 使最小化 \Vert B-QA\Vert_F。因為 \Vert A\Vert_F=\hbox{trace}(A^TA)=\hbox{trace}(AA^T)=\Vert A^T\Vert_F,推得 \Vert B-QA\Vert_F=\Vert B^T-A^TQ^T\Vert,其中 Q^T 是一個正交矩陣。
[3] 取自 G. H. Golub and C. F. Van Loan, Matrix Computations, 2nd edition, 1989, pp 582.
[4] 維基百科:排序不等式

相關閱讀:
Advertisements
本篇發表於 線性代數專欄, 二次型 並標籤為 , , , , , , 。將永久鏈結加入書籤。

8 則回應給 正交 Procrustes 問題

  1. 張盛東 說道:

    老師,這種矩陣近似一般在什麼場合下使用呢?

  2. 徐君 說道:

    老师,你在文中说:比較等號兩邊可得 QV=U,故 Q=UV^T 為一最佳解。那这个解是不是唯一的吗?或者说正交 Procrustes 問題是不是一个凸问题?

  3. 徐君 說道:

    老师,你在文中说: 設上式為零,可得 A^TB=(\Lambda+\Lambda^T)Q,或 A^TBQ^T=\Lambda+\Lambda^T。这里用到了QQ^T=I这个条件。可是,在Lagrangian 函數里还能用到之前有的约束条件(比如QQ^T=I)吗?

  4. 徐君 說道:

    正交 Procrustes 問題里的Q可以有更直观的解决方法,请看https://jianfengwang.files.wordpress.com/2015/07/proof_of_the_orthogonal_procrustes_problem.pdf

    • ccjou 說道:

      謝謝,確實是有更簡潔的推導法,有空我再補上來。

      (見“補充一個較快捷的推導方式”該段)

  5. 徐君 說道:

    谢谢老师!老师这篇帖子里提到的技巧帮助我证明了我论文中的一个定理。

發表迴響

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

WordPress.com Logo

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

Twitter picture

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

Facebook照片

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

Google+ photo

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

連結到 %s