avatar

主成分分析(PCA)

主成分分析(PCA)

主成分分析(Principal components analysis,简称PCA)是最重要的降维方法之一。在数据压缩消除冗余和数据噪音消除等领域都有广泛的应用。一般我们提到降维最容易想到的算法就是PCA。

PCA的思想

PCA顾名思义,就是找出数据里最主要的方面,用数据最主要的方面来替代原始的数据,具体的,加入我们的数据集是n维的,共有m个数据(x(1),x(2),,x(m))\left(x^{(1)}, x^{(2)}, \ldots, x^{(m)}\right),我们希望将这m个数据维度从n维降到n’ 维,希望m个n’ 维的数据集尽可能的代表原始的数据集。换句话来讲,就是将n维特征映射到n’维空间上(n’<n),这n’维特征是全新的正交特征,是重新构造出来的n’维特征,而不是简单地从n维特征中去除其余n−n’维特征。我们知道数据从n维降到n’ 维肯定会有损失,但是我们希望损失尽可能的小。

我们先看看最简单的情况,n=2,n’ = 1,也就是将数据从二维降到一维,数据如下图,我们希望找到某个维度的方向,他可以代表这两个维度的数据,图中,列了两个向量方向,u1u_{1}u2u_{2} ,那么那个向量可以更好的代表原始数据集呢?从直观上也可以看出,u1u_{1} 要比 u2u_{2}好。

在这里有两种解释,第一种解释是样本点到这个支线的距离足够的近,第二种解释是样本点在这个直线上的投影尽可能的分开。

假如我们把n’从一位推广到任意维,则我们的希望即降为的标准为:样本点到这个超平面的距离足够的近,或者说样本点在这个超平面上的投影尽可能的分开。

基于上边两种标准,我们可以得到PCA的两种等价的推导。

基于最小投影距离

m个n维数据(x(1),x(2),,x(m))\left(x^{(1)}, x^{(2)}, \ldots, x^{(m)}\right)(x(1),x(2),,x(m))\left(x^{(1)}, x^{(2)}, \ldots, x^{(m)}\right),都已经进行了去中心化,即i=1mx(i)=0\sum_{i=1}^{m} x^{(i)}=0,经过投影变换后得到的新坐标系为{w1,w2,,wn}\left\{w_{1}, w_{2}, \ldots, w_{n}\right\},其中w是标准正交基,即w2=1,wiTwj=0\|w\|_{2}=1, w_{i}^{T} w_{j}=0

如果我们将数据从n维降到n’ 维,即丢弃新坐标中的部分坐标,则新的坐标系为{w1,w2,,wn}\left\{w_{1}, w_{2}, \ldots, w_{n^{\prime}}\right\},样本点x(i)x^{(i)}在n’维坐标系中的投影为:z(i)=(z1(i),z2(i),,zn(i))Tz^{(i)}=\left(z_{1}^{(i)}, z_{2}^{(i)}, \ldots, z_{n^{\prime}}^{(i)}\right)^{T},其中zj(i)=wjTx(i)z_{j}^{(i)}=w_{j}^{T} x^{(i)}x(i)x^{(i)} 在低维坐标系里第 jj 维的坐标。

如果我们用 z(i)z^{(i)} 来恢复原始数据 x(i),x^{(i)}, 则得到的恢共頻据 xˉ(i)=j=1nzj(i)wj=Wz(i),\bar{x}^{(i)}=\sum_{j=1}^{n^{\prime}} z_{j}^{(i)} w_{j}=W z^{(i)}, 其中 ,, W为标准正交基組成的矩阵。

现在我们考虑整个样本集,我们希望所有的样本到这个平面的距离足够近,即最小化下面的L2:

i=1mxˉ(i)x(i)22\sum_{i=1}^{m}\left\|\bar{x}^{(i)}-x^{(i)}\right\|_{2}^{2}

对这个式子进行带入替换,打开平方公式,矩阵转置等操作:

i=1mxˉ(i)x(i)2i 译 =i=1mWz(i)x(i)22=i=1m(Wz(i))T(Wz(i))2i=1m(Wz(i))Tx(i)+i=1mx(i)Tx(i)=i=1mz(i)Tz(i)2i=1mz(i)TWTx(i)+i=1mx(i)Tx(i)=i=1mz(i)Tz(i)2i=1mz(i)Tz(i)+i=1mx(i)Tx(i)=i=1mz(i)Tz(i)+i=1mx(i)Tx(i)=tr(WT(i=1mx(i)x(i)T)W)+i=1mx(i)Tx(i)=tr(WTXXTW)+i=1mx(i)Tx(i)\begin{aligned} \sum_{i=1}^{m}\left\|\bar{x}^{(i)}-x^{(i)}\right\|_{2}^{i \text { 译 }} &=\sum_{i=1}^{m}\left\|W z^{(i)}-x^{(i)}\right\|_{2}^{2} \\ &=\sum_{i=1}^{m}\left(W z^{(i)}\right)^{T}\left(W z^{(i)}\right)-2 \sum_{i=1}^{m}\left(W z^{(i)}\right)^{T} x^{(i)}+\sum_{i=1}^{m} x^{(i) T} x^{(i)} \\ &=\sum_{i=1}^{m} z^{(i) T} z^{(i)}-2 \sum_{i=1}^{m} z^{(i) T} W^{T} x^{(i)}+\sum_{i=1}^{m} x^{(i) T} x^{(i)} \\ &=\sum_{i=1}^{m} z^{(i) T} z^{(i)}-2 \sum_{i=1}^{m} z^{(i) T} z^{(i)}+\sum_{i=1}^{m} x^{(i) T} x^{(i)} \\ &=-\sum_{i=1}^{m} z^{(i) T} z^{(i)}+\sum_{i=1}^{m} x^{(i) T} x^{(i)} \\ &=-t r\left(W^{T}\left(\sum_{i=1}^{m} x^{(i)} x^{(i) T}\right) W\right)+\sum_{i=1}^{m} x^{(i) T} x^{(i)} \\ &=-t r\left(W^{T} X X^{T} W\right)+\sum_{i=1}^{m} x^{(i) T} x^{(i)} \end{aligned}

其中,注意到 i=1mx(i)x(i)T\sum_{i=1}^{m} x^{(i)} x^{(i) T} 是数据集的协方差矩阵, \quad W的每一个向量 wj2w_{j^{2}} 是标准正交基。而 i=1mx(i)Tx(i)\sum_{i=1}^{m} x^{(i) T} x^{(i)} 是一个常量。最小化上式等价于

argminWtr(WTXXTW)\underbrace{\arg \min }_{W}-\operatorname{tr}\left(W^{T} X X^{T} W\right) s.t. WTW=IW^{T} W=I

这个最小化不难,直接观察也可以发现最小值对应的W由协方差矩阵 XXTX X^{T} 最大的n’个特征值对应的特征向量组成。当然用数学推导也很容易。利用拉格
朗日函数可以得到

J(W)=tr(WTXXTW+λ(WTWI))J(W)=-t r\left(W^{T} X X^{T} W+\lambda\left(W^{T} W-I\right)\right)

对W求导有- XXTW+λW=0,X X^{T} W+\lambda W=0, 整理下即为:

XXTW=λWX X^{T} W=\lambda W

这样可以更清楚的看出, W为 XXTX X^{T} 的n’个特征向量组成的矩阵, 而 λ\lambdaXXTX X^{T} 的若干特征值组成的矩阵, 特征值在主对角线上, 其余位置为0。当我们
将数据集从n维降到n’维时,需要找到最大的n’个特征值对应的特征向量。这n’个特征向量组成的矩阵W即为我们需要的矩阵。对于原始数据集,我们只需要用 z(i)=WTx(i)z^{(i)}=W^{T} x^{(i)},就可以把原始此据集降维到最小投影距离的n’维数据集。

基于最大投影方差

假设m个n维数据 (x(1),x(2),,x(m))\left(x^{(1)}, x^{(2)}, \ldots, x^{(m)}\right) 都已经进行了中心化, \quadi=1mx(i)=0.\sum_{i=1}^{m} x^{(i)}=0 . 经过投影变换后得到的新坐标系为 {w1,w2,,wn},\left\{w_{1}, w_{2}, \ldots, w_{n}\right\}, 其中 w是标准正交 w_{\text {是标准正交 }} 基, 即 ||w2=1,wiTwj=0w||_{2}=1, w_{i}^{T} w_{j}=0
如果我们将数据从n维降到n’维, 即丟弃新坐标系中的部分坐标, 则新的坐标系为 {w1,w2,,wn},\left\{w_{1}, w_{2}, \ldots, w_{n^{\prime}}\right\}, 样本点 x(i)x^{(i)} 在n’维坐标系中的投影为
[z{(i)}=\left(z_{1}{(i)}, z_{2}^{(i)}, \ldots, z_{n{\prime}}{(i)}\right)^{T} \text { .其中, } z_{j}{(i)}=w_{j}{T} x^{(i)} \text { 是 } x^{(i)} \text { 在低维坐标系里第 } j \text { 维的坐标。 }]
对于任意一个样本 x(i),x^{(i)}, 在新的坐标系中的投影为 WTx(i)W^{T} x^{(i)},在新坐标系中的投影方差为 x(i)TWWTx(i),x^{(i) T} W W^{T} x^{(i)}, 要使所有的样本的投影方差和最大, 也就是 最大化 i=1mWTx(i)x(i)TW\sum_{i=1}^{m} W^{T} x^{(i)} x^{(i) T} W 的迹,即:

argmaxWtr(WTXXTW)\underbrace{\arg \max }_{W} \operatorname{tr}\left(W^{T} X X^{T} W\right) s.t. WTW=IW^{T} W=I

利用拉格朗日函数可以得到:

J(W)=tr(WTXXTW+λ(WTWI))J(W)=\operatorname{tr}\left(W^{T} X X^{T} W+\lambda\left(W^{T} W-I\right)\right)

对W求导有 XXTW+λW=0,X X^{T} W+\lambda W=0, 整理下即为:
[X X^{T} W=(-\lambda) W]

和上面一样可以看出, W为X XTX^{T} 的n’个特征向量组成的矩阵, 而一 λ\lambdaXXTX X^{T} 的若干特征值组成的矩阵, 特征值在主对角线上, 其余位置为0。当我们
将数据集从n维降到n’维时,需要找到最大的n’个特征值对应的特征向量。这n’个特征向量组成的矩阵W即为我们需要的矩阵。对于原始数据集,我们只需要用
z(i)=WTx(i)z^{(i)}=W^{T} x^{(i)},就可以把原始数据集降维到最小投影距离的n’维数据集。

PCA算法流程

Input:n维样本集D,要降维的维数n’

output:降为后的样本集D’

1)对所有的样本进行中心化:x(i)=x(i)1mj=1mx(j)x^{(i)}=x^{(i)}-\frac{1}{m} \sum_{j=1}^{m} x^{(j)}

2)计算样本的协方差矩阵=XXT=X X^{T}

3)对矩阵=XXT=X X^{T}进行特征值分解

4)取出最大的n’个特征值向量=(w1,w2,,wn)=\left(w_{1}, w_{2}, \ldots, w_{n^{\prime}}\right),将所有的特征向量标准化后,组成特征向量矩阵W

5)对样本几种的每一个样本x(i)x^{(i)}转化为新的样本z(i)=WTx(i)z^{(i)}=W^{T} x^{(i)}

6)得到输出的样本集D’

有时候,我们不指定降维后的n’的值, 而是换种方式, 指定一个降维到的主成分比重闽值t。这个闽值t在 (0,1]之间。假如我们的n个特征值为 λ1λ2λn,\lambda_{1} \geq \lambda_{2} \geq \ldots \geq \lambda_{n}, 则n’可以通过下式得到:

i=1nλii=1nλit\frac{\sum_{i=1}^{n^{\prime}} \lambda_{i}}{\sum_{i=1}^{n} \lambda_{i}} \geq t

简单案例实现

#主成分分析
import numpy as np
#去中心化
def Mean(data_matrix):
#求均值
mean_val = np.mean(data_matrix,axis = 0)
#去中心化
new_data = data_matrix - mean_val
return new_data,mean_val

def PCA(data_matrix,n):
#得到去中心化后的矩阵及平均值
new_data ,mean_val = Mean(data_matrix)
#求协方差,rowvar = 0,传入的数据一行代表一个样本,rowvar非0,一列代表一个样本
cov_matrix = np.cov(new_data,rowvar = 0)
#计算特征值和特征矩阵,利用numpy.linalg
eigenvalues,eigenVects = np.linalg.eig(np.mat(cov_matrix))

#保留主要成分前n个,组成了新的特征空间的一组基n_eigVects
eigValIndice = np.argsort(eigenvalues)#对特征值从小到大排序

n_eigValIndice = eigValIndice[-1:-(n-1):-1]#取最大的n个特征值的下标

n_eigVect = eigenVects[:,n_eigValIndice]#最大n个特征值的特征向量

low_data_mat = new_data * n_eigVect#低维空间数据
re_mat = (low_data_mat * n_eigVect.T) + mean_val#重构数据
return low_data_mat,re_mat
#根据阈值求n
def percentage2n(eigenvalues,percentage):
sortArray = np.sort(eigenvalues)
sortArray = sortArray[-1::-1]
arraySum = sum(sortArray)
tmp = 0
num = 0
for i in sortArray:
tmp+=i
num+=1
if tmpSum >= arraySum*percentage:
return num
#指定一个降维到的主成分比重阈值percentage=0.99
def p_PCA(data_matrix,percentage=0.99):
#得到去中心化后的矩阵及平均值
new_data ,mean_val = Mean(data_matrix)
#求协方差,rowvar = 0,传入的数据一行代表一个样本,rowvar非0,一列代表一个样本
cov_matrix = np.cov(new_data,rowvar = 0)
#计算特征值和特征矩阵,利用numpy.linalg
eigenvalues,eigenVects = np.linalg.eig(np.mat(cov_matrix))
#计算所降的维度n
n = percentage2n(eigenvalues,percentage)
#保留主要成分前n个,组成了新的特征空间的一组基n_eigVects
egiValIndice = np.argsort(eigenvalues)#对特征值从小到大排序
n_eigValIndice = egiValIndice[-1:-(n-1):-1]#取最大的n个特征值的下标
n_eigVect = eigenVects[:,n_eigValIndice]#最大n个特征值的特征向量

low_data_mat = new_data * n_eigVect#低维空间数据
re_mat = (low_data_mat * n_eigVect.T) + mean_val#重构数据

return low_data_mat,re_mat

if __name__ == "__main__":
#data
data_matrix = [[-0.2314987 , 0.08387106],
[ 0.66963671 , 1.38535319],
[ 0.97355908 , 0.04134561],
[ 0.38108224 , 0.9434845 ],
[ 0.11467758 ,-0.72613803]]
#PCA方法一
low_data,re_mat = PCA(data_matrix,1)
print(low_data,'\n',re_mat)
#PCA方法二
low_data1,re_mat1 =p_PCA(data_matrix,percentage=0.99)
print(low_data1,'\n',re_mat1)

最后

作为一个非监督学习的降维方法,它只需要特征值分解,就可以对数据进行压缩,去噪。为了克服PCA的一些缺点,出现了很多PCA的变种,比如解决非线性降维的KPCA,还有解决内存限制的增量PCA方法Incremental PCA,以及解决稀疏数据降维的PCA方法Sparse PCA等。

PCA算法的主要优点有:

1)仅仅需要以方差衡量信息量,不受数据集以外的因素影响。

2)各主成分之间正交,可消除原始数据成分间的相互影响的因素。

3)计算方法简单,主要运算是特征值分解,易于实现。

PCA算法的主要缺点有:

1)主成分各个特征维度的含义具有一定的模糊性,不如原始样本特征的解释性强。

2)方差小的非主成分也可能含有对样本差异的重要信息,因降维丢弃可能对后续数据处理有影响。

Author: Hui Ning
Link: https://angelni.github.io/PCA/
Copyright Notice: All articles in this blog are licensed under CC BY-NC-SA 4.0 unless stating additionally.
Donate
  • 微信
    微信
  • 支付宝
    支付宝

Comment