主成分分析(PCA)
主成分分析(Principal Component Analysis)是目前为止最流行的降维算法。首先它找到接近数据集分布的超平面,然后将所有的数据都投影到这个超平面上。
在将训练集投影到较低维超平面之前,您首先需要选择正确的超平面。例如图 8-7 左侧是一个简单的二维数据集,以及三个不同的轴(即一维超平面)。图右边是将数据集投影到每个轴上的结果。正如你所看到的,投影到实线上保留了最大方差,而在点线上的投影只保留了非常小的方差,投影到虚线上保留的方差则处于上述两者之间。
图 8-7 选择投射到哪一个子空间
选择保持最大方差的轴看起来是合理的,因为它很可能比其他投影损失更少的信息。证明这种选择的另一种方法是,选择这个轴使得将原始数据集投影到该轴上的均方距离最小。这是就 PCA 背后的思想,相当简单。
PCA 寻找训练集中可获得最大方差的轴。在图 8-7 中,它是一条实线。它还发现了一个与第一个轴正交的第二个轴,选择它可以获得最大的残差。在这个 2D 例子中,没有选择:就只有这条点线。但如果在一个更高维的数据集中,PCA 也可以找到与前两个轴正交的第三个轴,以及与数据集中维数相同的第四个轴,第五个轴等。定义第i
个轴的单位矢量被称为第i
个主成分(PC)。在图 8-7 中,第一个 PC 是c1
,第二个 PC 是c2
。在图 8-2 中,前两个 PC 用平面中的正交箭头表示,第三个 PC 与上述 PC 形成的平面正交(指向上或下)。
概述: 主成分的方向不稳定:如果您稍微打乱一下训练集并再次运行 PCA,则某些新 PC 可能会指向与原始 PC 方向相反。但是,它们通常仍位于同一轴线上。在某些情况下,一对 PC 甚至可能会旋转或交换,但它们定义的平面通常保持不变。
那么如何找到训练集的主成分呢?幸运的是,有一种称为奇异值分解(SVD)的标准矩阵分解技术,可以将训练集矩阵X
分解为三个矩阵U·Σ·V^T
的点积,其中V^T
包含我们想要的所有主成分,如公式 8-1 所示。公式 8-1 主成分矩阵
下面的 Python 代码使用了 Numpy 提供的svd()
函数获得训练集的所有主成分,然后提取前两个 PC:
X_centered=X-X.mean(axis=0)U,s,V=np.linalg.svd(X_centered)c1=V.T[:,0]c2=V.T[:,1]
警告:PCA 假定数据集以原点为中心。正如我们将看到的,Scikit-Learn 的
PCA
类负责为您的数据集中心化处理。但是,如果您自己实现 PCA(如前面的示例所示),或者如果您使用其他库,不要忘记首先要先对数据做中心化处理。
d
维空间 一旦确定了所有的主成分,你就可以通过将数据集投影到由前d
个主成分构成的超平面上,从而将数据集的维数降至d
维。选择这个超平面可以确保投影将保留尽可能多的方差。例如,在图 8-2 中,3D 数据集被投影到由前两个主成分定义的 2D 平面,保留了大部分数据集的方差。因此,2D 投影看起来非常像原始 3D 数据集。
为了将训练集投影到超平面上,可以简单地通过计算训练集矩阵X
和Wd
的点积,Wd
定义为包含前d
个主成分的矩阵(即由V^T
的前d
列组成的矩阵),如公式 8-2 所示。公式 8-2 将训练集投影到d
维空间
下面的 Python 代码将训练集投影到由前两个主成分定义的超平面上:
W2=V.T[:,:2]X2D=X_centered.dot(W2)
好了你已经知道这个东西了!你现在已经知道如何给任何一个数据集降维而又能尽可能的保留原数据集的方差了。
Scikit-Learn 的 PCA 类使用 SVD 分解来实现,就像我们之前做的那样。以下代码应用 PCA 将数据集的维度降至两维(请注意,它会自动处理数据的中心化):
from sklearn.decomposition import PCApca=PCA(n_components=2)X2D=pca.fit_transform(X)
将 PCA 转化器应用于数据集后,可以使用components_
访问每一个主成分(注意,它返回以 PC 作为水平向量的矩阵,因此,如果我们想要获得第一个主成分则可以写成pca.components_.T[:,0]
)。
另一个非常有用的信息是每个主成分的方差解释率,可通过explained_variance_ratio_
变量获得。它表示位于每个主成分轴上的数据集方差的比例。例如,让我们看一下图 8-2 中表示的三维数据集前两个分量的方差解释率:
>>> print(pca.explained_variance_ratio_)array([0.84248607, 0.14631839])
这表明,84.2% 的数据集方差位于第一轴,14.6% 的方差位于第二轴。第三轴的这一比例不到1.2%,因此可以认为它可能没有包含什么信息。
通常我们倾向于选择加起来到方差解释率能够达到足够占比(例如 95%)的维度的数量,而不是任意选择要降低到的维度数量。当然,除非您正在为数据可视化而降低维度 -- 在这种情况下,您通常希望将维度降低到 2 或 3。
下面的代码在不降维的情况下进行 PCA,然后计算出保留训练集方差 95% 所需的最小维数:
pca=PCA()pac.fit(X)cumsum=np.cumsum(pca.explained_variance_ratio_)d=np.argmax(cumsum>=0.95)+1
你可以设置n_components = d
并再次运行 PCA。但是,有一个更好的选择:不指定你想要保留的主成分个数,而是将n_components
设置为 0.0 到 1.0 之间的浮点数,表明您希望保留的方差比率:
pca=PCA(n_components=0.95)X_reduced=pca.fit_transform(X)
另一种选择是画出方差解释率关于维数的函数(简单地绘制cumsum
;参见图 8-8)。曲线中通常会有一个肘部,方差解释率停止快速增长。您可以将其视为数据集的真正的维度。在这种情况下,您可以看到将维度降低到大约100个维度不会失去太多的可解释方差。
图 8-8 可解释方差关于维数的函数
显然,在降维之后,训练集占用的空间要少得多。例如,尝试将 PCA 应用于 MNIST 数据集,同时保留 95% 的方差。你应该发现每个实例只有 150 多个特征,而不是原来的 784 个特征。因此,尽管大部分方差都保留下来,但数据集现在还不到其原始大小的 20%!这是一个合理的压缩比率,您可以看到这可以如何极大地加快分类算法(如 SVM 分类器)的速度。
通过应用 PCA 投影的逆变换,也可以将缩小的数据集解压缩回 784 维。当然这并不会返回给你最原始的数据,因为投影丢失了一些信息(在5%的方差内),但它可能非常接近原始数据。原始数据和重构数据之间的均方距离(压缩然后解压缩)被称为重构误差(reconstruction error)。例如,下面的代码将 MNIST 数据集压缩到 154 维,然后使用inverse_transform()
方法将其解压缩回 784 维。图 8-9 显示了原始训练集(左侧)的几位数字在压缩并解压缩后(右侧)的对应数字。您可以看到有轻微的图像质量降低,但数字仍然大部分完好无损。
pca=PCA(n_components=154)X_mnist_reduced=pca.fit_transform(X_mnist)X_mnist_recovered=pca.inverse_transform(X_mnist_reduced)
图 8-9 MNIST 保留 95 方差的压缩
逆变换的公式如公式 8-3 所示公式 8-3 PCA逆变换,回退到原来的数据维度
先前 PCA 实现的一个问题是它需要在内存中处理整个训练集以便 SVD 算法运行。幸运的是,我们已经开发了增量 PCA(IPCA)算法:您可以将训练集分批,并一次只对一个批量使用 IPCA 算法。这对大型训练集非常有用,并且可以在线应用 PCA(即在新实例到达时即时运行)。
下面的代码将 MNIST 数据集分成 100 个小批量(使用 NumPy 的array_split()
函数),并将它们提供给 Scikit-Learn 的IncrementalPCA
类,以将 MNIST 数据集的维度降低到 154 维(就像以前一样)。请注意,您必须对每个最小批次调用partial_fit()
方法,而不是对整个训练集使用fit()
方法:
from sklearn.decomposition import IncrementalPCAn_batches=100inc_pca=IncrementalPCA(n_components=154)for X_batch in np.array_spplit(X_mnist,n_batches): inc_pca.partial_fit(X_batch)X_mnist_reduced=inc_pca.transform(X_mnist)
或者,您可以使用 NumPy 的memmap
类,它允许您操作存储在磁盘上二进制文件中的大型数组,就好像它完全在内存中;该类仅在需要时加载内存中所需的数据。由于增量 PCA 类在任何时间内仅使用数组的一小部分,因此内存使用量仍受到控制。这可以调用通常的fit()
方法,如下面的代码所示:
X_mm=np.memmap(filename,dtype='float32',mode='readonly',shape=(m,n))batch_size=m//n_batchesinc_pca=IncrementalPCA(n_components=154,batch_size=batch_size)inc_pca.fit(X_mm)
Scikit-Learn 提供了另一种执行 PCA 的选择,称为随机 PCA。这是一种随机算法,可以快速找到前d个主成分的近似值。它的计算复杂度是O(m × d^2) + O(d^3)
,而不是O(m × n^2) + O(n^3)
,所以当d
远小于n
时,它比之前的算法快得多。
rnd_pca=PCA(n_components=154,svd_solver='randomized')X_reduced=rnd_pca.fit_transform(X_mnist)
在第 5 章中,我们讨论了核技巧,一种将实例隐式映射到非常高维空间(称为特征空间)的数学技术,让支持向量机可以应用于非线性分类和回归。回想一下,高维特征空间中的线性决策边界对应于原始空间中的复杂非线性决策边界。
事实证明,同样的技巧可以应用于 PCA,从而可以执行复杂的非线性投影来降低维度。这就是所谓的核 PCA(kPCA)。它通常能够很好地保留投影后的簇,有时甚至可以展开分布近似于扭曲流形的数据集。
例如,下面的代码使用 Scikit-Learn 的KernelPCA
类来执行带有 RBF 核的 kPCA(有关 RBF 核和其他核的更多详细信息,请参阅第 5 章):
from sklearn.decomposition import KernelPCArbf_pca=KernelPCA(n_components=2,kernel='rbf',gamma=0.04)X_reduced=rbf_pca.fit_transform(X)
图 8-10 展示了使用线性核(等同于简单的使用 PCA 类),RBF 核,sigmoid 核(Logistic)将瑞士卷降到 2 维。
图 8-10 使用不同核的 kPCA 将瑞士卷降到 2 维
由于 kPCA 是无监督学习算法,因此没有明显的性能指标可以帮助您选择最佳的核方法和超参数值。但是,降维通常是监督学习任务(例如分类)的准备步骤,因此您可以简单地使用网格搜索来选择可以让该任务达到最佳表现的核方法和超参数。例如,下面的代码创建了一个两步的流水线,首先使用 kPCA 将维度降至两维,然后应用 Logistic 回归进行分类。然后它使用Grid SearchCV
为 kPCA 找到最佳的核和gamma
值,以便在最后获得最佳的分类准确性:
from sklearn.model_selection import GridSearchCV from sklearn.linear_model import LogisticRegression from sklearn.pipeline import Pipelineclf = Pipeline([ ("kpca", KernelPCA(n_components=2)), ("log_reg", LogisticRegression())])param_grid = [{ "kpca__gamma": np.linspace(0.03, 0.05, 10), "kpca__kernel": ["rbf", "sigmoid"] }]grid_search = GridSearchCV(clf, param_grid, cv=3)grid_search.fit(X, y)
你可以通过调用best_params_
变量来查看使模型效果最好的核和超参数:
>>> print(grid_search.best_params_){'kpca__gamma': 0.043333333333333335, 'kpca__kernel': 'rbf'}
另一种完全为非监督的方法,是选择产生最低重建误差的核和超参数。但是,重建并不像线性 PCA 那样容易。这里是原因:图 8-11 显示了原始瑞士卷 3D 数据集(左上角),并且使用 RBF 核应用 kPCA 后生成的二维数据集(右上角)。由于核技巧,这在数学上等同于使用特征映射φ
将训练集映射到无限维特征空间(右下),然后使用线性 PCA 将变换的训练集投影到 2D。请注意,如果我们可以在缩减空间中对给定实例实现反向线性 PCA 步骤,则重构点将位于特征空间中,而不是位于原始空间中(例如,如图中由x
表示的那样)。由于特征空间是无限维的,我们不能找出重建点,因此我们无法计算真实的重建误差。幸运的是,可以在原始空间中找到一个贴近重建点的点。这被称为重建前图像(reconstruction pre-image)。一旦你有这个前图像,你就可以测量其与原始实例的平方距离。然后,您可以选择最小化重建前图像错误的核和超参数。
图 8-11 核 PCA 和重建前图像误差
您可能想知道如何进行这种重建。一种解决方案是训练一个监督回归模型,将预计实例作为训练集,并将原始实例作为训练目标。如果您设置了fit_inverse_transform = True
,Scikit-Learn 将自动执行此操作,代码如下所示:
rbf_pca = KernelPCA(n_components = 2, kernel="rbf", gamma=0.0433,fit_inverse_transform=True)X_reduced = rbf_pca.fit_transform(X)X_preimage = rbf_pca.inverse_transform(X_reduced)
概述:默认条件下,
fit_inverse_transform = False
并且KernelPCA
没有inverse_tranfrom()
方法。这种方法仅仅当fit_inverse_transform = True
的情况下才会创建。
你可以计算重建前图像误差:
>>> from sklearn.metrics import mean_squared_error>>> mean_squared_error(X, X_preimage) 32.786308795766132
现在你可以使用交叉验证的方格搜索来寻找可以最小化重建前图像误差的核方法和超参数。
局部线性嵌入(Locally Linear Embedding)是另一种非常有效的非线性降维(NLDR)方法。这是一种流形学习技术,不依赖于像以前算法那样的投影。简而言之,LLE 首先测量每个训练实例与其最近邻(c.n.)之间的线性关系,然后寻找能最好地保留这些局部关系的训练集的低维表示(稍后会详细介绍) 。这使得它特别擅长展开扭曲的流形,尤其是在没有太多噪音的情况下。
例如,以下代码使用 Scikit-Learn 的LocallyLinearEmbedding
类来展开瑞士卷。得到的二维数据集如图 8-12 所示。正如您所看到的,瑞士卷被完全展开,实例之间的距离保存得很好。但是,距离不能在较大范围内保留的很好:展开的瑞士卷的左侧被挤压,而右侧的部分被拉长。尽管如此,LLE 在对流形建模方面做得非常好。
from sklearn.manifold import LocallyLinearEmbeddinglle=LocallyLinearEmbedding(n_components=2,n_neighbors=10)X_reduced=lle.fit_transform(X)
图 8-12 使用 LLE 展开瑞士卷
这是LLE的工作原理:首先,对于每个训练,该算法识别其最近的k
个邻居(在前面的代码中k = 10
中),然后尝试重构为这些邻居的线性函数。更特殊的,它假设如果 不是的k
个最近邻之一,就找到权重从而使和之间的平方距离尽可能的小。因此,LLE 的第一步是方程 8-4 中描述的约束优化问题,其中W
是包含所有权重的权重矩阵。第二个约束简单地对每个训练实例的权重进行归一化。
公式 8-4 LLE 第一步:对局部关系进行线性建模
在这步之后,权重矩阵(包含权重对训练实例的线形关系进行编码。现在第二步是将训练实例投影到一个d
维空间(d < n
)中去,同时尽可能的保留这些局部关系。如果是在这个d
维空间的图像,那么我们想要和之间的平方距离尽可能的小。这个想法让我们提出了公式8-5中的非限制性优化问题。它看起来与第一步非常相似,但我们要做的不是保持实例固定并找到最佳权重,而是恰相反:保持权重不变,并在低维空间中找到实例图像的最佳位置。请注意,Z
是包含所有的矩阵。
公式 8-5 LLE 第二步:保持关系的同时进行降维
Scikit-Learn 的 LLE 实现具有如下的计算复杂度:查找k
个最近邻为O(m log(m) n log(k))
,优化权重为O(m n k^3)
,建立低维表示为O(d m^2)
。不幸的是,最后一项m^2
使得这个算法在处理大数据集的时候表现较差。
还有很多其他的降维方法,Scikit-Learn 支持其中的好几种。这里是其中最流行的:
**缩放(MDS)在尝试保持实例之间距离的同时降低了维度(参见图 8-13)
Isomap 通过将每个实例连接到最近的邻居来创建图形,然后在尝试保持实例之间的测地距离时降低维度。
t-分布随机邻域嵌入(t-Distributed Stochastic Neighbor Embedding,t-SNE)可以用于降低维度,同时试图保持相似的实例临近并将不相似的实例分开。它主要用于可视化,尤其是用于可视化高维空间中的实例(例如,可以将MNIST图像降维到 2D 可视化)。
线性判别分析(Linear Discriminant Analysis,LDA)实际上是一种分类算法,但在训练过程中,它会学习类之间最有区别的轴,然后使用这些轴来定义用于投影数据的超平面。LDA 的好处是投影会尽可能地保持各个类之间距离,所以在运行另一种分类算法(如 SVM 分类器)之前,LDA 是很好的降维技术。
图 8-13 使用不同的技术将瑞士卷降维至 2D
减少数据集维度的主要动机是什么?主要缺点是什么?
什么是维度爆炸?
一旦对某数据集降维,我们可能恢复它吗?如果可以,怎样做才能恢复?如果不可以,为什么?
PCA 可以用于降低一个高度非线性对数据集吗?
假设你对一个 1000 维的数据集应用 PCA,同时设置方差解释率为 95%,你的最终数据集将会有多少维?
在什么情况下你会使用普通的 PCA,增量 PCA,随机 PCA 和核 PCA?
你该如何评价你的降维算法在你数据集上的表现?
将两个不同的降维算法串联使用有意义吗?
加载 MNIST 数据集(在第 3 章中介绍),并将其分成一个训练集和一个测试集(将前 60,000 个实例用于训练,其余 10,000 个用于测试)。在数据集上训练一个随机森林分类器,并记录了花费多长时间,然后在测试集上评估模型。接下来,使用 PCA 降低数据集的维度,设置方差解释率为 95%。在降维后的数据集上训练一个新的随机森林分类器,并查看需要多长时间。训练速度更快?接下来评估测试集上的分类器:它与以前的分类器比较起来如何?
使用 t-SNE 将 MNIST 数据集缩减到二维,并使用 Matplotlib 绘制结果图。您可以使用 10 种不同颜色的散点图来表示每个图像的目标类别。或者,您可以在每个实例的位置写入彩色数字,甚至可以绘制数字图像本身的降维版本(如果绘制所有数字,则可视化可能会过于混乱,因此您应该绘制随机样本或只在周围没有其他实例被绘制的情况下绘制)。你将会得到一个分隔良好的的可视化数字集群。尝试使用其他降维算法,如 PCA,LLE 或 MDS,并比较可视化结果。