本文主要思路如下：

1 PCA優化目標

PCA（主成分分析）是一種數據降維的方法，即用較少特徵地數據表達較多特徵地數據（數據壓縮，PCA屬於有損壓縮）。PCA推導有兩種主要思路：

  1. 最大化數據投影后的的方差（讓數據更分散）
  2. 最小化投影造成的損失

本文採用第一種思路完成推導過程，下圖中旋轉的是新坐標軸，每個數據點在改坐標軸上垂直投影，最佳的坐標軸為數據投影后的數據之間距離最大。

圖1 數據投影到新坐標軸

要完成PCA推導過程，需要如下第 2 章部分的理論依據

2 理論依據

2.1 矩陣換基底

坐標變換地目標是，找到一組新的正交單位向量，替換原來的正交單位向量。下面通過具體例子說明。

假設存在向量 vec{a}=egin{bmatrix} 3\2end{bmatrix} ，要變換導以 vec{u}=egin{bmatrix}{1oversqrt{2}}\ {1oversqrt{2}}end{bmatrix}, vec{v}=egin{bmatrix}-{1oversqrt{2}}\ {1oversqrt{2}}end{bmatrix} 為新基底地坐標上，求在心坐標系中的坐標

ecause 向量 vec{a} 在向量 vec{u} 上的投影距離 s：

s=||vec{a}||cdot cos heta = {vec{a}cdotvec{u}over ||vec{u}||}=vec{a}cdotvec{u}

其中： heta 表示兩個向量之間的夾角

herefore a_u=vec{u}^Tcdot vec{a}, a_v=vec{v}^Tcdot vec{a}

herefore 向量 vec{a} 在新坐標系中的坐標可以表示為：

vec{a}_{new} =egin{bmatrix}vec{u}&vec{v}end{bmatrix}^Tcdot vec{a} =egin{bmatrix} vec{u}^Tcdotvec{a}\ vec{v}^Tcdotvec{a}\ end{bmatrix}

如果矩陣 A 的列向量分別表示原來坐標系中的點，那麼在新坐標系中的坐標為：

A_{new} =egin{bmatrix}vec{u}&vec{v}end{bmatrix}^Tcdot A

如果 vec{a}_{center} 表示一系列數據點的中心，那麼可以證明：

vec{a}_{newcenter}=egin{bmatrix}vec{u}&vec{v}end{bmatrix}^Tcdot vec{a}_{center}

經過上面的變換之後，新坐標系相比原坐標系順時針旋轉了45度； vec{a} 相對新坐標系位置和相對原坐標系位置發生了逆時針旋轉45度。即：上述變換過程為向量的旋轉過程，旋轉的角度=-坐標系旋轉角度

如果 ||vec{u}|| eq1, ||vec{v}|| eq 1 ，那麼：

vec{u}^Tcdot vec{a}=scdot ||vec{u}||\ vec{v}^Tcdot vec{a}=scdot ||vec{v}||

即： vec{a}_{new} 相比 vec{a} ，2個坐標分別放大了 ||vec{u}|| 倍和 ||vec{v}|| 倍。即向量發生了伸縮。

2.2 拉格朗日乘子法

拉格朗日乘子法主要提供了一種求解函數在約束條件下極值的方法。下面還是通過一個例子說明。

假設存在一個函數 f(x,y) ，求該函數在 g(x,y)=c 下的極值（可以是極大，也可以極小）

通過觀察我們發現，在極值點的時候兩個函數必然相切，即此時各自的導數成正比，從而：

{partial foverpartial x}=lambda{partial goverpartial x}\ {partial foverpartial y}=lambda{partial goverpartial y}\ g(x,y)=c

通過聯立上述三個公式，既可以求出最終結果。拉格朗日運算元的主要思路同上，不過他假設了一個新的函數：

F(x,y,lambda)=f(x,y)+lambda[c-g(x,y)]

然後分解求：

{partial Foverpartial x}=0\ {partial Foverpartial y}=0\ {partial Foverpartial lambda}=0

從而完成求解過程

2.3 協方差矩陣

假設有一組數據：

egin{array}{c|cccc} 樣本編號&變數x（如發傳單數量）&變數y（如購買數量）&變數z（如購買總價）\ hline 1& 1 &2&3\ 2&35&25&55\ cdots&cdots&cdots&cdots\ end{array}

協方差研究的目的是變數（特徵）之間的關係，也就是上表中的發傳單數量、購買數量、購買總額之間的相關情況

上表數據用矩陣表示為：

X=egin{bmatrix} 1&35&cdots\ 2&25&cdots\ 3&55&cdots\ end{bmatrix}

那麼兩兩變數之間的關係：

cov(x,y)=E[(1-E(x))(2-E(y))+(35-E(x))(25-E(y))+cdots]\ cov(x,z)=E[(1-E(x))(3-E(z))+(35-E(x))(55-E(z))+cdots]\ cdots

如果E(x)=E(y)=E(z)=0（可以通過數據初始化實現），那麼上述的協方差關係可以用如下矩陣乘法表示：

cov(X)={1over m}XX^T=egin{bmatrix} cov(x,x)&cov(x,y)&cov(x,z)\ cov(y,x)&cov(y,y)&cov(y,z)\ cov(z,x)&cov(z,y)&cov(z,z)\ end{bmatrix}

如果把對角線上的數據加起來會發現：

cov(x,x)+cov(y,y)+cov(z,z)=E[([(1-E(x)]^2+[2-E(y)]^2+[3-E(z)]^2)+cdots]={1over m}sum_{i=1}^m||X_i-X_{center}||^2

也就是說每個樣本點到樣本中心距離的平方和的平均 = 樣本各個特徵方差和（自身協方差）= ? sum diag({1over m}XX^T) ，即樣本的方差

2.4 特徵向量和奇異值分解

2.4.1 特徵向量

參考：

特徵值和特徵向量?

www.jianshu.com圖標

假設：左側矩形由 egin{bmatrix} vec{i}&vec{j} end{bmatrix}=egin{bmatrix} 1&0\ 0&1\ end{bmatrix} 定義，右側矩形由 egin{bmatrix} vec{i}&vec{j} end{bmatrix}=egin{bmatrix} 2&0\ 0&0.5\ end{bmatrix} 定義。

根據 2.1 矩陣拉伸變換的結果，變換矩陣 A=egin{bmatrix} vec{u}^T\ vec{v}^T\ end{bmatrix}=egin{bmatrix} 2&0\ 0&0.5\ end{bmatrix} ，即：

Acdotegin{bmatrix} vec{i}&vec{j} end{bmatrix}=egin{bmatrix} vec{i}&vec{j}\ end{bmatrix}

在應用變換矩陣變換時，我們發現存在與上圖中紅色向量平行的向量 vec{a} ，他們總滿足：

Acdot vec{a} // vec{a}

即：

Acdot vec{a}=lambdacdot vec{a}

所以：紅色的特徵向量不受變換矩陣的影響，仍保持原來的方向，我們稱這類向量為變換矩陣A的特徵向量，對應的 lambda 為特徵值。又因為特徵向量有很多個，即：

Acdot vec{a_i}=lambda_icdot vec{a_i}

所以：

Acdotegin{bmatrix} vec{a_1}&vec{a_2}&cdots end{bmatrix}=egin{bmatrix} vec{a_1}&vec{a_2}&cdots end{bmatrix}cdotegin{bmatrix} lambda_1&&\ &lambda_2&\ &&ddots end{bmatrix}Rightarrow A=Qcdot Sigmacdot Q^{-1}

其中：Q的列向量都是A變換矩陣的特徵向量

另外，在做旋轉變換時，要求變換前後的坐標維度不發生改變，即A須為方陣

綜上：如果方陣A滿足 A=Qcdot Sigmacdot Q^{-1} ，那麼Q為特徵向量， Sigma 為對應的特徵值

2.4.2 奇異值分解

奇異值分解（svd: singular value decomposition）定義：對於任意的矩陣A，存在：

A_{m imes n}=U_{m imes m}cdotSigma_{m imes n}cdot V_{n imes n}^T

其中：

U^Tcdot U=I_m\ V^Tcdot V=I_n

即：U的列向量兩兩正交且模為1，V列向量兩兩正交且模為1，即：

U^T=U^{-1}\ V^T=V^{-1}

2.4.3 特徵向量和奇異值分解的關係

對於任意矩陣 A ，對A做svd有：

AA^T=USigma V^Tcdot VSigma U^T=USigma^2U^{-1}

令 Sigma=Sigma^2 ，則：

AA^T=USigma U^{-1}\ 滿足A=QSigma Q^{-1}特徵向量定義

所以 AA^T 能實現特徵分解，又因為：

AA^T=underbrace{USigma V^T}_{svd}

所以：

U=U\ Sigma=Sigma\ U^{-1}=V^TRightarrow U=V

因此：對 AA^T 做SVD，那麼得到的U列向量為特徵向量（對應A的U矩陣）， Sigma 為特徵值對角陣

同理：對 A^TA做SVD，那麼得到的U列向量為特徵向量（對應A的V矩陣）， Sigma 為特徵值對角矩陣

3 PCA

3.1 PCA推導

PCA的目標是找到一組新的正交基 {u_1,u_2,cdots,u_k} （從n維下降到k維），使得數據點在該正交基構成的平面上投影后，數據間的距離最大，即數據間的方差最大。如果數據在每個正交基上投影后的方差最大，那麼同樣滿足在正交基所構成的平面上投影距離最大。

根據2.1，設正交基 u_j ，數據點 x_i 在該基底上的投影距離為 x_i^Tcdot u_j ，所以所有數據在該基底上投影的方差 J_j 為：

J_j={1over m}sum_{i=1}^m(x_i^Tu_j-x_{center}^Tu_j)^2

其中：m為樣本數量，在數據運算之前對數據 x 進行0均值初始化，即 x_{center}=0 ，從而：

J_j={1over m}sum_{i=1}^m(x_i^Tu_j)^2 ={1over m}sum_{i=1}^m(u_j^Tx_icdot x_i^T u_j) =u_j^Tcdot {1over m}sum_{i=1}^m(x_ix_i^T)cdot u_j

所以：

J_j =u_j^Tcdot {1over m}(x_1x_1^T+x_2x_2^T+cdots+x_mx_m^T)cdot u_j =u_j^Tcdot {1over m}(egin{bmatrix}x_1&cdots&x_mend{bmatrix}cdotegin{bmatrix}x_1\ vdots \x_mend{bmatrix})cdot u_j =={1over m}u_j^TXX^Tu_j

由於 {1over m}XX^T 為常數，這裡假設 S={1over m}XX^T ，則： J_j=u_j^Tcdot Scdot u_j ，根據PCA目標，我們需要求解 J_j 最大時對應的 u_j

根據 2.2 中的拉格朗日運算元（求極值）求解：

J_{j}=u_j^TSu_j\ s.t. u_j^Tu_j=1

則構造函數：

F(u_j)=u_j^TSu_j+lambda_j(1-u_j^Tu_j)

求解 {partial Foverpartial u_j}=0 ，得：

2Scdot u_j-2lambda_jcdot u_j=0 Rightarrow Su_j=lambda_ju_j

結合2.4.1則：當 u_j、lambda_j 分別為S矩陣的特徵向量、特徵值時， J_j 有極值，把上述結果帶回公式得：

{J_j}_m=u_j^Tlambda_j u_j=lambda_j

所以對於任意滿足條件的正交基，對應的數據在上面投影后的方差值為S矩陣的特徵向量，從而：

J_{max}=sum_{j=1}^k lambda_j，lambda從大到小排序

所以投影正交基為S的特徵向量中的前k個最大特徵值對應的特徵向量。

接下來對S進行特徵分解，根據2.4.3特徵向量和svd的關係結論，S的特徵向量集合：

U = Uof svd(S)=U of svd({1over m}XX^T)

另外，由於 S={1over m}XX^T 由於X已0均值處理，根據2.3 協方差矩陣定義：S為數據集X的協方差矩陣。

綜上，即可得到滿足投影后數據距離最大的新的正交基 {u_1,u_2,cdots,u_k}

因此：

X_{new_{k imes m}}=egin{bmatrix} u_1^T\ u_2^T\ vdots\ u_k^T\ end{bmatrix}_{k imes n}cdot X_{n imes m}

3.2 PCA過程總結

PCA流程如下：

  1. 初始化X，使得所有樣本之間的特徵值均值為0，同時應用feature scaling，縮放到-0.5～0.5 ;
  2. 計算X的協方差矩陣S;
  3. 對S進行SVD分解，U即我們要求的新坐標系集合， Sigma 為特徵值集合（計算時特徵值都會大於0，且結果會從小到大排列）;
  4. 按照特徵值從大到小排序，要降低為k維，那麼取前k個特徵值對應的特徵向量，就是新的k個坐標軸
  5. 把X映射到新的坐標系中，完整降維操作;

根據之前的公式，做PCA投影后，投影數據的方差：

Var_{X_{project}}=sum_{j=1}^kJ_j=sum_{j=1}^klambda_j

又因為：數據從n維投影新的n維的坐標系，方差不會發生改變（向量的模長度相等且為1，可以用2D坐標系投影到45-135度坐標系驗證），即：

Var_X=Var_{X_{project}}=sum_{j=1}^nJ_j=sum_{j=1}^nlambda_j

即：X的協方差矩陣的特徵值和對應X的方差

3.3 主成份數量的選擇

PCA使得數據從n維降低為k維度，接下來介紹如何選擇合適的k。一般選擇標準為：投影前後方差比例值，作為k值的選擇標準。距離來說，我們期望：

{Var_{X_{project}}over Var_X}geq q

其中q一般選擇0.99。根據PCA總結中特徵協方差矩陣和X方差的關係得：

{Var_{X_{project}}over Var_X}={sum_{j=1}^klambda_jover sum_{j=1}^nlambda_j}geq 0.99

因此主成份數量k根據上述公式求得滿足條件的最小k

本文同時發佈於CSDN博客：

詳細推導PCA演算法（包括演算法推導必備的知識） - 怪獸 - CSDN博客?

blog.csdn.net圖標
推薦閱讀：
相关文章