《菜菜的机器学习 sklearn 课堂》降维算法 PCA 和 SVD
PCA 作为矩阵分解算法的核心算法,其实没有太多参数,但每个参数的意义和运用都很难,因为几乎每个参数都涉及到高深的数学原理。为了参数的运用和意义变得明朗,我们来看一组简单的二维数据的降维。
我们现在有一组简单的数据,有特征 x1 和 x2,三个样本数据的坐标点分别为(1,1),(2,2),(3,3)。我们可以让 x1 和 x2 分别作为两个特征向量,很轻松地用一个二维平面来描述这组数据。这组数据现在每个特征的均值都为 2,方差则等于:
x 1 _ v a r = x 2 _ v a r = ( 1 ? 2 ) 2 + ( 2 ? 2 ) 2 + ( 3 ? 2 ) 2 2 = 1 x_1\_{var} = x_2\_var = \frac {(1-2)^2 + (2-2)^2 + (3-2)^2} 2 = 1 x1?_var=x2?_var=2(1?2)2+(2?2)2+(3?2)2?=1
每个特征的数据一模一样,因此方差也都为 1,数据的方差总和是 2。
现在我们的目标是:只用一个特征向量来描述这组数据,即将二维数据降为一维数据,并且尽可能地保留信息量,即让数据的总方差尽量靠近 2。于是,我们将原本的直角坐标系逆时针旋转 45°,形成了新的特征向量 x1*和 x2*组成的新平面,在这个新平面中,三个样本数据的坐标点可以表示为: ( 2 , 0 ) (\sqrt2,0) (2 ?,0), ( 2 2 , 0 ) (2\sqrt2,0) (22 ?,0), ( 3 2 , 0 ) (3\sqrt2,0) (32 ?,0)。可以注意到,x2*上的数值此时都变成了 0,因此 x2*明显不带有任何有效信息了(此时 x2*的方差也为 0 了)。此时,x1*特征上的数据均值是 2 2 2\sqrt2 22 ? ,而方差则可表示成:
x 2 ? _ v a r = ( 2 ? 2 2 ) 2 + ( 2 2 ? 2 2 ) 2 + ( 3 2 ? 2 2 ) 2 2 = 2 x_2^*\_var = \frac {(\sqrt2 - 2\sqrt2)^2 + (2\sqrt2 -2\sqrt2)^2 + (3\sqrt2 - 2\sqrt2)^2} 2 = 2 x2??_var=2(2 ??22 ?)2+(22 ??22 ?)2+(32 ??22 ?)2?=2
x1* 上的数据均值为 0,方差也为 0。
此时,我们根据信息含量的排序,取信息含量最大的一个特征,因为我们想要的是一维数据。所以我们可以将 x2*删除,同时也删除图中的 x2*特征向量,剩下的 x1*就代表了曾经需要两个特征来代表的三个样本点。通过旋转原有特征向量组成的坐标轴来找到新特征向量和新坐标平面,我们将三个样本点的信息压缩到了一条直线上,实现了二维变一维,并且尽量保留原始数据的信息。一个成功的降维就实现了。
在这个降维过程中,有几个重要的步骤,同时将 2 维延伸到 n 维:
在步骤 3 当中,我们用来找出 n 个新特征向量,让数据能够被压缩到少数特征上并且总信息量不损失太多的技术就是矩阵分解。
PCA 和 SVD 是两种不同的降维算法,但他们都遵从上面的过程来实现降维,只是两种算法中矩阵分解的方法不同,信息量的衡量指标不同罢了。
PCA 使用方差作为信息量的衡量指标,并且特征值分解来找出空间 V。
降维时,它会通过一系列数学的神秘操作(比如说,产生协方差矩阵 1 n X X T \frac 1 n XX^T n1?XXT)将特征矩阵 X 分解为以下三个矩阵,其中 Q Q Q 和 Q ? 1 Q^{-1} Q?1 是辅助的矩阵, Σ \Sigma Σ是一个对角矩阵(即除了对角线上有值,其他位置都是 0 的矩阵),其对角线上的元素就是方差。
降维完成之后,PCA 找到的每个新特征向量就叫做“主成分”,而被丢弃的特征向量被认为信息量很少,这些信息很可能就是噪音。
X → 数 学 神 秘 的 宇 宙 → Q Σ Q ? 1 X \rightarrow 数学神秘的宇宙 \rightarrow Q \Sigma Q^{-1} X→数学神秘的宇宙→QΣQ?1
SVD 使用奇异值分解来找出空间 V,其中 Σ \Sigma Σ也是一个对角矩阵,不过它对角线上的元素是奇异值,这也是 SVD 中用来衡量特征上的信息量的指标。
U 和 V T V^{T} VT 分别是左奇异矩阵和右奇异矩阵,也都是辅助矩阵。
X → 另 一 个 数 学 神 秘 宇 宙 → U Σ V T X \rightarrow 另一个数学神秘宇宙 \rightarrow U\Sigma V^T X→另一个数学神秘宇宙→UΣVT
在数学原理中,无论是 PCA 和 SVD 都需要遍历所有的特征和样本来计算信息量指标。并且在矩阵分解的过程之中,会产生比原来的特征矩阵更大的矩阵,比如原数据的结构是(m,n),在矩阵分解中为了找出最佳新特征空间 V,可能需要产生(n,n),(m,m)大小的矩阵,还需要产生协方差矩阵去计算更多的信息。而现在无论是 Python 还是 R,或者其他的任何语言,在大型矩阵运算上都不是特别擅长,无论代码如何简化,我们不可避免地要等待计算机去完成这个非常庞大的数学计算过程。因此,降维算法的计算量很大,运行比较缓慢,但无论如何,它们的功能无可替代,它们依然是机器学习领域的宠儿。
思考:PCA 和特征选择技术都是特征工程的一部分,它们有什么不同?
特征工程中有三种方式:特征提取,特征创造和特征选择。
特征选择是从已存在的特征中选取携带信息最多的,选完之后的特征依然具有可解释性,我们依然知道这个特征在原数据的哪个位置,代表着原数据上的什么含义。
而 PCA,是将已存在的特征进行压缩,降维完毕后的特征不是原本的特征矩阵中的任何一个特征,而是通过某些方式组合起来的新特征。通常来说,在新的特征矩阵生成之前,我们无法知晓 PCA 都建立了怎样的新特征向量,新特征矩阵生成之后也不具有可读性,我们无法判断新特征矩阵的特征是从原数据中的什么特征组合而来,新特征虽然带有原始数据的信息,却已经不是原数据上代表着的含义了。以 PCA 为代表的降维算法因此是特征创造(feature creation,或 feature construction)的一种。
可以想见,PCA 一般不适用于探索特征和标签之间的关系的模型(如线性回归),因为无法解释的新特征和标签之间的关系不具有意义。在线性回归模型中,我们使用特征选择。
========================================================================================
n_components 是我们降维后需要的维度,即降维后需要保留的特征数量,一般输入[0, min(X.shape)]范围中的整数。
如果我们希望可视化一组数据来观察数据分布,我们往往将数据降到三维以下,很多时候是二维,即 n_components 的取值为 2。
…
================================================================================
1.调用库和模块
import matplotlib.pyplot as plt
from sklearn.datasets import load_iris # 鸢尾花数据集
from sklearn.decomposition import PCA
2.提取数据集
iris = load_iris()
y = iris.target # 分类
iris.target_names # 分类的名字
X = iris.data # 数据
作为数组,X 是几维?
X.shape #作为数组,它有 150 行 4 列,是 2 维
#作为数据表或特征矩阵,X 是几维?
import pandas as pd
pd.DataFrame(X) #作为特征矩阵,它有 4 个特征,所以是 4 维
3.建模
#调用 PCA
pca = PCA(n_components=2) #实例化
pca = pca.fit(X) #拟合模型
X_dr = pca.transform(X) #获取降维后的新特征矩阵
fit_transform 可以代替上面三行代码
X_dr = PCA(2).fit_transform(X)
X_dr.shape # 150 行,2 列
pd.DataFrame(X_dr) # 降维后特征矩阵有 2 个特征,所以是二维
4.可视化
要将三种鸢尾花的数据分布显示在二维平面坐标系中,
对应的两个坐标(两个特征向量)应该是三种鸢尾花降维后的 x1 和 x2,
怎样才能取出三种鸢尾花下不同的 x1 和 x2 呢?
y # 结果显示为 3 分类
"""
array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2])
"""
X_dr[y == 0, 0] # 布尔索引,表示取分类为 0 的第一列
X_dr[y == 1, 0] # 取分类为 1 的第一列
要展示三种分类的分布
,需要对三种鸢尾花分别绘图
可以写成三行代码,也可以写成 for 循环
iris.target_names #分类的名字
画法 1
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()
画法 2
colors = ['red', 'black', 'orange']
plt.figure()
for i in [0, 1, 2]:
plt.scatter(X_dr[y == i, 0],
X_dr[y == i, 1],
alpha=.7, #点的透明程度 70%
c=colors[i],
label=iris.target_names[i])
plt.legend()
plt.title('PCA of IRIS dataset')
plt.show()
鸢尾花的分布被展现在我们眼前了,明显这是一个分簇的分布,并且每个簇之间的分布相对比较明显,也许 versicolor 和 virginia 这两种花之间会有一些分类错误,但 setosa 肯定不会被分错。
这样的数据很容易分类,可以遇见:KNN,随机森林,神经网络,朴素贝叶斯,Adaboost 这些分类器在鸢尾花数据集上,未调整的时候都可以有 95%上下的准确率。
5.探索降维后的数据
#属性 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()
#降维后保留了 97.768%的信息
0.977685206318795
6.选择最好的 n_components:累积可解释方差贡献率曲线
当参数 n_components 中不填写任何值,则默认返回min(X.shape)
个特征,一般样本量都会大于特征数目,所以什么都不填就相当于转换了新特征空间,但没有减少特征的个数。
一般来说,不会使用这种输入方式。但我们却可以使用这种输入方式来画出累计可解释方差贡献率曲线,以此选择最好的 n_components 的整数取值。
累积可解释方差贡献率曲线:
横坐标:降维后保留的特征个数
纵坐标:降维后新特征矩阵捕捉到的可解释方差贡献率
累积可解释方差贡献率曲线能够帮助我们决定 n_components 最好的取值。
import numpy as np
pca_line = PCA().fit(X)
pca_line.transform(X).shape # (150, 4) 特征值个数没有变
降维后新特征矩阵捕捉到的可解释方差贡献率
pca_line.explained_variance_ratio_
"""
array([0.92461872, 0.05306648, 0.01710261, 0.00521218])
"""
[1,2,3,4]是设置横坐标,默认是 0,1,2,3
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 ratio")
plt.show()
=====================================================================
评论