1. 背景
在研究多变量的数据时,变量太多不仅导致数据的复杂性,而且不同的变量之间可能存在一定的相关关系,即两个变量之间的信息有一定的重叠,这无疑为分析问题增加了难度。为了提取数据中的主要信息成分,提出了一个新的降维方法PCA。主成分分析(Principal Component Analysis,简称PCA)是最常用的一种降维方法。它是将原先提出的所有变量中关系紧密的变量删除,建立尽可能少的新变量,使得这些新变量两两不相关,而且这些新变量尽可能的包含原有的信息。
2. 基础知识
协方差与散度矩阵
样本均值 x ˉ = 1 n ∑ i = 1 N x i \bar x=\frac{1}{n}\sum_{i=1}^Nx_i xˉ=n1i=1∑Nxi样本方差 S 2 = 1 n − 1 ∑ i = 1 N ( x i − x ˉ ) 2 S^2=\frac{1}{n-1}\sum_{i=1}^N(x_i-\bar x)^2 S2=n−11i=1∑N(xi−xˉ)2样本X和样本Y的协方差 c o v ( X , Y ) = E [ ( X − E ( X ) ) ( Y − E ( Y ) ) ] = 1 n − 1 ∑ i = 1 n ( x i − x ˉ ) ( y i − y ˉ ) \begin {aligned}cov(X,Y)&=E[(X-E(X))(Y-E(Y))] \\ & = \frac{1}{n-1}\sum_{i=1}^n(x_i-\bar x)(y_i-\bar y)\end {aligned} cov(X,Y)=E[(X−E(X))(Y−E(Y))]=n−11i=1∑n(xi−xˉ)(yi−yˉ)假如样本X和Y的维度为2维。那么它们之间的协方差为: c o v ( X , Y ) = [ c o v ( x , x ) c o v ( x , y ) c o v ( y , x ) c o v ( y , y ) ] cov(X,Y)=\begin{bmatrix} cov(x,x)& cov(x,y) \\ cov(y,x) & cov(y,y) \end{bmatrix} cov(X,Y)=[cov(x,x)cov(y,x)cov(x,y)cov(y,y)]如果协方差为正,说明X和Y是正相关关系;协方差为负,说明X和Y是负相关关系;协方差为0,说明X与Y不相关。
散度矩阵 S = ∑ i = 1 n ( x i − m ) ( x i − m ) T S=\sum_{i=1}^n(\boldsymbol {x_i}-\boldsymbol{m})(\boldsymbol{x_i}-\boldsymbol{m})^T S=i=1∑n(xi−m)(xi−m)T m = 1 n ∑ i = 1 n x i \boldsymbol{m}=\frac{1}{n}\sum_{i=1}^n\boldsymbol{x_i} m=n1i=1∑nxi数据X的散度矩阵为 X X T XX^T XXT。散度矩阵等于协方差矩阵乘以(总数据量-1)。因此,它们的特征值与特征向量是相同的。
3. PCA
假设数据集中有m个数据{ x 1 , x 2 , . . . , x m x_1,x_2,...,x_m x1,x2,...,xm},每一个数据是d维数据。目标是在尽可能保证原有数据的信息条件下,将数据的维度从 d d d维降到 d ′ d^{'} d′维。虽然对数据降维会导致原有数据信息的损失,但是希望是使得损失尽可能的小。那么,如何将 d d d维数据降到 d ′ d^{'} d′维度呢?
举一个例子,如下图是将二维数据降到一维。希望是找到一个维度方向,可以代表这两个维度的数据。图中,有两个维度方向,
u
1
u_1
u1和
u
2
u_2
u2。那么,哪一个向量可以更好的代表原有数据呢?
直观来讲,
u
1
u_1
u1比
u
2
u_2
u2更适合作为一维方向。原因可从以下两个方面来解释:
- 最近重构性:样本点到这个直线的距离都足够近;
- 最大可分性:样本点在这个直线上的投影尽可能分开;
3.1 PCA的推导方法:最近重构性
假设m个d维数据{
x
1
,
x
2
,
.
.
.
,
x
m
x_1,x_2,...,x_m
x1,x2,...,xm},首先进行数据中心化,即
∑
i
=
1
m
x
i
=
0
\sum_{i=1}^mx_i=0
∑i=1mxi=0,坐标系为{
w
1
,
w
2
,
.
.
.
,
w
d
w_1,w_2,...,w_d
w1,w2,...,wd},其中w是标准正交基,即满足
∣
∣
w
∣
∣
2
=
1
,
w
i
T
w
j
=
0
||w||_2=1,w_i^Tw_j=0
∣∣w∣∣2=1,wiTwj=0。如果将数据从
d
d
d维降到
d
′
d^{'}
d′维,那么在新坐标系中的投影为
z
i
=
(
z
i
1
,
z
i
2
,
.
.
.
,
z
i
d
′
)
z_i=(z_{i1},z_{i2},...,z_{id^{'}})
zi=(zi1,zi2,...,zid′),其中
z
i
j
=
w
j
T
x
i
z_{ij}=\mathbf w_j^Tx_i
zij=wjTxi是
x
i
x_i
xi在低维坐标系下第j维的坐标。如果用
z
i
z_i
zi来重构数据
x
i
x_i
xi,那么恢复数据为
x
^
i
=
∑
j
=
1
d
′
z
i
j
w
j
=
W
z
i
\hat {x}_i=\sum_{j=1}^{d^{'}}z_{ij}\mathbf w_j=Wz_i
x^i=∑j=1d′zijwj=Wzi,其中W为标准正交基组成的矩阵。
考虑到整个数据集,希望所有样本到这个超平面的距离尽可能的小,即
∑
i
=
1
m
∥
∑
j
=
1
d
′
z
i
j
w
j
−
x
i
∥
2
2
=
∑
i
=
1
m
∣
∣
W
z
i
−
x
i
∣
∣
2
2
=
∑
i
=
1
m
z
i
T
W
T
W
z
i
−
2
∑
i
=
1
m
(
W
z
i
)
T
x
i
+
∑
i
=
1
m
x
i
T
x
i
=
∑
i
=
1
z
i
T
z
i
−
2
∑
i
=
1
m
z
i
T
z
i
+
∑
i
=
1
m
x
i
T
x
i
=
−
∑
i
=
1
m
z
i
T
z
i
+
∑
i
=
1
m
x
i
T
x
i
=
−
t
r
(
Z
T
Z
)
+
∑
i
=
1
m
x
i
T
x
i
=
−
t
r
(
W
T
X
X
T
W
)
+
∑
i
=
1
m
x
i
T
x
i
\begin {aligned}\sum_{i=1}^m\left \| \sum_{j=1}^{d^{'}}z_{ij}w_j-x_i \right \|_2^2 &=\sum_{i=1}^m||Wz_i-x_i||^2_2 \\ &=\sum_{i=1}^mz_i^TW^TWz_i-2\sum_{i=1}^m(Wz_i)^Tx_i + \sum_{i=1}^mx_i^Tx_i \\ &=\sum_{i=1}z_i^Tz_i-2\sum_{i=1}^mz_i^Tz_i+\sum_{i=1}^mx_i^Tx_i \\ &=-\sum_{i=1}^mz_i^Tz_i+\sum_{i=1}^mx_i^Tx_i \\ &=-tr(Z^TZ)+\sum_{i=1}^mx_i^Tx_i \\ &=-tr(W^TXX^TW)+\sum_{i=1}^mx_i^Tx_i\end {aligned}
i=1∑m∥∥∥∥∥∥j=1∑d′zijwj−xi∥∥∥∥∥∥22=i=1∑m∣∣Wzi−xi∣∣22=i=1∑mziTWTWzi−2i=1∑m(Wzi)Txi+i=1∑mxiTxi=i=1∑ziTzi−2i=1∑mziTzi+i=1∑mxiTxi=−i=1∑mziTzi+i=1∑mxiTxi=−tr(ZTZ)+i=1∑mxiTxi=−tr(WTXXTW)+i=1∑mxiTxi其中,用到的公式有
(
A
B
)
T
=
B
T
A
T
,
W
T
W
=
I
,
z
i
=
W
T
x
i
(AB)^T=B^TA^T,W^TW=I,z_i=W^Tx_i
(AB)T=BTAT,WTW=I,zi=WTxi。由于
∑
i
=
1
m
x
i
T
x
i
\sum_{i=1}^mx_i^Tx_i
∑i=1mxiTxi是数据集的协方差矩阵,为常数。W的每一个向量
w
j
\mathbf w_j
wj是标准正交基,因此上式可转化为:
a
r
g
min
W
−
t
r
(
W
T
X
X
T
W
)
s
.
t
.
W
T
W
=
I
arg\min_{W} \quad -tr(W^TXX^TW) \quad s.t. \quad W^TW=I
argWmin−tr(WTXXTW)s.t.WTW=I通过拉格朗日函数推导得到:
J
(
W
)
=
−
t
r
(
W
T
X
X
T
W
+
λ
(
W
T
W
−
I
)
)
J(W)=-tr(W^TXX^TW+\lambda(W^TW-I))
J(W)=−tr(WTXXTW+λ(WTW−I))对W求导得到:
X
X
T
W
=
λ
W
XX^TW=\lambda W
XXTW=λW此时可以看出,W为
X
T
X
X^TX
XTX的
d
′
d^{'}
d′个特征向量组成的矩阵,
λ
\lambda
λ为
X
T
X
X^TX
XTX矩阵的特征值组成的矩阵,特征值在主对角线上,非对角线上的元素为0。当数据集从
d
d
d维降到
d
′
d^{'}
d′时,只需按照特征值的由大到小顺序找到前
d
′
d^{'}
d′个特征值及对应的特征向量,这些特征向量组成的矩阵就是我们要求得矩阵。再通过矩阵变换
z
i
=
W
T
x
i
z_i=W^Tx_i
zi=WTxi将原始数据映射到
d
′
d^{'}
d′维空间中。
3.2 PCA的另一推导方法:最大可分性
假设m个d维数据{ x 1 , x 2 , . . . , x m x_1,x_2,...,x_m x1,x2,...,xm},首先进行数据中心化,即 ∑ i = 1 m x i = 0 \sum_{i=1}^mx_i=0 ∑i=1mxi=0,坐标系为{ w 1 , w 2 , . . . , w d w_1,w_2,...,w_d w1,w2,...,wd},其中w是标准正交基,即满足 ∣ ∣ w ∣ ∣ 2 = 1 , w i T w j = 0 ||w||_2=1,w_i^Tw_j=0 ∣∣w∣∣2=1,wiTwj=0。如果将数据从 d d d维降到 d ′ d^{'} d′维,那么在新坐标系中的投影为 z i = ( z i 1 , z i 2 , . . . , z i d ′ ) z_i=(z_{i1},z_{i2},...,z_{id^{'}}) zi=(zi1,zi2,...,zid′),其中 z i j = w j T x i z_{ij}=\mathbf w_j^Tx_i zij=wjTxi是 x i x_i xi在低维坐标系下第j维的坐标。
从最大可分性出发,对于任意一个样本 x i x_i xi,在新空间超平面上的投影为 W T x i W^Tx_i WTxi,在新坐标系中的投影方差为 W T x i x i T W W^Tx_ix_i^TW WTxixiTW。若所有样本的投影尽可能分开,则应该使得投影后的样本方差最大化。即最大化 ∑ i = 1 m W T x i x i T W \sum_{i=1}^mW^Tx_ix_i^TW ∑i=1mWTxixiTW a r g max W t r ( W T X X T W ) s . t . W T W = I arg\max_W \quad tr(W^TXX^TW) \quad s.t. \quad W^TW=I argWmaxtr(WTXXTW)s.t.WTW=I和上一个方法相同,用拉格朗日函数得到: J ( W ) = t r ( W T X X T W + λ ( W T W − I ) ) J(W)=tr(W^TXX^TW+\lambda(W^TW-I)) J(W)=tr(WTXXTW+λ(WTW−I))对W求导得到: X X T W = − λ W XX^TW=-\lambda W XXTW=−λW类似,W为 X T X X^TX XTX的 d ′ d^{'} d′个特征向量组成的矩阵, λ \lambda λ为 X T X X^TX XTX矩阵的特征值组成的矩阵,特征值在主对角线上,非对角线上的元素为0。当数据集从 d d d维降到 d ′ d^{'} d′时,只需按照特征值的由大到小顺序找到前 d ′ d^{'} d′个特征值及对应的特征向量,这些特征向量组成的矩阵就是我们要求得矩阵。再通过矩阵变换 z i = W T x i z_i=W^Tx_i zi=WTxi将原始数据映射到 d ′ d^{'} d′维空间中。
4. PCA算法
输入:
d
d
d维样本集
D
=
{
x
1
,
x
2
,
.
.
.
,
x
m
}
D=\left\{x_1,x_2,...,x_m\right\}
D={x1,x2,...,xm}
输出:降维后的样本集
D
′
D^{'}
D′
- 对所有样本进行中心化: x i = x i − 1 m ∑ i = 1 m x i x_i=x_i-\frac{1}{m}\sum_{i=1}^mx_i xi=xi−m1∑i=1mxi
- 计算样本的协方差矩阵 X X T XX^T XXT
- 对协方差矩阵 X X T XX^T XXT进行特征值分解
- 取最大的 d ′ d^{'} d′个特征值所对应的特征向量 { w 1 , w 2 , . . . , w d ′ } \left\{w_1,w_2,...,w_{d^{'}}\right\} {w1,w2,...,wd′},对所有特征向量进行标准化后,组成特征向量矩阵W
- 对每一个样本 x i x_i xi,投影得到新的样本 z i = W T x i z_i=W^Tx_i zi=WTxi
- 得到输出样本集 D ′ = { z 1 , z 2 , . . . , z m } D^{'}=\left\{z_1,z_2,...,z_m\right\} D′={z1,z2,...,zm}
对于降维后的维数 d ′ d^{'} d′可由用户事先指定,也可通过设置一个主成分比重阈值 t ∈ ( 0 , 1 ] t\in (0,1] t∈(0,1]。假设特征值为 λ 1 ≥ λ 2 ≥ . . . ≥ λ n \lambda_1\ge \lambda_2\ge ...\ge \lambda_n λ1≥λ2≥...≥λn则选取使下式成立的最小 d ′ d^{'} d′值: ∑ i = 1 d ′ λ i ∑ i = 1 d λ i ≥ t \frac{\sum_{i=1}^{d^{'}}\lambda_i}{\sum_{i=1}^d \lambda_i}\ge t ∑i=1dλi∑i=1d′λi≥t
5. 核化线性降维
如果数据是线性的,线性降维方法可将数据从高维空间映射到低维空间。然而,一般情况下,数据是非线性的,不能使用PCA线性降维。因此,可先将n维数据映射到更高维N的空间,然后再从N维空间降维到低维空间n’。
非线性降维是一种常用的方法,核函数的主成分分析称为核主成分分析KPCA。假设n维空间到高维空间的映射函数为
ϕ
\phi
ϕ,n维空间的特征分解为:
∑
i
=
1
m
x
i
x
i
T
W
=
λ
W
\sum_{i=1}^mx_ix_i^TW=\lambda W
i=1∑mxixiTW=λW高维空间的特征分解为:
∑
i
=
1
m
ϕ
(
x
i
)
ϕ
(
x
i
)
T
W
=
∑
i
=
1
m
k
(
x
i
,
x
i
)
W
=
λ
W
\sum_{i=1}^m\phi(x_i)\phi(x_i)^TW=\sum_{i=1}^m\boldsymbol k(x_i,x_i)W=\lambda W
i=1∑mϕ(xi)ϕ(xi)TW=i=1∑mk(xi,xi)W=λW