【數學】主成分分析（PCA）的詳細深度推導過程

Based on Deep Learning (2017, MIT) book.

本文基於Deep Learning (2017, MIT)，推導過程補全了所涉及的知識及書中推導過程中跳躍和省略的部分。
blog

1 概述

現代資料集，如網路索引、高解析度影像、氣象學、實驗測量等，通常包含高維特徵，高緯度的資料可能不清晰、冗餘，甚至具有誤導性。資料視覺化和解釋變數之間的關係很困難，而使用這種高維資料訓練的神經網路模型往往容易出現過擬合（維度詛咒）。
主成分分析（PCA）是一種簡單而強大的無監督機器學習技術，用於資料降維。它旨在從大型變數集中提取一個較小的資料集，同時儘可能保留原始資訊和特徵（有失真壓縮）。PCA有助於識別資料集中最顯著和有意義的特徵，使資料易於視覺化。應用場景包括：統計學、去噪和為機器學習演算法預處理資料。

• 主成分是什麼？
  主成分是構建為原始變數的線性組合的新變數。這些新變數是不相關的，並且包含原始資料中大部分的資訊。

---

2 背景數學知識

這些知識對下一節的推導很重要。

• 正交向量和矩陣：
  • 如果兩個向量垂直，則它們是正交的。即兩個向量的點積為零。
  • 正交矩陣是一個方陣，其行和列是相互正交的單位向量；每兩行和兩列的點積為零，每一行和每一列的大小為1。
  • 如果\(A^T=A^{-1}\)或\(AA^T=A^TA=I\)，則\(A\)是正交矩陣。
  • 在機器人學中，旋轉矩陣通常是一個\(3\times3\)的正交矩陣，在空間變換中它會旋轉向量的方向但保持原始向量的大小。
• 矩陣、向量乘法規則：
  • \((AB)^T=B^TA^T\)，兩個矩陣的乘積的轉置。
  • \(\vec{a}^T\vec{b}=\vec{b}^T\vec{a}\)，兩個結果都是標量，標量的轉置是相同的。
  • \((A + B)C = AC + BC\)，乘法是可分配的。
  • \(AB \neq{} BA\)，乘法一般不滿足交換律。
  • \(A(BC)=(AB)C\)，乘法滿足結合律。
• 對稱矩陣：
  • \(A=A^T\)，\(A\)是對稱矩陣。
  • \(X^TX\)是對稱矩陣，因為\((X^TX)^T=X^TX\)。
• 向量導數規則（\(B\)是常量矩陣）：
  • \(d(x^TB)/dx=B\)
  • \(d(x^Tx)/dx=2x\)
  • \(d(x^TBx)/dx=2Bx\)
• 矩陣跡規則：
  • \(Tr(A)=Tr(A^T)\)
  • \(Tr(AB)=Tr(BA)\)
  • \(Tr(A)=\sum_i{\lambda_i}\)，其中\(\lambda\)是\(A\)的特徵值。
  • 跡在迴圈移位下不變：\(Tr(ABCD)=Tr(BCDA)=Tr(CDAB)=Tr(DABC)\)
• 向量和矩陣範數：
  • 向量的\(L^2\)範數，也稱為歐幾里得範數：\(||x||_2=\sqrt{\sum_i|x_i|^2}\)。
  • 通常使用平方的\(L^2\)範數來衡量向量的大小，可以計算為\(x^Tx\)。
  • Frobenius範數用於衡量矩陣的大小：\(||A||_F=\sqrt{\sum_{i,j}A^2_{i,j}}\)
  • Frobenius範數是所有矩陣元素的絕對平方和的平方根。
  • Frobenius範數是矩陣版本的歐幾里得範數。
• 特徵值分解和特徵值：
  • 方陣\(A\)的特徵向量是一個非零向量\(v\)，使得\(A\)的乘法僅改變\(v\)的比例：\(Av=\lambda v\)。\(\lambda\)是特徵值，\(v\)是特徵向量。
  • 假設矩陣\(A\)有\(n\)個線性無關的特徵向量\(v^{(i)}\)，我們可以將所有特徵向量連線起來形成一個矩陣\(V=[v^{(1)},\ldots,v^{(n)}]\)，並透過連線所有特徵值\(\lambda=[\lambda_1,\ldots,\lambda_n]^T\)形成一個向量，那麼\(A\)的特徵分解是\(A=Vdiag(\lambda)V^{-1}\)
  • 每個實對稱矩陣都可以分解為\(A=Q\Lambda Q^T\)，其中\(Q\)是由\(A\)的特徵向量組成的正交矩陣，\(\Lambda\)（讀作'lambda'）是一個對角矩陣。
• 拉格朗日乘數法：
  • 拉格朗日乘數法是一種在方程約束下尋找函式區域性最大值和最小值的策略。
  • 一般形式：\(\mathcal{L}(x,\lambda)=f(x)+\lambda\cdot g(x)\)，\(\lambda\)稱為拉格朗日乘子。

---

3 詳細PCA推導

需求描述

我們有\(m\)個點的輸入資料，表示為\({x^{(1)},...,x^{(m)}}\)在\(\mathbb{R}^{n}\)的實數集中。因此，每個點\(x^{(i)}\)是一個列向量，具有\(n\)維特徵。

需要對輸入資料進行有失真壓縮，將這些點編碼以表示它們的較低維度版本。換句話說，我們想要找到編碼向量\(c^{(i)}\in \mathbb{R}^{l}\)，\((l<n)\)來表示每個輸入點\(x^{(i)}\)。我們的目標是找到產生輸入的編碼向量的編碼函式\(f(x)=c\)，以及相應的重構（解碼）函式\(x\approx g(f(x))\)，根據編碼向量\(c\)計算原始輸入。

解碼的\(g(f(x))\)是一組新的點（變數），因此它與原始\(x\)是近似的。儲存\(c^{(i)}\)和解碼函式比儲存\(x^{(i)}\)更節省空間，因為\(c^{(i)}\)的維度較低。

解碼矩陣

我們選擇使用矩陣\(D\)作為解碼矩陣，將編碼向量\(c^{(i)}\)對映回\(\mathbb{R}^{n}\)，因此\(g(c)=Dc\)，其中\(D\in \mathbb{R}^{n\times l}\)。為了簡化編碼問題，PCA將\(D\)的列約束為彼此正交。

衡量重構的表現

在繼續之前，我們需要弄清楚如何生成最優的編碼點\(c^{*}\)，我們可以測量輸入點\(x\)與其重構\(g(c^*)\)之間的距離，使用\(L^2\)範數（或歐幾里得範數）：\(c^{*}=\arg\min_c||x-g(c)||_2\)。由於\(L^2\)範數是非負的，並且平方操作是單調遞增的，所以我們可以轉而使用平方的\(L^2\)範數：

\[c^{*}={\arg\min}_c||x-g(c)||_2^2 \]

向量的\(L^2\)範數是其分量的平方和，它等於向量與自身的點積，例如\(||x||_2=\sqrt{\sum|x_i|^2}=\sqrt{x^Tx}\)，因此平方的\(L^2\)範數可以寫成以下形式：

\[||x-g(c)||_2^2 = (x-g(c))^T(x-g(c)) \]

由分配率：

\[=(x^T-g(c)^T)(x-g(c))=x^Tx-x^Tg(c)-g(c)^Tx+g(c)^Tg(c) \]

由於\(x^Tg(c)\)和\(g(c)^Tx\)是標量，標量等於其轉置，\((g(c)^Tx)^T=x^Tg(c)\)，所以：

\[=x^Tx-2x^Tg(c)+g(c)^Tg(c) \]

為了找到使上述函式最小化的\(c\)，第一項可以省略，因為它不依賴於\(c\)，所以：

\[c^*={\arg\min}_c-2x^Tg(c)+g(c)^Tg(c) \]

然後用\(g(c)\)的定義\(Dc\)進行替換：

\[={\arg\min}_c-2x^TDc+c^TD^TDc \]

由於\(D\)的正交性和單位範數約束：

\[c^*={\arg\min}_c-2x^TDc+c^TI_lc \]

\[= {\arg\min}_c-2x^TDc+c^Tc \]

目標函式

現在目標函式是\(-2x^TDc+c^Tc\)，我們需要找到\(c^*\)來最小化目標函式。使用向量微積分，並令其導數等於0：

\[\nabla_c(-2x^TDc+c^Tc)=0 \]

根據向量導數規則：

\[-2D^Tx+2c=0 \Rightarrow c=D^Tx \]

找到編碼矩陣 \(D\)

所以編碼器函式是 \(f(x)=D^Tx\)。因此我們可以定義 PCA 重構操作為 \(r(x)=g(f(x))=D(D^Tx)=DD^Tx\)。

因此編碼矩陣 \(D\) 也被重構過程使用。我們需要找到最優的 \(D\) 來最小化重構誤差，即輸入和重構之間所有維度特徵的距離。這裡使用 Frobenius 範數（矩陣範數）定義目標函式：

\[D^*={\arg\min}_D\sqrt{\sum_{i,j}(x_j^{(i)}-r(x^{i})_j)^2},\quad D^TD=I_l \]

從考慮 \(l=1\) 的情況開始（這也是第一個主成分），\(D\) 是一個單一向量 \(d\)，並使用平方 \(L^2\) 範數形式：

\[d^*={\arg\min}_d{\sum_{i}||(x^{(i)}-r(x^{i}))}||_2^2, ||d||_2=1 \]

\[= {\arg\min}_d{\sum_{i}||(x^{(i)}-dd^Tx^{(i)})||_2^2}, ||d||_2=1 \]

\(d^Tx^{(i)}\) 是一個標量：

\[= {\arg\min}_d{\sum_{i}||(x^{(i)}-d^Tx^{(i)}d)}||_2^2, ||d||_2=1 \]

標量等於其自身的轉置：

\[d^*= {\arg\min}_d{\sum_{i}||(x^{(i)}-x^{(i)T}dd)}||_2^2, ||d||_2=1 \]

使用矩陣形式表示

令 \(X\in \mathbb{R}^{m\times n}\) 表示所有描述點的向量堆疊，即 \(\{x^{(1)^T}, x^{(2)^T}, \ldots, x^{(i)^T}, \ldots, x^{(m)^T}\}\)，使得 \(X_{i,:}=x^{(i)^T}\)。

\[X = \begin{bmatrix} x^{(1)^T}\\ x^{(2)^T}\\ \ldots\\ x^{(m)^T} \end{bmatrix} \Rightarrow Xd = \begin{bmatrix} x^{(1)^T}d\\ x^{(2)^T}d\\ \ldots\\ x^{(m)^T}d \end{bmatrix} \]

\[\Rightarrow Xdd^T = \begin{bmatrix} x^{(1)^T}dd^T\\ x^{(2)^T}dd^T\\ \ldots\\ x^{(m)^T}dd^T\\ \end{bmatrix} \]

\[\Rightarrow X-Xdd^T = \begin{bmatrix} x^{(1)^T}-x^{(1)^T}dd^T\\ x^{(2)^T}-x^{(2)^T}dd^T\\ \ldots\\ x^{(m)^T}-x^{(m)^T}dd^T\\ \end{bmatrix} \]

矩陣中的一行的轉置：

\[(x^{(i)^T}-x^{(i)^T}dd^T)^T=x^{(i)}-dd^Tx^{(i)} \]

由於 \(d^Tx^{(i)}\) 是標量：

\[=x^{(i)}-d^Tx^{(i)}d=x^{(i)}-x^{(i)^T}dd \]

所以我們知道 \(X\) 的第 \(i\) 行的 \(L^2\) 範數與原始形式相同，因此我們可以使用矩陣重寫問題，並省略求和符號：

\[d^*={\arg\min}_{d}||X-Xdd^T||_F^2, \quad d^Td=1 \]

利用矩陣跡規則簡化 Frobenius 範數部分如下：

\[{\arg\min}_{d}||X-Xdd^T||_F^2 \]

\[={\arg\min}_{d}Tr((X-Xdd^T)^T(X-Xdd^T)) \]

\[={\arg\min}_{d}-Tr(X^TXdd^T)-Tr(dd^TX^TX)+Tr(dd^TX^TXdd^T) \]

\[={\arg\min}_{d}-2Tr(X^TXdd^T)+Tr(X^TXdd^Tdd^T) \]

由於 \(d^Td=1\)：

\[={\arg\min}_{d}-2Tr(X^TXdd^T)+Tr(X^TXdd^T) \]

\[={\arg\min}_{d}-Tr(X^TXdd^T) \]

\[={\arg\max}_{d}Tr(X^TXdd^T) \]

由於跡是迴圈置換不變的，將方程重寫為：

\[d^*={\arg\max}_{d}Tr(d^TX^TXd), \quad d^Td=1 \]

由於 \(d^TX^TXd\) 是實數，因此跡符號可以省略：

\[d^*={\arg\max}_{d}d^TX^TXd,\quad d^Td=1 \]

尋找最優的 \(d\)

現在的問題是找到最優的 \(d\) 來最大化 \(d^TX^TXd\)，並且有約束條件 \(d^Td=1\)。

使用拉格朗日乘子法來將問題描述為關於 \(d\) 的形式：

\[\mathcal{L}(d,\lambda)=d^TX^TXd+\lambda(d^Td-1) \]

對 \(d\) 求導數（向量導數規則）：

\[\nabla_d\mathcal{L}(d,\lambda)=2X^TXd+2\lambda d \]

令導數等於0，\(d\) 將是最優的：

\[2X^TXd+2\lambda d=0 \]

\[X^TXd=-\lambda d \]

\[X^TXd=\lambda' d,\quad(\lambda'=-\lambda) \]

這個方程是典型的矩陣特徵值分解形式，\(d\) 是矩陣 \(X^TX\) 的特徵向量，\(\lambda'\) 是對應的特徵值。

利用上述結果，讓我們重新審視原方程：

\[d^*={\arg\max}_{d}d^TX^TXd, \quad d^Td=1 \]

\[={\arg\max}_{d}d^T\lambda' d \]

\[={\arg\max}_{d}\lambda'd^T d \]

\[={\arg\max}_{d}\lambda' \]

現在問題已經變的非常清楚了，\(X^TX\) 的最大特徵值會最大化原方程的結果，因此最優的 \(d\) 是矩陣 \(X^TX\) 對應最大特徵值的特徵向量。

這個推導是針對 \(l=1\) 的情況，只包含第一個主成分。當 \(l>1\) 時，\(D=[d_1, d_2, \ldots]\)，第一個主成分 \(d_1\) 是矩陣 \(X^TX\) 對應最大特徵值的特徵向量，第二個主成分 \(d_2\) 是對應第二大特徵值的特徵向量，以此類推。

---

4 總結

我們有一個資料集，包含 \(m\) 個點，記為 \({x^{(1)},...,x^{(m)}}\)。
令 \(X\in \mathbb{R}^{m\times n}\) 為將所有這些點堆疊而成的矩陣：\([x^{(1)^T}, x^{(2)^T}, \ldots, x^{(i)^T}, \ldots, x^{(m)^T}]\)。

主成分分析（PCA）編碼函式表示為 \(f(x)=D^Tx\)，重構函式表示為 \(x\approx g(c)=Dc\)，其中 \(D=[d_1, d_2, \ldots]\) 的列是 \(X^TX\) 的特徵向量，特徵向量對應的特徵值大小為降序排列。\(D^Tx\)即是降維度之後的資料。