主成分分析PCA(Principal Component Analysis)在sklearn中的应用及部分源码分析

最近太忙,又有一段时间没写东西了。

pca是机器学习中一个重要的降维技术,是特征提取的代表。关于pca的实现原理,在此不做过多赘述,相关参考书和各大神牛的博客都已经有各种各样的详细介绍。 如需学习相关数学理论,请移驾。T_T

简单说一下pca的实现,首先对于一个矩阵X,我们计算X·XT,显然这个一个半正定矩阵,可以做特征值分解,然后取出k个最大的特征值及其对应的特征向量就可以表达整个原矩阵。若X·XT=p-1Λp,因为p是单位矩阵,所以p-1=pT,即X·XT=p-1·Λ1/2·(p-1·Λ1/2)T,也就是降维的X后用来p-1·Λ1/2表示。

其实从SVD的角度来理解也是一样的,若X=UΣVT,则X·XT=UΣ2UT,同样我们用来UΣ来表示原X。

当我看sklearn的文档时,文档并没有具体解释它的方法得到的结果在数学上的表示什么,钻研了半天,看了源码后才知道。

sklearn的方法是通过SVD来实现的。这里着重介绍sklearn的pca类中的一个属性(components_)和两个方法(fit,transform)。

首先,给定一个矩阵,设置参数后,通过调用fit方法得到降维模型,也就是一个基矩阵。我们看一下fit中的部分关键代码。

 

...
self.mean_ = np.mean(X, axis=0)
X -= self.mean_
U, S, V = linalg.svd(X, full_matrices=False)
U, V = svd_flip(U, V)
components_ = V
...
self.components_ = components_[:n_components]
...

 

首先做的工作就是对数据进行按列中心化,然后做svd分解,然后把V的前k个向量保存为模型,模型的关键内容就是components_。(看似几百行代码,T_T)

接下来看看transfrom的部分关键源码。

 

...
X = np.dot(X, self.components_.T)
return X

 

我只贴这两句,T_T,因为真的只有这两句特别重要。也就是说我们要把一个新矩阵用到训练好的pca模型中时,其实只是做了一次矩阵乘法而已。

 

怎么来理解作者的做法呢?

其实就是我们上面的提到的,用svd方式来实现pca时,我们实际上用UΣ来表示降维后的数据。综合svd公式,可以看成XV=UΣ,也就是把原矩阵与V做乘法,实际上这里理解为投影,把X投影到V的单位正交基所表示的子空间中去得到X的低维表示。对于一个新的矩阵Y,同样用YV来表示其在V子空间降维后的结果,这也就是为什么transform方法为什么最关键的步骤只有一步乘法了。回过头看,fit训练模型就是要得到V,然后在transform降维时只需要一步乘法就可以了。

我们用下面的代码做个小实验。

 

import numpy as np
from sklearn.decomposition import PCA
from sklearn import preprocessing
from sklearn.utils.extmath import svd_flip
 
svd = np.linalg.svd
 

A = np.random.randint(0,5, (5,3))

print('A=')
print(A)

pca = PCA(n_components=2, svd_solver='full')

# print(A)
model = pca.fit(A)
print('model.components_=')
print(model.components_)
X = model.transform(A)
print('pca(A)=')
print(X)

A = A - np.mean(A, axis=0)
u,s,v = svd(A, full_matrices=False)

print('V=')
print(v[:2])
print('pva_by_svd=')
print(np.dot(A, v[:2].T))

 

运行得到的结果如下

A=
[[3 1 4]
 [4 1 1]
 [1 2 0]
 [1 4 2]
 [1 3 2]]
model.components_=
[[ 0.70734192 -0.61231721  0.35317848]
 [ 0.19273774 -0.31363722 -0.92977624]]
pca(A)=
[[ 2.21911523 -1.47640531]
 [ 1.86692171  1.50566115]
 [-1.22059974  1.54358693]
 [-1.73887721 -0.94323999]
 [-1.12656    -0.62960277]]
V=
[[-0.70734192  0.61231721 -0.35317848]
 [ 0.19273774 -0.31363722 -0.92977624]]
pva_by_svd=
[[-2.21911523 -1.47640531]
 [-1.86692171  1.50566115]
 [ 1.22059974  1.54358693]
 [ 1.73887721 -0.94323999]
 [ 1.12656    -0.62960277]]

咦,结果跟之前的分析有点小小的差别,通过svd和sklearn的结果对比,V的第一行和降维结果的第一列正负相反了。

多运行几次,这个现象并不一定出现。其实是通过上面svd_flip函数来实现的。

sklearn对奇异分解结果进行了一个处理,因为uii*vi=(-ui)*σi*(-vi),也就是u和v同时取反得到的结果是一样的,而这会导致通过pca降维得到不一样的结果(虽然都是正确的)。为了追求唯一的表示,首先定位ui向量中绝对值最大的元素位置,如果它为负数,则ui和vi取反,否则不变。这部分的源码在这里

posted @ 2017-06-13 16:17  092000  阅读(7806)  评论(0编辑  收藏  举报