sklearn中的降维算法PCA和SVD_在scikit-learn中有sklearn.decomposition.kernelpca()函数

1 概述

1.1 什么叫“维度”

对于数组和Series来说，维度就是功能shape返回的结果，shape中返回了几个数字，就是几维。索引以外的数据，不分行列的叫一维（此时shape返回唯一的维度上的数据个数），有行列之分叫二维（shape返回行*列），也称为表。一张表最多二维，复数的表构成了更高的维度。但一个数组中存在2张3行4列的表时，shape返回的是（2,3,4）。
数组中的每一张表，都可以是一个特征矩阵（矩阵可以是任意正数维，但特征矩阵特指二维）或一个DataFrame，这些结构永远只有一张表，所以一定有行列，其中行是样本，列是特征。针对每一张表，维度指的是样本的数量或特征的数量，一般无特别说明，指的都是特征的数量。除了索引之外，n个特征就是n维特征矩阵。
对图像来说，维度就是图像中特征向量的数量。特征向量可以理解为坐标轴，一个特征向量定义一条直线，是一维，两个相互垂直的特征向量定义一个平面，即一个直角坐标系，就是二维，三个相互垂直的特征向量定义一个空间，即一个立体直角坐标系，就是三维。三个以上的特征向量相互垂直，定义是人眼无法看见的，也无法想象的高维空间。
降维算法中的“降维”，指的是降低特征矩阵中特征的数量。降维的目的是为了让算法运算更快，效果更好。但还有另一种需求：数据可视化。

1.2 sklearn中的降维算法

sklearn中降维算法都被包括在模块decomposition中，这个模块本质是一个矩阵分解模块。

• 主成分分析
  decomposition.PCA——主成分分析（PCA）
  decomposition.IncrementalPCA——增量主成分分析（IPCA）
  decomposition.KernelPCA——核主成分分析（KPCA）
  decomposition.MiniBathchSparsePCA——小批量稀疏主成分分析
  decomposition.SparsePCA——稀疏主成分分析（SparsePCA）
  decomposition.TruncatedSVD——截断的SVD（aka LSA）
• 因子分析
  decomposition.FactorAnalysis——因子分析（FA）
• 独立成分分布
  decomposition.FastICA——独立成分分析的快速算法
• 字典学习
  decomposition.DictionaryLearning——字典学习
  decomposition.MiniBatchDictionaryLearning——小批量字典学习
  decomposition.dict_learning——字典学习用于矩阵分解
  decomposition.dict_learning_online——在线字典学习用于矩阵分解
• 高级矩阵分解
  decomposition.LatterDirichletAllocation——具有在线变分贝叶斯算法的隐含狄利克雷分布（LDA）
  decomposition.NMF——非负矩阵分解（NMF）
• 其他矩阵分解
  decomposition.SparseCoder——稀疏编码
  SVD和主成分分析PCA都属于矩阵分解算法中的入门算法，都是通过分解特征矩阵来进行降维。

2 PCA和SVD

在降维过程中，会减少特征的数量，这意味着删除数据，数据量变少则表示模型可以获取的信息会变少，模型的表现可能会因此受到影响。同时，在高维数据中，必然有一些特征是不带有有效的信息的（如噪音），或者有些特征带有的信息和其他一些特征是重复的（如多重共线性）。希望能够找到一种办法来衡量特征上所带的信息量，使得在降维过程中，能够既减少特征的数量，又保留大部分有效信息。将那些带有重复信息的特征合并，并删除哪些带无效信息的特征等等，创造出能够代表原特征矩阵大部分信息的、特征更少的新特征矩阵。
特征工程中有一种重要的特征选择方法：方差过滤。如果一个特征的方差很小，则意味着这个特征很可能有大量取值都想要，那这个特征的取值对样本而言就没有区分度，这种特征就不带有有效信息。从方差的这一应用可以推断出，如果一个特征的方差很大，则说明这个特征上带有大量的信息。因此，在降维中，PCA使用的信息量衡量指标，就是样本方差，又称可解释性方差，方差越大，特征所带的信息量越多。
$Var=\frac{1}{n-1}\sum (x_i-\hat{x}{)}^2$

2.1 降维究竟怎样实现

sklearn.decomposition.PCA(n_estimators = None, copy = True, whiten = False, svd_solver = ‘auto’, tol = 0.0, iterated_power = ‘auto’, random_state = None)
n维特征矩阵降维的步骤：

• 输入原数据，结构为(m,n)，找出原本的n个特征向量构成的n维空间V
• 决定降维后的特征数量k
• 通过某种变化，找出n个新的特征向量，以及他们构成的新n维空间V
• 找出原始数据在新特征空间V中的n个新特征向量上对应的值，即“将数据映射到新空间中”
• 选取前k个信息量最大的特征，删除没有被选中的特征，成功将n维空间V降为k维
  在步骤3中，用来找出n个新特征向量，让数据能够被压缩到少数特征上并且总信息量不损失太多的技术就是矩阵分解。PCA和SVD是两种不同的降维算法，但两者都遵从上述步骤来实现降维，只是两种算法中矩阵分解的方法不同，信息量的衡量指标不同。
  PCA使用方差作为信息量的衡量指标，并且通过特征值分解来找出空间V。降维时将特征矩阵$X$分解为$Q\sum Q^{-1}$，其中$Q$和$Q^{-1}$是辅助的矩阵，$\sum$是一个对角矩阵，其对角线上的元素就是方差。降维完成之后，PCA找到的每个新特征向量就叫做“主成分”，而被丢弃的特征向量被认为信息量很少，这些信息很可能就是噪音。
  SVD使用奇异值分解来找出空间V，将特征矩阵$X$分解为$U\sum V^T$，其中$\sum$也是一个对角矩阵，不过其对角线上的元素是奇异值，这也是SVD中用来衡量特征上的信息量的指标。$U$和$V^T$分别是左奇异矩阵和右奇异矩阵，也都是辅助矩阵。
  在数学原理中，无论是PCA和SVD都需要遍历所有的特征和样本来计算信息量指标。并且在矩阵分解的过程之中，会产生比原来的特征矩阵更大的矩阵（如原数据的结构是(m,n)，在矩阵分解中为了找出最佳新特征空间V，可能需要产生(n,n)，(m,m)大小的矩阵，还需要产生协方差矩阵去计算更多的信息）。而无论是Python还是R，或者其他任何语言，在大型矩阵运算上都不是特别擅长，无论代码如何简化，都需要等待计算机去完成这个非常庞大的数学运算过程，因此，降维算法的计算量很大，运行比较缓慢。
  **思考：PCA和特征选择技术都是特征工程的一部分，它们有什么不同？**特征工程中有是三种方式：特征提取、特征创造和特征选择。特征选择是从已存在的特征中选取携带信息最多的，选完之后的特征依然具有可解释性，依然知道选出的特征在原数据的哪个位置，代表着原数据上的什么含义。而PCA是将已存在的特征进行压缩，降维完毕就的特征不是原本的特征矩阵中的任何一个特征，而是通过某些方式组合起来的新特征。通常来说，在新的特征矩阵生成之前，无法知晓PCA都建立了怎样的新特征向量，新特征矩阵生成之后也不具有可读性，无法判断新特征矩阵的特征是从原数据中的什么特征组合而成，新特征虽然带有原始数据的信息，却已经不是原数据上代表着的含义了。以PCA为代表的降维算法因此是特征创造（feature creation或feature construction）的一种。因此，PCA一般不适用于探索特征和标签之间的关系的模型（如线性回归），无法解释的新特征和标签之间的关系不具有意义。在线性回归模型中，使用特征选择。

2.2 重要参数n_components

n_components是降维后需要的维度，即降维后需要保留的特征数量。降维流程步骤2中需要确定的k值，一般输入$[0,min(X.shape)]$范围中的整数。这是一个需要认为去确认的超参数，并且设定的数值会影响到模型的表现。如果留下的特征太多，就达不到降维的效果，如果留下的特征太少，那么新特征向量可能无法包含原始数据集中的大部分信息。
从降维的目标而言，如果希望可视化一组数据来观察数据分布，往往将数据降到三维以下 ，很多时候是二维，即n_components的取值为2。

2.2.1 迷你案例：高维数据的可视化

#调用库和模块
import matplotlib.pyplot as plt
from sklearn.datasets import load_iris
from sklearn.decomposition import PCA
#提取数据集
iris = load_iris()
y = iris.target
x = iris.data
x.shape
import pandas as pd
pd.DataFrame(x)
#建模
#调用PCA
pca = PCA(n_components = 2)#实例化
pca = pca.fit(x)#拟合模型
x_dr = pca.transform(x)#获取新矩阵
x_dr
#也可以使用fit_transform一步到位
#x_dr = PCA(2).fit_transform(x)
#要展示三种分类的分布，需要对三种鸢尾花分别绘图
#可以写成三行代码，也可以写成for循环
#plt.figure()#初始化画布
#plt.scatter(x_dr[y == 0,0],x_dr[y == 0,1],c = 'red',label = iris.target_names[0])
#plt.scatter(x_dr[y == 1,0],x_dr[y == 1,1],c = 'black',label = iris.target_names[1])
#plt.scatter(x_dr[y == 2,0],x_dr[y == 2,1],c = 'orange',label = iris.target_names[2])
#plt.legend()
#plt.title("PCA of IRIS dataset")
#plt.show()
colors = ["red","black","orange"]
iris.target_names
plt.figure()
for i in [0,1,2]:
    plt.scatter(x_dr[y == i,0]
               ,x_dr[y == i,1]
               ,alpha = .7
               ,c = colors[i]    #透明度
               ,label = iris.target_names[i]
               )
plt.legend()
plt.title("PCA of IRIS dataset")
plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41

#探索降维后的数据
#属性explained_variance_，查看降维后每个新特征向量上所带的信息量大小（可解释性方差的大小）
pca.explained_variance_
#array([4.22824171, 0.24267075])
#属性explained_variance_ratio_，查看降维后每个新特征向量所占的信息量占原始数据总信息量的百分比,又叫做可解释方差贡献率
pca.explained_variance_ratio_
#array([0.92461872, 0.05306648])
#大部分信息都被有效地集中在了第一个特征上
pca.explained_variance_ratio_.sum()
#0.9776852063187949
1
2
3
4
5
6
7
8
9
10

当参数n_components中不填写任何值，则默认返回min(X.shape)个特征，一般来说，样本量都会大于特征数目，所以什么都不填就相当于转换了新特征空间，但没有减少特征的个数。一般来说，不会使用这种输入方式，但却可以使用这种输入方式来画出累积可解释方差贡献率曲线，以此来选择最好的n_components的整数取值。
累积可解释方差贡献率曲线是一条以降维后保留的特征个数为横坐标，降维后新特征矩阵捕捉到的可解释方差贡献率为纵坐标的曲线，能够帮助确定n_components最好的取值。

#选择最好的n_components：累积可解释方差贡献率曲线
import numpy as np
pca_line = PCA().fit(x)
plt.plot([1,2,3,4],np.cumsum(pca_line.explained_variance_ratio_))
plt.xticks([1,2,3,4])#这是为了限制坐标轴显示为整数
plt.xlabel("number of components after dimension reduction")
plt.ylabel("cumulative explained variance")
plt.show()
1
2
3
4
5
6
7
8

2.2.2 最大似然估计自选超参数

让PCA用最大似然估计（maximum likelihood estimation）自选超参数的方法，输入“mle”作为n_components的参数输入，就可以调用这个方法。

pca_mle = PCA(n_components = 'mle')
pca_mle = pca_mle.fit(x)
x_mle = pca_mle.transform(x)
x_mle
#可以发现，mle自动选择了3个特征
pca_mle.explained_variance_ratio_.sum()
1
2
3
4
5
6

2.2.3 按信息量占比选超参数

输入$[0,1]$之间的浮点数，并且让参数svd_solver==‘full’，表示希望降维后的总解释性方差占比大于n_components指定的百分比，即希望保留百分之多少的信息量。（如希望保留97%的信息量，就可以输入n_components=0.97）

pca_f = PCA(n_components = 0.97,svd_solver = 'full')
pca_f = pca_f.fit(x)
x_f = pca_f.transform(x)
pca_f.explained_variance_ratio_
#array([0.92461872, 0.05306648])
1
2
3
4
5

2.3 PCA中的SVD

2.3.1 PCA中的SVD哪里来？

SVD可以跳过未知的数学变换，不计算协方差矩阵，直接找出一个新特征向量组成的n维空间，而这个n维空间就是奇异值分解后的右矩阵$V^T$。
右奇异矩阵$V^T$有这一性质：$X_{dr}=X\ast V[:k]$
k就是n_components，是降维后希望得到的维度。若$X$为(m,n)的特征矩阵，$V^T$就是结构为(n,n)的矩阵，取这个矩阵的前k行（进行切片），即将$V$转换为结构为(k,n)的矩阵。而$V_{(k,n)}$与原特征矩阵$X$相乘，即可得到降维后的特征矩阵X_dr。也就是说，奇异值分解可以不计算协方差矩阵等等结构复杂计算冗长的矩阵，就直接求出新特征空间和降维后的特征矩阵。
简而言之，SVD在矩阵分解中的过程比PCA简单快速，虽然两个算法都走一样的分解流程，但SVD可以直接算出$V$，但SVD的信息量衡量指标比较复杂，要理解“奇异值”远不如理解“方差”来得容易，因此，sklearn将降维流程拆成了两个部分：一部分是计算特征空间$V$，由奇异值分解完成，另一部分是映射数据和求解新特征矩阵，由主成分分析完成，实现了用SVD的性质减少计算量，却让信息量的评估指标是方差。

通过SVD和PCA的合作，sklearn实现了一种计算更快更简单，但效果却更好的的“合作降维”。很多人理解SVD，是将SVD当做PCA的一种求解方法，其实指的就是在矩阵分解时不使用PCA本身的特征值分解，而使用奇异值分解来减少计算量。这种方法确实存在，但在sklearn中，矩阵$\sum$和$U$虽然会被计算出来，但完全不会被用到，也无法调取查看或者使用，因此可以认为$\sum$和$U$在fit过后就被遗弃了。奇异值分解追求的仅仅是$V$，只要计算出$V$，就可以得出降维后的特征矩阵。在transform过程之后，fit中奇异值分解的结果除了$V(k,n)$以外，都会被舍弃，而$V(k,n)$会被保存在属性components_当中，可以调用查看。

PCA(2).fit(x).components_
PCA(2).fit(x).components_.shape#(k,n)
1
2

2.3.2 重要参数svd_solver与random_state

参数svd_solver是在降维过程中，用来控制矩阵分解的一些细节的参数，有四种模式可选：“auto”，“full”，“arpack”，“randomized”，默认为“auto”。

• “auto”：基于X.shape和n_components的默认策略来选择分解器：如果输入数据的尺寸大于500*500且要提取的特征数小于数据最小维度min(X,shape)的80%，就启用效率更高的“randomized”方法。否则，精确完成的SVD将被计算，截断将会在矩阵被分解完成后有选择的发生
• “full”：从scipy.linalg.svd中调用标准的LAPACK分解器来生成精确完整的SVD，适合数据量比较适中，计算时间充足的情况，生成的精确完整的SVD的结构为：$U_{(m,m)},{\sum }_{(m,n)},V_{(n,n)}^T$
• “arpack”：从scipy.sparse.linalg.svds调用ARPACK分解器来运行截断奇异值分解（SVD truncated），分解时就将 特征数量降到n_components中输入的数值k，可以加快运算速度，适合特征矩阵很大的时候，但一般用于特征矩阵为稀疏矩阵的情况，此过程包含一定的随机性。截断后的SVD分解出的结构为：$U_{(m,k)},{\sum }_{(k,k)},V_{(n,n)}^T$
• “randomized”：通过Halko等人(Finding structure with randomness: Probabilistic algorithms for constructiong approximate matrix decompositions, 2011)的随机方法进行随机SVD。在“full”方法中，分解器会根据原始数据和输入的n_components值去计算和寻找符合需求的新特征向量，但是在“randomized”方法中，分解器会先生成多个随机向量，然后一一去检测这些随机向量中是否有任何一个符合分解需求，如果符合，就保留这个随机向量，并基于这个随机向量来构建后续的向量空间。这个方法已经被Halko等人证明，比“full”模式下计算快很多，并且能够保证模型运行效果。适合特征矩阵巨大，计算量庞大的情况。
  而参数random_state在参数svd_solver的值为“arpack”或者“randomized”的时候生效，可以控制这两种SVD模式中的随机模式。通常就选用“auto”，不必对这个参数纠结太多。

2.3.3 重要属性components_

$V(k,n)$是新特征空间，是要将原始数据进行映射的那些新特征向量组成的矩阵，用它来计算新的特征矩阵。$V(k,n)$保存在components_这个属性中。
PCA与特征选择的区别：特征选择后特征矩阵是可解读的，而PCA降维后的特征矩阵是不可解读的。PCA是将已存在的特征进行压缩，降维完毕后的特征不是原本的特征矩阵中的任何一个特征，而是通过某些方式组合起来的新特征。通常而言，在新的特征矩阵生成之前，无法知晓PCA都建立了怎样的新特征向量，新特征矩阵生成之后也不具有可读性，无法判断新特征矩阵的特征是从原数据的什么特征组合而来，新特征虽然带有原始数据的信息，却已经不是原数据上代表着的含义了。
但其实，在矩阵分解时，PCA是有目标的：在原有特征的基础上，找出能够让信息尽量聚集的新特征向量。在sklearn使用的PCA和SVD联合的降维方法中，这些新特征向量组成的新特征空间其实就是$V(k,n)$。当$V(k,n)$是数字的时候，无法判断$V(k,n)$和原有特征之间有怎样的联系。但是如果原特征矩阵是图像，$V(k,n)$这个空间矩阵也可以别可视化，通过图像比较，就可以发现新特征空间究竟从原数据中提取了什么信息。

#导入需要的库和模块
from sklearn.datasets import fetch_lfw_people
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
import numpy as np
#实例化数据集，探索数据
faces = fetch_lfw_people(min_faces_per_person = 60)
faces.images.shape#数据的维度
#1384是矩阵中图像的个数#62是每个图像的特征矩阵的行#47是每个图像的特征矩阵的列
#结果：(1348, 62, 47)
faces.data.shape#行是样本#列是样本对应的所有特征
#结果：(1348, 2914)
X = faces.data
#看看图像什么样，将原特征矩阵进行可视化
#数据本身是图像，和数据本身只是数字，使用的可视化方法不同
#创建画布和子图对象
fig, axes = plt.subplots(4,5
                        ,figsize = (8,4)
                        ,subplot_kw = {"xticks":[],"yticks":[]}#不要显示各自的坐标轴
                        )
#fig
#axes
axes[0][0].imshow(faces.images[0,:,:])
#不难发现，axes中的一个对象对应fig中的一个空格
#在每一个子图对象中填充图像（共20张图），因此需要写一个在子图对象中遍历的循环
axes.shape
#二维结构可以有两种循环方式，一种是使用索引，循环一次同时生成一列上的三个图，另一种是把数据拉成一维，循环一次只生成一个图
#究竟使用哪种循环方式，要看画得图的信息，储存在一个怎样的结构中
#使用  子图对象.imshow  来将图像填充到空白画布上
#而imshow要求的数据格式必须是一个(m,n)格式的矩阵，即每个数据都是一张单独的图
#因此需要遍历的是faces.images，其结构是(1348, 62, 47)
#要从一个数据集中取出20个图，明显是一次性的循环切片[i,:,:]更加便利
#因此要把axes的结构拉成一维来循环
[*axes.flat]
axes.flat#图像的集合
#填充图像
for i, ax in enumerate(axes.flat):
    ax.imshow(faces.images[i,:,:]
             ,cmap = "gray"#选择色彩的模式
             )
fig, axes = plt.subplots(4,5
                        ,figsize = (8,4)
                        ,subplot_kw = {"xticks":[],"yticks":[]}#不要显示各自的坐标轴
                        )
for i, ax in enumerate(axes.flat):
    ax.imshow(faces.images[i,:,:]
             ,cmap = "gray"#选择色彩的模式
             )
#色彩模式可以选择，参考：https://matplotlib.org/tutorials/colors/colormaps.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49

#建模降维，提取新特征空间矩阵
#原本有2900维，现在来降到150维
pca = PCA(150).fit(X)
v = pca.components_
v.shape
#结果：(150, 2914)
#将新特征空间矩阵可视化
fig, axes = plt.subplots(4,5,figsize = (8,4),subplot_kw = {'xticks':[],'yticks':[]})
for i, ax in enumerate(axes.flat):
    ax.imshow(v[i,:].reshape(62,47),cmap = 'gray')
1
2
3
4
5
6
7
8
9
10

比起降维前的数据，新特征空间可视化后的人脸非常模糊，这是因为原始数据还没有被映射到特征空间中。

2.4 重要接口inverse_transform

inverse_transform可以将归一化、标准化、甚至做过哑变量的特征矩阵还原回原始数据中的特征矩阵。是否可以认为，任何有inverse_transform接口的过程都是可逆的。PCA也是如此，在sklearn中，通过让原特征矩阵$X$右乘新特征空间矩阵$V_{(k,n)}$来生成新空间特征矩阵X_dr。那理论上来说，只要让新特征矩阵X_dr右乘$V_{(k,n)}$的逆矩阵$V_{(k,n)}^{-1}$，就可以将新特征矩阵X_dr还原为原本的特征矩阵。

2.4.1 迷你案例：用人脸识别看PCA降维后的信息保存量

人脸识别是最有效探索inverse_transform功能的数据。先调用一组人脸数据$X_{(m,n)}$，对人脸图像进行绘制，然后对人脸数据进行降维得到X_dr，之后再使用inverse_transform(X_dr)返回X_inverse，并对这个新矩阵中的人脸图像也进行绘制。如果PCA的降维过程是可逆的，应该期待X和X_inverse返回一模一样的图像，即携带一模一样的信息。

#导入需要的库和模块
from sklearn.datasets import fetch_lfw_people
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
import numpy as np
#导入数据，探索数据
#实例化数据集，探索数据
faces = fetch_lfw_people(min_faces_per_person = 60)
faces.images.shape
#结果：(1348, 62, 47)
faces.data.shape
#结果：(1348, 2914)
X = faces.data
#建模降维，获取降维后特征矩阵X_dr
pca = PCA(150)#实例化
x_dr = pca.fit_transform(X)#拟合+提取结果
x_dr.shape
#结果：(1348, 150)
#将降维后矩阵用inverse_transform返回原空间
x_inverse = pca.inverse_transform(x_dr)
x_inverse.shape
#结果：(1348, 2914)
#将特征矩阵X和x_inverse可视化
fig, ax = plt.subplots(2,10,figsize = (10,2.5)
                        ,subplot_kw = {'xticks':[],'yticks':[]}
                        )
#对子图对象进行遍历的循环，来将图像填入子图中
#由于ax中是2行10列，第一行是原数据，第二行是inverse_transform后返回的数据
#所以需要同时循环两份数据，即一次循环画一列上的两张图，而不是把ax拉平
for i in range(10):
    ax[0,i].imshow(faces.images[i,:,:],cmap = 'binary_r')
    ax[1,i].imshow(x_inverse[i].reshape(62,47),cmap = 'binary_r')
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32

可以明显看出，两组数据可视化后，由降维后再通过inverse_transform转换回原维度的数据画出的图像和原始数据画出图像大致相似，但原数据的图像更加清晰。这说明，inverse_transform并没有实现数据的完全逆转。这是因为，在降维的时候，部分信息被舍弃了，X_dr中往往不会包含原数据100%的信息，所以在逆转的时候，即便维度升高，原数据中已经被舍弃的信息也不可能再回来了。所以，降维不是完全可逆的。
inverse_transform的功能，是基于X_dr中的数据进行升维，将数据重新映射到原数据所在的特征空间中，而并非恢复所有原有的功能。但同也可以看出，降维到300以后的数据，的确保留了原数据的大部分信息，所以图像看起来，才会和原数据高度相似，只是稍稍模糊罢了。

2.4.2 迷你案例：用PCA做噪音过滤

降维的目的之一就是希望抛弃掉对模型带来负面影响的特征，带有效信息的特征的方差应该是远远大于噪音的。所以相比噪音，有效的特征所带来的信息应该不会在PCA过程中被大量抛弃。inverse_transform功能能够在不恢复原始数据的情况下，将降维后的数据返回到原本的高维空间，即是说能够实现“保证维度，但去掉方差很小的特征所带的信息”。利用inverse_transform这个性质，能够实现噪音过滤。

#导入所需要的库和模块
from sklearn.datasets import load_digits
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
import numpy as np
#导入数据，探索数据
digits = load_digits()
digits.data.shape
#结果：(1797, 64)
#定义画图函数
def plot_digits(data):
    fig, axes = plt.subplots(4,20,figsize = (10,4)
                            ,subplot_kw = {'xticks':[],'yticks':[]}
                            )
    for i,ax in enumerate(axes.flat):
        ax.imshow(data[i].reshape(8,8),cmap = 'binary')
plot_digits(digits.data)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

#为数据加上噪音
rng = np.random.RandomState(42)
#在指定的数据集中，随机抽取服从正态分布的数据
#两个参数，分别是指定的数据集，和抽取出来的正态分布的方差
noisy = rng.normal(digits.data,2)
plot_digits(noisy)
1
2
3
4
5
6

#降维
pca = PCA(0.5).fit(noisy)#保留原数据50%的信息
x_dr = pca.transform(noisy)
x_dr.shape
#结果：(1797, 6)
#逆转降维结果，实现降噪
without_noise = pca.inverse_transform(x_dr)
plot_digits(without_noise)
1
2
3
4
5
6
7
8

2.5 重要接口、参数和属性总结

重要参数：n_components、svd_solver、random_state
三个重要属性：components_、explained_variance_、explained_variance_ratio_
重要接口：inverse_transform、fit_transform

3 案例：PCA对手写数字数据集的降维

使用kaggle手写数字的数据集

#导入需要的模块和库
from sklearn.decomposition import PCA
from sklearn.ensemble import RandomForestClassifier as RFC
from sklearn.model_selection import cross_val_score
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
#导入数据，探索数据
data = pd.read_csv(r"C:\Users\86188\Desktop\coding\Digit Recognizer_官方数据集\train.csv")
x = data.iloc[:,1:]
y = data.iloc[:,0]
#画累计方差贡献率曲线，找最佳降维后维度的范围
pca_line = PCA().fit(x)
plt.figure(figsize = [20,5])
plt.plot(np.cumsum(pca_line.explained_variance_ratio_))
plt.xlabel("number of components after dimension reduction")
plt.ylabel("cumulative explained variance")
plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18

#降维后维度的学习曲线，继续缩小最佳维度的范围
score=[]
for i in range(1,101,10):
    x_dr = PCA(i).fit_transform(x)
    once = cross_val_score(RFC(n_estimators = 10,random_state = 0)
                          ,x_dr,y,cv = 5).mean()
    score.append(once)
plt.figure(figsize = [20,5])
plt.plot(range(1,101,10),score)
plt.show()
1
2
3
4
5
6
7
8
9
10

#细化学习曲线，找出降维后的最佳维度
score = []
for i in range(10,25):
    x_dr = PCA(i).fit_transform(x)
    once = cross_val_score(RFC(n_estimators = 10,random_state = 0),x_dr,y,cv=5)
    score.append(once)
plt.figure(figsize = [20,5])
plt.plot(range(10,25),score)
plt.show()
1
2
3
4
5
6
7
8
9

#导入找出的最佳维度进行降维，查看模型效果
x_dr = PCA(21).fit_transform(x)
cross_val_score(RFC(n_estimators = 10,random_state = 0),x_dr,y,cv=5).mean()
#结果：0.9177857142857142
cross_val_score(RFC(n_estimators = 100,random_state = 0),x_dr,y,cv=5).mean()
#结果：0.9441190476190476
1
2
3
4
5
6

模型的效果达到了94.41%的水平，但还是没有使用嵌入法特征选择后的96%高，还能怎样提高模型的效果呢？之前因为计算量太大，一直使用随机森林，但事实上，KNN的效果会比随机森林更好，KNN在未调参的状况下已经达到96%的准确率，而随机森林在未调参前只能达到93%，这是模型本身的限制带来的，这个数据使用KNN的效果就是会更好。

#特征数量已经不足原来的3%，换成KNN模型
from sklearn.neighbors import KNeighborsClassifier as KNN
cross_val_score(KNN(),x_dr,y,cv=5).mean()
#结果：0.9676428571428571
#KNN的k值学习曲线
score=[]
for i in range(10):
    x_dr = PCA(21).fit_transform(x)
    once = cross_val_score(KNN(i+1),x_dr,y,cv=5).mean()
    score.append(once)
plt.figure(figsize = [20,5])
plt.plot(range(10),score)
plt.show()
1
2
3
4
5
6
7
8
9
10
11
12
13

#定下超参数后，模型的效果如何，模型的运行时间如何
cross_val_score(KNN(3),x_dr,y,cv=5).mean()
#结果：0.9682619047619048
%%timeit
cross_val_score(KNN(3),x_dr,y,cv=5).mean()
#结果：19.4 s ± 123 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
1
2
3
4
5
6

参考资料：https://www.bilibili.com/video/BV1h34y1m7TU?p=53