我们要讨论的第二种机器学习算法是无监督学习算法。无监督学习包括没有已知输出、没有老师指导学习算法的各种机器学习。在无监督学习中,学习算法只有输入数据,并需要从这些数据中提取知识。
3.1 无监督学习的类型
本章将研究两种类型的无监督学习:数据集变换与聚类。
数据集的无监督变换(unsupervised transformation)是创建数据新的表示的算法,与数据的原始表示相比,新的表示可能更容易被人或其他机器学习算法所理解。无监督变换的一个常见应用是降维(dimensionality reduction),它接受包含许多特征的数据的高维表示,并找到表示该数据的一种新方法,用较少的特征就可以概括其重要特性。降维的一个常见应用是为了可视化将数据降为二维。
无监督变换的另一个应用是找到“构成”数据的各个组成部分。这方面的一个例子就是对文本文档集合进行主题提取。这里的任务是找到每个文档中讨论的未知主题,并学习每个文档中出现了哪些主题。这可以用于追踪社交媒体上的话题讨论,比如选举、枪支管制或流行歌手等话题。
与之相反, 聚类算法(clustering algorithm)将数据划分成不同的组,每组包含相似的物项。思考向社交媒体网站上传照片的例子。为了方便你整理照片,网站可能想要将同一个人的照片分在一组。但网站并不知道每张照片是谁,也不知道你的照片集中出现了多少个人。明智的做法是提取所有的人脸,并将看起来相似的人脸分在一组。但愿这些人脸对应同一个人,这样图片的分组也就完成了。
3.2 无监督学习的挑战
无监督学习的一个主要挑战就是评估算法是否学到了有用的东西。无监督学习算法一般用于不包含任何标签信息的数据,所以我们不知道正确的输出应该是什么。因此很难判断一个模型是否“表现很好”。例如,假设我们的聚类算法已经将所有的侧脸照片和所有的正面照片进行分组。这肯定是人脸照片集合的一种可能的划分方法,但并不是我们想要的那种方法。然而,我们没有办法“告诉”算法我们要的是什么,通常来说,评估无监督算法结果的唯一方法就是人工检查。
因此,如果数据科学家想要更好地理解数据,那么无监督算法通常可用于探索性的目的,而不是作为大型自动化系统的一部分。无监督算法的另一个常见应用是作为监督算法的预处理步骤。学习数据的一种新表示,有时可以提高监督算法的精度,或者可以减少内存占用和时间开销。
在开始学习“真正的”无监督算法之前,我们先简要讨论几种简单又常用的预处理方法。虽然预处理和缩放通常与监督学习算法一起使用,但缩放方法并没有用到与“监督”有关的信息,所以它是无监督的。
3.3 预处理与缩放
上一章我们学到,一些算法(如神经网络和 SVM)对数据缩放非常敏感。因此,通常的做法是对特征进行调节,使数据表示更适合于这些算法。通常来说,这是对数据的一种简单的按特征的缩放和移动。下面的代码(图 3-1)给出了一个简单的例子:
In[2]:
mglearn.plots.plot_scaling()
3.3.1 不同类型的预处理
在图 3-1 中,第一张图显示的是一个模拟的有两个特征的二分类数据集。第一个特征(x轴)位于 10 到 15 之间。第二个特征(y 轴)大约位于 1 到 9 之间。
接下来的 4 张图展示了 4 种数据变换方法,都生成了更加标准的范围。 scikit-learn 中的 StandardScaler 确保每个特征的平均值为 0、方差为 1,使所有特征都位于同一量级。但这种缩放不能保证特征任何特定的最大值和最小值。 RobustScaler 的工作原理与StandardScaler 类似,确保每个特征的统计属性都位于同一范围。但 RobustScaler 使用的是中位数和四分位数 1,而不是平均值和方差。这样 RobustScaler 会忽略与其他点有很大不同的数据点(比如测量误差)。这些与众不同的数据点也叫异常值(outlier),可能会给其他缩放方法造成麻烦。
与之相反, MinMaxScaler 移动数据,使所有特征都刚好位于 0 到 1 之间。对于二维数据集来说,所有的数据都包含在 x 轴 0 到 1 与 y 轴 0 到 1 组成的矩形中。最后, Normalizer 用到一种完全不同的缩放方法。它对每个数据点进行缩放,使得特征向量的欧式长度等于 1。换句话说,它将一个数据点投射到半径为 1 的圆上(对于更高维度的情况,是球面)。这意味着每个数据点的缩放比例都不相同(乘以其长度的倒数)。如果只有数据的方向(或角度)是重要的,而特征向量的长度无关紧要,那么通常会使用这种归一化。
3.3.2 应用数据变换
前面我们已经看到不同类型的变换的作用,下面利用 scikit-learn 来应用这些变换。我们将使用第 2 章见过的 cancer 数据集。通常在应用监督学习算法之前使用预处理方法(比如缩放)。举个例子,比如我们想要将核 SVM(SVC)应用在 cancer 数据集上,并使用MinMaxScaler 来预处理数据。首先加载数据集并将其分为训练集和测试集(我们需要分开的训练集和数据集来对预处理后构建的监督模型进行评估):
In[3]:
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split
cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target,
random_state=1)
print(X_train.shape)
print(X_test.shape)
Out[3]:
(426, 30)
(143, 30)
提醒一下,这个数据集包含 569 个数据点,每个数据点由 30 个测量值表示。我们将数据集分成包含 426 个样本的训练集与包含 143 个样本的测试集。
与之前构建的监督模型一样,我们首先导入实现预处理的类,然后将其实例化:
In[4]:
from sklearn.preprocessing import MinMaxScaler
scaler = MinMaxScaler()
然后,使用 fit 方法拟合缩放器(scaler),并将其应用于训练数据。对于 MinMaxScaler 来说, fit 方法计算训练集中每个特征的最大值和最小值。与第 2 章中的分类器和回归器(regressor)不同,在对缩放器调用 fit 时只提供了 X_train,而不用 y_train:
In[5]:
scaler.fit(X_train)
Out[5]:
MinMaxScaler(copy=True, feature_range=(0, 1))
为了应用刚刚学习的变换(即对训练数据进行实际缩放),我们使用缩放器的 transform方法。在 scikit-learn 中,每当模型返回数据的一种新表示时,都可以使用 transform方法:
In[6]:
# 变换数据
X_train_scaled = scaler.transform(X_train)
# 在缩放之前和之后分别打印数据集属性
print("transformed shape: {}".format(X_train_scaled.shape))
print("per-feature minimum before scaling:\n {}".format(X_train.min(axis=0)))
print("per-feature maximum before scaling:\n {}".format(X_train.max(axis=0)))
print("per-feature minimum after scaling:\n {}".format(
X_train_scaled.min(axis=0)))
print("per-feature maximum after scaling:\n {}".format(
X_train_scaled.max(axis=0)))
Out[6]:
transformed shape: (426, 30)
per-feature minimum before scaling:
[ 6.98 9.71 43.79 143.50 0.05 0.02 0. 0. 0.11
0.05 0.12 0.36 0.76 6.80 0. 0. 0. 0.
0.01 0. 7.93 12.02 50.41 185.20 0.07 0.03 0.
0. 0.16 0.06]
per-feature maximum before scaling:
[ 28.11 39.28 188.5 2501.0 0.16 0.29 0.43 0.2
0.300 0.100 2.87 4.88 21.98 542.20 0.03 0.14
0.400 0.050 0.06 0.03 36.04 49.54 251.20 4254.00
0.220 0.940 1.17 0.29 0.58 0.15]
per-feature minimum after scaling:
[ 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.]
per-feature maximum after scaling:
[ 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.]
变换后的数据形状与原始数据相同,特征只是发生了移动和缩放。你可以看到,现在所有特征都位于 0 到 1 之间,这也符合我们的预期。
为了将 SVM 应用到缩放后的数据上,还需要对测试集进行变换。这可以通过对 X_test 调用 transform 方法来完成:
In[7]:
# 对测试数据进行变换
X_test_scaled = scaler.transform(X_test)
# 在缩放之后打印测试数据的属性
print("per-feature minimum after scaling:\n{}".format(X_test_scaled.min(axis=0)))
print("per-feature maximum after scaling:\n{}".format(X_test_scaled.max(axis=0)))
Out[7]:
per-feature minimum after scaling:
[ 0.034 0.023 0.031 0.011 0.141 0.044 0. 0. 0.154 -0.006
-0.001 0.006 0.004 0.001 0.039 0.011 0. 0. -0.032 0.007
0.027 0.058 0.02 0.009 0.109 0.026 0. 0. -0. -0.002]
per-feature maximum after scaling:
[ 0.958 0.815 0.956 0.894 0.811 1.22 0.88 0.933 0.932 1.037
0.427 0.498 0.441 0.284 0.487 0.739 0.767 0.629 1.337 0.391
0.896 0.793 0.849 0.745 0.915 1.132 1.07 0.924 1.205 1.631]
你可以发现,对测试集缩放后的最大值和最小值不是 1 和 0,这或许有些出乎意料。有些特征甚至在 0~1 的范围之外!对此的解释是, MinMaxScaler(以及其他所有缩放器)总是对训练集和测试集应用完全相同的变换。也就是说, transform 方法总是减去训练集的最小值,然后除以训练集的范围,而这两个值可能与测试集的最小值和范围并不相同。
3.3.3 对训练数据和测试数据进行相同的缩放
为了让监督模型能够在测试集上运行,对训练集和测试集应用完全相同的变换是很重要的。如果我们使用测试集的最小值和范围,下面这个例子(图 3-2)展示了会发生什么:
In[8]:
from sklearn.datasets import make_blobs
# 构造数据
X, _ = make_blobs(n_samples=50, centers=5, random_state=4, cluster_std=2)
# 将其分为训练集和测试集
X_train, X_test = train_test_split(X, random_state=5, test_size=.1)
# 绘制训练集和测试集
fig, axes = plt.subplots(1, 3, figsize=(13, 4))
axes[0].scatter(X_train[:, 0], X_train[:, 1],
c=mglearn.cm2(0), label="Training set", s=60)
axes[0].scatter(X_test[:, 0], X_test[:, 1], marker='^',
c=mglearn.cm2(1), label="Test set", s=60)
axes[0].legend(loc='upper left')
axes[0].set_title("Original Data")
# 利用MinMaxScaler缩放数据
scaler = MinMaxScaler()
scaler.fit(X_train)
X_train_scaled = scaler.transform(X_train)
X_test_scaled = scaler.transform(X_test)
# 将正确缩放的数据可视化
axes[1].scatter(X_train_scaled[:, 0], X_train_scaled[:, 1],
c=mglearn.cm2(0), label="Training set", s=60)
axes[1].scatter(X_test_scaled[:, 0], X_test_scaled[:, 1], marker='^',
c=mglearn.cm2(1), label="Test set", s=60)
axes[1].set_title("Scaled Data")
# 单独对测试集进行缩放
# 使得测试集的最小值为0,最大值为1
# 千万不要这么做!这里只是为了举例
test_scaler = MinMaxScaler()
test_scaler.fit(X_test)
X_test_scaled_badly = test_scaler.transform(X_test)
# 将错误缩放的数据可视化
axes[2].scatter(X_train_scaled[:, 0], X_train_scaled[:, 1],
c=mglearn.cm2(0), label="training set", s=60)
axes[2].scatter(X_test_scaled_badly[:, 0], X_test_scaled_badly[:, 1],
marker='^', c=mglearn.cm2(1), label="test set", s=60)
axes[2].set_title("Improperly Scaled Data")
for ax in axes:
ax.set_xlabel("Feature 0")
ax.set_ylabel("Feature 1")
第一张图是未缩放的二维数据集,其中训练集用圆形表示,测试集用三角形表示。第二张图中是同样的数据,但使用 MinMaxScaler 缩放。这里我们调用 fit 作用在训练集上,然后调用 transform 作用在训练集和测试集上。你可以发现,第二张图中的数据集看起来与第一张图中的完全相同,只是坐标轴刻度发生了变化。现在所有特征都位于 0 到 1 之间。你还可以发现,测试数据(三角形)的特征最大值和最小值并不是 1 和 0。
第三张图展示了如果我们对训练集和测试集分别进行缩放会发生什么。在这种情况下,对训练集和测试集而言,特征的最大值和最小值都是 1 和 0。但现在数据集看起来不一样。测试集相对训练集的移动不一致,因为它们分别做了不同的缩放。我们随意改变了数据的排列。这显然不是我们想要做的事情。
再换一种思考方式,想象你的测试集只有一个点。对于一个点而言,无法将其正确地缩放以满足 MinMaxScaler 的最大值和最小值的要求。但是,测试集的大小不应该对你的处理方式有影响。
3.3.4 预处理对监督学习的作用
现在我们回到 cancer 数据集,观察使用 MinMaxScaler 对学习 SVC 的作用(这是一种不同的方法,实现了与第 2 章中相同的缩放)。首先,为了对比,我们再次在原始数据上拟合 SVC:
In[10]:
from sklearn.svm import SVC
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target,
random_state=0)
svm = SVC(C=100)
svm.fit(X_train, y_train)
print("Test set accuracy: {:.2f}".format(svm.score(X_test, y_test)))
Out[10]:
Test set accuracy: 0.63
下面先用 MinMaxScaler 对数据进行缩放,然后再拟合 SVC:
In[11]:
# 使用0-1缩放进行预处理
scaler = MinMaxScaler()
scaler.fit(X_train)
X_train_scaled = scaler.transform(X_train)
X_test_scaled = scaler.transform(X_test)
# 在缩放后的训练数据上学习SVM
svm.fit(X_train_scaled, y_train)
# 在缩放后的测试集上计算分数
print("Scaled test set accuracy: {:.2f}".format(
svm.score(X_test_scaled, y_test)))
Out[11]:
Scaled test set accuracy: 0.97
正如我们上面所见,数据缩放的作用非常显著。虽然数据缩放不涉及任何复杂的数学,但良好的做法仍然是使用 scikit-learn 提供的缩放机制,而不是自己重新实现它们,因为即使在这些简单的计算中也容易犯错。
你也可以通过改变使用的类将一种预处理算法轻松替换成另一种,因为所有的预处理类都具有相同的接口,都包含 fit 和 transform 方法:
In[12]:
# 利用零均值和单位方差的缩放方法进行预处理
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
scaler.fit(X_train)
X_train_scaled = scaler.transform(X_train)
X_test_scaled = scaler.transform(X_test)
# 在缩放后的训练数据上学习SVM
svm.fit(X_train_scaled, y_train)
# 在缩放后的测试集上计算分数
print("SVM test accuracy: {:.2f}".format(svm.score(X_test_scaled, y_test)))
Out[12]:
SVM test accuracy: 0.96
前面我们已经看到了用于预处理的简单数据变换的工作原理,下面继续学习利用无监督学习进行更有趣的变换。
3.4 降维、 特征提取与流形学习
前面讨论过,利用无监督学习进行数据变换可能有很多种目的。最常见的目的就是可视化、压缩数据,以及寻找信息量更大的数据表示以用于进一步的处理。
为了实现这些目的,最简单也最常用的一种算法就是主成分分析。我们也将学习另外两种算法:非负矩阵分解(NMF)和 t-SNE,前者通常用于特征提取,后者通常用于二维散点图的可视化。
3.4.1 主成分分析
主成分分析(principal component analysis, PCA)是一种旋转数据集的方法,旋转后的特征在统计上不相关。在做完这种旋转之后,通常是根据新特征对解释数据的重要性来选择它的一个子集。下面的例子(图 3-3)展示了 PCA 对一个模拟二维数据集的作用:
In[13]:
mglearn.plots.plot_pca_illustration()
第一张图(左上)显示的是原始数据点,用不同颜色加以区分。算法首先找到方差最大的方向,将其标记为“成分 1”(Component 1)。这是数据中包含最多信息的方向(或向量),换句话说,沿着这个方向的特征之间最为相关。然后,算法找到与第一个方向正交(成直角)且包含最多信息的方向。在二维空间中,只有一个成直角的方向,但在更高维的空间中会有(无穷)多的正交方向。虽然这两个成分都画成箭头,但其头尾的位置并不重要。我们也可以将第一个成分画成从中心指向左上,而不是指向右下。利用这一过程找到的方向被称为主成分(principal component),因为它们是数据方差的主要方向。一般来说,主成分的个数与原始特征相同。
第二张图(右上)显示的是同样的数据,但现在将其旋转,使得第一主成分与 x 轴平行且第二主成分与 y 轴平行。在旋转之前,从数据中减去平均值,使得变换后的数据以零为中心。在 PCA 找到的旋转表示中,两个坐标轴是不相关的,也就是说,对于这种数据表示,除了对角线,相关矩阵全部为零。
我们可以通过仅保留一部分主成分来使用 PCA 进行降维。在这个例子中,我们可以仅保留第一个主成分,正如图 3-3 中第三张图所示(左下)。这将数据从二维数据集降为一维数据集。但要注意,我们没有保留原始特征之一,而是找到了最有趣的方向(第一张图中从左上到右下)并保留这一方向,即第一主成分。
最后,我们可以反向旋转并将平均值重新加到数据中。这样会得到图 3-3 最后一张图中的数据。这些数据点位于原始特征空间中,但我们仅保留了第一主成分中包含的信息。这种变换有时用于去除数据中的噪声影响,或者将主成分中保留的那部分信息可视化。
1. 将 PCA 应用于 cancer 数据集并可视化
PCA 最常见的应用之一就是将高维数据集可视化。正如第 1 章中所说,对于有两个以上特征的数据,很难绘制散点图。对于 Iris(鸢尾花)数据集,我们可以创建散点图矩阵(见第 1 章图 1-3),通过展示特征所有可能的两两组合来展示数据的局部图像。但如果我们想要查看乳腺癌数据集,即便用散点图矩阵也很困难。这个数据集包含 30 个特征,这就导致需要绘制 30 * 14 = 420 张散点图!我们永远不可能仔细观察所有这些图像,更不用说试图理解它们了。
不过我们可以使用一种更简单的可视化方法——对每个特征分别计算两个类别(良性肿瘤和恶性肿瘤)的直方图(见图 3-4)。
In[14]:
fig, axes = plt.subplots(15, 2, figsize=(10, 20))
malignant = cancer.data[cancer.target == 0]
benign = cancer.data[cancer.target == 1]
ax = axes.ravel()
for i in range(30):
_, bins = np.histogram(cancer.data[:, i], bins=50)
ax[i].hist(malignant[:, i], bins=bins, color=mglearn.cm3(0), alpha=.5)
ax[i].hist(benign[:, i], bins=bins, color=mglearn.cm3(2), alpha=.5)
ax[i].set_title(cancer.feature_names[i])
ax[i].set_yticks(())
ax[0].set_xlabel("Feature magnitude")
ax[0].set_ylabel("Frequency")
ax[0].legend(["malignant", "benign"], loc="best")
fig.tight_layout()
这里我们为每个特征创建一个直方图,计算具有某一特征的数据点在特定范围内(叫作bin)的出现频率。每张图都包含两个直方图,一个是良性类别的所有点(蓝色),一个是恶性类别的所有点(红色)。这样我们可以了解每个特征在两个类别中的分布情况,也可以猜测哪些特征能够更好地区分良性样本和恶性样本。例如,“smoothness error”特征似乎没有什么信息量,因为两个直方图大部分都重叠在一起,而“worst concave points”特征看起来信息量相当大,因为两个直方图的交集很小。
但是,这种图无法向我们展示变量之间的相互作用以及这种相互作用与类别之间的关系。
利用 PCA,我们可以获取到主要的相互作用,并得到稍为完整的图像。我们可以找到前两个主成分,并在这个新的二维空间中用散点图将数据可视化。在应用 PCA 之前,我们利用 StandardScaler 缩放数据,使每个特征的方差均为 1:
In[15]:
from sklearn.datasets import load_breast_cancer
cancer = load_breast_cancer()
scaler = StandardScaler()
scaler.fit(cancer.data)
X_scaled = scaler.transform(cancer.data)
学习并应用 PCA 变换与应用预处理变换一样简单。我们将 PCA 对象实例化,调用 fit 方法找到主成分,然后调用 transform 来旋转并降维。默认情况下, PCA 仅旋转(并移动)数据,但保留所有的主成分。为了降低数据的维度,我们需要在创建 PCA 对象时指定想要保留的主成分个数:
In[16]:
from sklearn.decomposition import PCA
# 保留数据的前两个主成分
pca = PCA(n_components=2)
# 对乳腺癌数据拟合PCA模型
pca.fit(X_scaled)
# 将数据变换到前两个主成分的方向上
X_pca = pca.transform(X_scaled)
print("Original shape: {}".format(str(X_scaled.shape)))
print("Reduced shape: {}".format(str(X_pca.shape)))
Out[16]:
Original shape: (569, 30)
Reduced shape: (569, 2)
现在我们可以对前两个主成分作图(图 3-5):
In[17]:
# 对第一个和第二个主成分作图,按类别着色
plt.figure(figsize=(8, 8))
mglearn.discrete_scatter(X_pca[:, 0], X_pca[:, 1], cancer.target)
plt.legend(cancer.target_names, loc="best")
plt.gca().set_aspect("equal")
plt.xlabel("First principal component")
plt.ylabel("Second principal component")
重要的是要注意, PCA 是一种无监督方法,在寻找旋转方向时没有用到任何类别信息。它只是观察数据中的相关性。对于这里所示的散点图,我们绘制了第一主成分与第二主成分的关系,然后利用类别信息对数据点进行着色。你可以看到,在这个二维空间中两个类别被很好地分离。这让我们相信,即使是线性分类器(在这个空间中学习一条直线)也可以在区分这个两个类别时表现得相当不错。我们还可以看到,恶性点比良性点更加分散,这一点也可以在图 3-4 的直方图中看出来。
PCA 的一个缺点在于,通常不容易对图中的两个轴做出解释。主成分对应于原始数据中的方向,所以它们是原始特征的组合。但这些组合往往非常复杂,这一点我们很快就会看到。在拟合过程中,主成分被保存在 PCA 对象的 components_ 属性中:
In[18]:
print("PCA component shape: {}".format(pca.components_.shape))
Out[18]:
PCA component shape: (2, 30)
components_ 中的每一行对应于一个主成分,它们按重要性排序(第一主成分排在首位,以此类推)。列对应于 PCA 的原始特征属性,在本例中即为“mean radius”“mean texture”等。我们来看一下 components_ 的内容:
In[19]:
print("PCA components:\n{}".format(pca.components_))
Out[19]:
PCA components:
[[ 0.219 0.104 0.228 0.221 0.143 0.239 0.258 0.261 0.138 0.064
0.206 0.017 0.211 0.203 0.015 0.17 0.154 0.183 0.042 0.103
0.228 0.104 0.237 0.225 0.128 0.21 0.229 0.251 0.123 0.132]
[-0.234 -0.06 -0.215 -0.231 0.186 0.152 0.06 -0.035 0.19 0.367
-0.106 0.09 -0.089 -0.152 0.204 0.233 0.197 0.13 0.184 0.28
-0.22 -0.045 -0.2 -0.219 0.172 0.144 0.098 -0.008 0.142 0.275]]
我们还可以用热图将系数可视化(图 3-6),这可能更容易理解:
In[20]:
plt.matshow(pca.components_, cmap='viridis')
plt.yticks([0, 1], ["First component", "Second component"])
plt.colorbar()
plt.xticks(range(len(cancer.feature_names)),
cancer.feature_names, rotation=60, ha='left')
plt.xlabel("Feature")
plt.ylabel("Principal components")
你可以看到,在第一个主成分中,所有特征的符号相同(均为正,但前面我们提到过,箭头指向哪个方向无关紧要)。这意味着在所有特征之间存在普遍的相关性。如果一个测量值较大的话,其他的测量值可能也较大。第二个主成分的符号有正有负,而且两个主成分都包含所有 30 个特征。这种所有特征的混合使得解释图 3-6 中的坐标轴变得十分困难。
2. 特征提取的特征脸
前面提到过, PCA 的另一个应用是特征提取。特征提取背后的思想是,可以找到一种数据表示,比给定的原始表示更适合于分析。特征提取很有用,它的一个很好的应用实例就是图像。图像由像素组成,通常存储为红绿蓝(RGB)强度。图像中的对象通常由上千个像素组成,它们只有放在一起才有意义。
我们将给出用 PCA 对图像做特征提取的一个简单应用,即处理 Wild 数据集 Labeled Faces(标记人脸)中的人脸图像。这一数据集包含从互联网下载的名人脸部图像,它包含从 21世纪初开始的政治家、歌手、演员和运动员的人脸图像。我们使用这些图像的灰度版本,并将它们按比例缩小以加快处理速度。你可以在图 3-7 中看到其中一些图像:
from sklearn.datasets import fetch_lfw_people
people = fetch_lfw_people(min_faces_per_person=20, resize=0.7)
image_shape = people.images[0].shape
fig, axes = plt.subplots(2, 5, figsize=(15, 8),
subplot_kw={'xticks': (), 'yticks': ()})
for target, image, ax in zip(people.target, people.images, axes.ravel()):
ax.imshow(image)
ax.set_title(people.target_names[target])
一共有 3023 张图像,每张大小为 87 像素 × 65 像素,分别属于 62 个不同的人:
In[22]:
print("people.images.shape: {}".format(people.images.shape))
print("Number of classes: {}".format(len(people.target_names)))
Out[22]:
people.images.shape: (3023, 87, 65)
Number of classes: 62
但这个数据集有些偏斜,其中包含 George W. Bush(小布什)和 Colin Powell(科林 • 鲍威尔)的大量图像,正如你在下面所见:
In[23]:
# 计算每个目标出现的次数
counts = np.bincount(people.target)
# 将次数与目标名称一起打印出来
for i, (count, name) in enumerate(zip(counts, people.target_names)):
print("{0:25} {1:3}".format(name, count), end=' ')
if (i + 1) % 3 == 0:
print()
为了降低数据偏斜,我们对每个人最多只取 50 张图像(否则,特征提取将会被 George W.Bush 的可能性大大影响):
In[24]:
mask = np.zeros(people.target.shape, dtype=np.bool)
for target in np.unique(people.target):
mask[np.where(people.target == target)[0][:50]] = 1
X_people = people.data[mask]
y_people = people.target[mask]
# 将灰度值缩放到0到1之间,而不是在0到255之间
# 以得到更好的数据稳定性
X_people = X_people / 255.
人脸识别的一个常见任务就是看某个前所未见的人脸是否属于数据库中的某个已知人物。
这在照片收集、社交媒体和安全应用中都有应用。解决这个问题的方法之一就是构建一个分类器,每个人都是一个单独的类别。但人脸数据库中通常有许多不同的人,而同一个人的图像很少(也就是说,每个类别的训练样例很少)。这使得大多数分类器的训练都很困难。另外,通常你还想要能够轻松添加新的人物,不需要重新训练一个大型模型。
一种简单的解决方法是使用单一最近邻分类器,寻找与你要分类的人脸最为相似的人脸。这个分类器原则上可以处理每个类别只有一个训练样例的情况。下面看一下KNeighborsClassifier 的表现如何:
In[25]:
from sklearn.neighbors import KNeighborsClassifier
# 将数据分为训练集和测试集
X_train, X_test, y_train, y_test = train_test_split(
X_people, y_people, stratify=y_people, random_state=0)
# 使用一个邻居构建KNeighborsClassifier
knn = KNeighborsClassifier(n_neighbors=1)
knn.fit(X_train, y_train)
print("Test set score of 1-nn: {:.2f}".format(knn.score(X_test, y_test)))
Out[25]:
Test set score of 1-nn: 0.27
我们得到的精度为 26.6%。对于包含 62 个类别的分类问题来说,这实际上不算太差(随机猜测的精度约为 1/62=1.5%),但也不算好。我们每识别四次仅正确识别了一个人。
这里就可以用到 PCA。想要度量人脸的相似度,计算原始像素空间中的距离是一种相当糟糕的方法。用像素表示来比较两张图像时,我们比较的是每个像素的灰度值与另一张图像对应位置的像素灰度值。这种表示与人们对人脸图像的解释方式有很大不同,使用这种原始表示很难获取到面部特征。例如,如果使用像素距离,那么将人脸向右移动一个像素将会发生巨大的变化,得到一个完全不同的表示。我们希望,使用沿着主成分方向的距离可以提高精度。这里我们启用 PCA 的白化(whitening)选项,它将主成分缩放到相同的尺度。变换后的结果与使用 StandardScaler 相同。再次使用图 3-3 中的数据,白化不仅对应于旋转数据,还对应于缩放数据使其形状是圆形而不是椭圆(参见图 3-8):
In[26]:
mglearn.plots.plot_pca_whitening()
我们对训练数据拟合 PCA 对象,并提取前 100 个主成分。然后对训练数据和测试数据进行变换:
In[27]:
pca = PCA(n_components=100, whiten=True, random_state=0).fit(X_train)
X_train_pca = pca.transform(X_train)
X_test_pca = pca.transform(X_test)
print("X_train_pca.shape: {}".format(X_train_pca.shape))
Out[27]:
X_train_pca.shape: (1547, 100)
新数据有 100 个特征,即前 100 个主成分。现在,可以对新表示使用单一最近邻分类器来将我们的图像分类:
In[28]:
knn = KNeighborsClassifier(n_neighbors=1)
knn.fit(X_train_pca, y_train)
print("Test set accuracy: {:.2f}".format(knn.score(X_test_pca, y_test)))
Out[28]:
Test set accuracy: 0.36
我们的精度有了相当显著的提高,从 26.6% 提升到 35.7%,这证实了我们的直觉,即主成分可能提供了一种更好的数据表示。
对于图像数据,我们还可以很容易地将找到的主成分可视化。请记住,成分对应于输入空间里的方向。这里的输入空间是 87 像素 × 65 像素的灰度图像,所以在这个空间中的方向也是 87 像素 × 65 像素的灰度图像。
我们来看一下前几个主成分(图 3-9):
In[29]:
print("pca.components_.shape: {}".format(pca.components_.shape))
Out[29]:
pca.components_.shape: (100, 5655)
In[30]:
fig, axes = plt.subplots(3, 5, figsize=(15, 12),
subplot_kw={'xticks': (), 'yticks': ()})
for i, (component, ax) in enumerate(zip(pca.components_, axes.ravel())):
ax.imshow(component.reshape(image_shape),
cmap='viridis')
ax.set_title("{}. component".format((i + 1)))
虽然我们肯定无法理解这些成分的所有内容,但可以猜测一些主成分捕捉到了人脸图像的哪些方面。第一个主成分似乎主要编码的是人脸与背景的对比,第二个主成分编码的是人脸左半部分和右半部分的明暗程度差异,如此等等。虽然这种表示比原始像素值的语义稍强,但它仍与人们感知人脸的方式相去甚远。由于 PCA 模型是基于像素的,因此人脸的相对位置(眼睛、下巴和鼻子的位置)和明暗程度都对两张图像在像素表示中的相似程度有很大影响。但人脸的相对位置和明暗程度可能并不是人们首先感知的内容。在要求人们评价人脸的相似度时,他们更可能会使用年龄、性别、面部表情和发型等属性,而这些属性很难从像素强度中推断出来。重要的是要记住,算法对数据(特别是视觉数据,比如人们非常熟悉的图像)的解释通常与人类的解释方式大不相同。
不过让我们回到 PCA 的具体案例。我们对 PCA 变换的介绍是:先旋转数据,然后删除方差较小的成分。另一种有用的解释是尝试找到一些数字(PCA 旋转后的新特征值),使我们可以将测试点表示为主成分的加权求和(见图 3-10)。
这里 x0、 x1 等是这个数据点的主成分的系数,换句话说,它们是图像在旋转后的空间中的表示。
我们还可以用另一种方法来理解 PCA 模型,就是仅使用一些成分对原始数据进行重建。在图 3-3 中,在去掉第二个成分并来到第三张图之后,我们反向旋转并重新加上平均值,这样就在原始空间中获得去掉第二个成分的新数据点,正如最后一张图所示。我们可以对人脸做类似的变换,将数据降维到只包含一些主成分,然后反向旋转回到原始空间。回到原始特征空间可以通过 inverse_transform 方法来实现。这里我们分别利用 10 个、 50 个、100 个和 500 个成分对一些人脸进行重建并将其可视化(图 3-11):
In[32]:
mglearn.plots.plot_pca_faces(X_train, X_test, image_shape)
可以看到,在仅使用前 10 个主成分时,仅捕捉到了图片的基本特点,比如人脸方向和明暗程度。随着使用的主成分越来越多,图像中也保留了越来越多的细节。这对应于图 3-10的求和中包含越来越多的项。如果使用的成分个数与像素个数相等,意味着我们在旋转后不会丢弃任何信息,可以完美重建图像。
我们还可以尝试使用 PCA 的前两个主成分,将数据集中的所有人脸在散点图中可视化(图 3-12),其类别在图中给出。这与我们对 cancer 数据集所做的类似:
In[33]:
mglearn.discrete_scatter(X_train_pca[:, 0], X_train_pca[:, 1], y_train)
plt.xlabel("First principal component")
plt.ylabel("Second principal component")
如你所见,如果我们只使用前两个主成分,整个数据只是一大团,看不到类别之间的分界。这并不意外,因为即使有 10 个成分(正如图 3-11 所示), PCA 也仅捕捉到人脸非常粗略的特征。
3.4.2 非负矩阵分解
非负矩阵分解(non-negative matrix factorization, NMF)是另一种无监督学习算法,其目的在于提取有用的特征。它的工作原理类似于 PCA,也可以用于降维。与 PCA 相同,我们试图将每个数据点写成一些分量的加权求和,正如图 3-10 所示。但在 PCA 中,我们想要的是正交分量,并且能够解释尽可能多的数据方差;而在 NMF 中,我们希望分量和系数均为非负,也就是说,我们希望分量和系数都大于或等于 0。因此,这种方法只能应用于每个特征都是非负的数据,因为非负分量的非负求和不可能变为负值。
将数据分解成非负加权求和的这个过程,对由多个独立源相加(或叠加)创建而成的数据特别有用,比如多人说话的音轨或包含多种乐器的音乐。在这种情况下, NMF 可以识别出组成合成数据的原始分量。总的来说,与 PCA 相比, NMF 得到的分量更容易解释,因为负的分量和系数可能会导致难以解释的抵消效应(cancellation effect)。举个例子,图3-9 中的特征脸同时包含正数和负数,我们在 PCA 的说明中也提到过,正负号实际上是任意的。在将 NMF 应用于人脸数据集之前,我们先来简要回顾一下模拟数据。
1. 将 NMF 应用于模拟数据
与使用 PCA 不同,我们需要保证数据是正的, NMF 能够对数据进行操作。这说明数据相对于原点 (0, 0) 的位置实际上对 NMF 很重要。因此,你可以将提取出来的非负分量看作是从 (0, 0) 到数据的方向。
下面的例子(图 3-13)给出了 NMF 在二维玩具数据上的结果:
In[34]:
mglearn.plots.plot_nmf_illustration()
对于两个分量的 NMF(如左图所示),显然所有数据点都可以写成这两个分量的正数组合。如果有足够多的分量能够完美地重建数据(分量个数与特征个数相同),那么算法会选择指向数据极值的方向。
如果我们仅使用一个分量,那么 NMF 会创建一个指向平均值的分量,因为指向这里可以对数据做出最好的解释。你可以看到,与 PCA 不同,减少分量个数不仅会删除一些方向,而且会创建一组完全不同的分量! NMF 的分量也没有按任何特定方法排序,所以不存在“第一非负分量”:所有分量的地位平等。
NMF 使用了随机初始化,根据随机种子的不同可能会产生不同的结果。在相对简单的情况下(比如两个分量的模拟数据),所有数据都可以被完美地解释,那么随机性的影响很小(虽然可能会影响分量的顺序或尺度)。在更加复杂的情况下,影响可能会很大。
2. 将 NMF 应用于人脸图像
现在我们将 NMF 应用于之前用过的 Wild 数据集中的 Labeled Faces。 NMF 的主要参数是我们想要提取的分量个数。通常来说,这个数字要小于输入特征的个数(否则的话,将每个像素作为单独的分量就可以对数据进行解释)。
首先,我们来观察分量个数如何影响 NMF 重建数据的好坏(图 3-14):
In[35]:
mglearn.plots.plot_nmf_faces(X_train, X_test, image_shape)
反向变换的数据质量与使用 PCA 时类似,但要稍差一些。这是符合预期的,因为 PCA 找到的是重建的最佳方向。 NMF 通常并不用于对数据进行重建或编码,而是用于在数据中寻找有趣的模式。
我们尝试仅提取一部分分量(比如 15 个),初步观察一下数据。其结果见图 3-15。
In[36]:
from sklearn.decomposition import NMF
nmf = NMF(n_components=15, random_state=0)
nmf.fit(X_train)
X_train_nmf = nmf.transform(X_train)
X_test_nmf = nmf.transform(X_test)
fig, axes = plt.subplots(3, 5, figsize=(15, 12),
subplot_kw={'xticks': (), 'yticks': ()})
for i, (component, ax) in enumerate(zip(nmf.components_, axes.ravel())):
ax.imshow(component.reshape(image_shape))
ax.set_title("{}. component".format(i))
这些分量都是正的,因此比图 3-9 所示的 PCA 分量更像人脸原型。例如,你可以清楚地看到,分量 3(component 3)显示了稍微向右转动的人脸,而分量 7(component 7)则显示了稍微向左转动的人脸。我们来看一下这两个分量特别大的那些图像,分别如图 3-16 和图3-17 所示。
In[37]:
compn = 3
# 按第3个分量排序,绘制前10张图像
inds = np.argsort(X_train_nmf[:, compn])[::-1]
fig, axes = plt.subplots(2, 5, figsize=(15, 8),
subplot_kw={'xticks': (), 'yticks': ()})
fig.suptitle("Large component 3")
for i, (ind, ax) in enumerate(zip(inds, axes.ravel())):
ax.imshow(X_train[ind].reshape(image_shape))
compn = 7
# 按第7个分量排序,绘制前10张图像
inds = np.argsort(X_train_nmf[:, compn])[::-1]
fig.suptitle("Large component 7")
fig, axes = plt.subplots(2, 5, figsize=(15, 8),
subplot_kw={'xticks': (), 'yticks': ()})
for i, (ind, ax) in enumerate(zip(inds, axes.ravel())):
ax.imshow(X_train[ind].reshape(image_shape))
正如所料,分量 3 系数较大的人脸都是向右看的人脸(图 3-16),而分量 7 系数较大的人脸都向左看(图 3-17)。如前所述,提取这样的模式最适合于具有叠加结构的数据,包括音频、基因表达和文本数据。我们通过一个模拟数据的例子来看一下这种用法。
假设我们对一个信号感兴趣,它是三个不同信号源合成的(图 3-18):
In[38]:
S = mglearn.datasets.make_signals()
plt.figure(figsize=(6, 1))
plt.plot(S, '-')
plt.xlabel("Time")
plt.ylabel("Signal")
不幸的是,我们无法观测到原始信号,只能观测到三个信号的叠加混合。我们想要将混合信号分解为原始分量。假设我们有许多种不同的方法来观测混合信号(比如有 100 台测量装置),每种方法都为我们提供了一系列测量结果。
In[39]:
# 将数据混合成100维的状态
A = np.random.RandomState(0).uniform(size=(100, 3))
X = np.dot(S, A.T)
print("Shape of measurements: {}".format(X.shape))
Out[39]:
Shape of measurements: (2000, 100)
我们可以用 NMF 来还原这三个信号:
In[40]:
nmf = NMF(n_components=3, random_state=42)
S_ = nmf.fit_transform(X)
print("Recovered signal shape: {}".format(S_.shape))
Out[40]:
Recovered signal shape: (2000, 3)
为了对比,我们也应用了 PCA:
In[41]:
pca = PCA(n_components=3)
H = pca.fit_transform(X)
图 3-19 给出了 NMF 和 PCA 发现的信号活动:
In[42]:
models = [X, S, S_, H]
names = ['Observations (first three measurements)',
'True sources',
'NMF recovered signals',
'PCA recovered signals']
fig, axes = plt.subplots(4, figsize=(8, 4), gridspec_kw={'hspace': .5},
subplot_kw={'xticks': (), 'yticks': ()})
for model, name, ax in zip(models, names, axes):
ax.set_title(name)
ax.plot(model[:, :3], '-')
图中包含来自 X 的 100 次测量中的 3 次,用于参考。可以看到, NMF 在发现原始信号源时得到了不错的结果,而 PCA 则失败了,仅使用第一个成分来解释数据中的大部分变化。
要记住, NMF 生成的分量是没有顺序的。在这个例子中, NMF 分量的顺序与原始信号完全相同(参见三条曲线的颜色),但这纯属偶然。
还有许多其他算法可用于将每个数据点分解为一系列固定分量的加权求和,正如 PCA 和NMF 所做的那样。讨论所有这些算法已超出了本书的范围,而且描述对分量和系数的约束通常要涉及概率论。如果你对这种类型的模式提取感兴趣,我们推荐你学习 scikitlearn 用户指南中关于独立成分分析(ICA)、因子分析(FA)和稀疏编码(字典学习)等的内容,所有这些内容都可以在关于分解方法的页面中找到(http://scikit-learn.org/stable/modules/decomposition.html)。
3.4.3 用t-SNE进行流形学习
虽然 PCA 通常是用于变换数据的首选方法,使你能够用散点图将其可视化,但这一方法的性质(先旋转然后减少方向)限制了其有效性,正如我们在 Wild 数据集 Labeled Faces的散点图中所看到的那样。有一类用于可视化的算法叫作流形学习算法(manifold learning algorithm),它允许进行更复杂的映射, 通常也可以给出更好的可视化。其中特别有用的一个就是 t-SNE 算法。
流形学习算法主要用于可视化,因此很少用来生成两个以上的新特征。其中一些算法(包括 t-SNE)计算训练数据的一种新表示,但不允许变换新数据。这意味着这些算法不能用于测试集:更确切地说,它们只能变换用于训练的数据。流形学习对探索性数据分析是很有用的,但如果最终目标是监督学习的话,则很少使用。 t-SNE 背后的思想是找到数据的一个二维表示,尽可能地保持数据点之间的距离。 t-SNE 首先给出每个数据点的随机二维表示,然后尝试让在原始特征空间中距离较近的点更加靠近,原始特征空间中相距较远的点更加远离。 t-SNE 重点关注距离较近的点,而不是保持距离较远的点之间的距离。换句话说,它试图保存那些表示哪些点比较靠近的信息。
我们将对 scikit-learn 包含的一个手写数字数据集 2 应用 t-SNE 流形学习算法。在这个数据集中,每个数据点都是 0 到 9 之间手写数字的一张 8× 8 灰度图像。图 3-20 给出了每个类别的一个例子。
In[43]:
from sklearn.datasets import load_digits
digits = load_digits()
fig, axes = plt.subplots(2, 5, figsize=(10, 5),
subplot_kw={'xticks':(), 'yticks': ()})
for ax, img in zip(axes.ravel(), digits.images):
ax.imshow(img)
我们用 PCA 将降到二维的数据可视化。我们对前两个主成分作图,并按类别对数据点着色(图 3-21):
In[44]:
# 构建一个PCA模型
pca = PCA(n_components=2)
pca.fit(digits.data)
# 将digits数据变换到前两个主成分的方向上
digits_pca = pca.transform(digits.data)
colors = ["#476A2A", "#7851B8", "#BD3430", "#4A2D4E", "#875525",
"#A83683", "#4E655E", "#853541", "#3A3120", "#535D8E"]
plt.figure(figsize=(10, 10))
plt.xlim(digits_pca[:, 0].min(), digits_pca[:, 0].max())
plt.ylim(digits_pca[:, 1].min(), digits_pca[:, 1].max())
for i in range(len(digits.data)):
# 将数据实际绘制成文本,而不是散点
plt.text(digits_pca[i, 0], digits_pca[i, 1], str(digits.target[i]),
color = colors[digits.target[i]],
fontdict={'weight': 'bold', 'size': 9})
plt.xlabel("First principal component")
plt.ylabel("Second principal component")
实际上,这里我们用每个类别对应的数字作为符号来显示每个类别的位置。利用前两个主成分可以将数字 0、 6 和 4 相对较好地分开,尽管仍有重叠。大部分其他数字都大量重叠在一起。
我们将 t-SNE 应用于同一个数据集,并对结果进行比较。由于 t-SNE 不支持变换新数据,所以 TSNE 类没有 transform 方法。我们可以调用 fit_transform 方法来代替,它会构建模型并立刻返回变换后的数据(见图 3-22):
In[45]:
from sklearn.manifold import TSNE
tsne = TSNE(random_state=42)
# 使用fit_transform而不是fit,因为TSNE没有transform方法
digits_tsne = tsne.fit_transform(digits.data)
In[46]:
plt.figure(figsize=(10, 10))
plt.xlim(digits_tsne[:, 0].min(), digits_tsne[:, 0].max() + 1)
plt.ylim(digits_tsne[:, 1].min(), digits_tsne[:, 1].max() + 1)
for i in range(len(digits.data)):
# 将数据实际绘制成文本,而不是散点
plt.text(digits_tsne[i, 0], digits_tsne[i, 1], str(digits.target[i]),
color = colors[digits.target[i]],
fontdict={'weight': 'bold', 'size': 9})
plt.xlabel("t-SNE feature 0")
plt.xlabel("t-SNE feature 1")
t-SNE 的结果非常棒。所有类别都被明确分开。数字 1 和 9 被分成几块,但大多数类别都形成一个密集的组。要记住,这种方法并不知道类别标签:它完全是无监督的。但它能够找到数据的一种二维表示,仅根据原始空间中数据点之间的靠近程度就能够将各个类别明确分开。
t-SNE 算法有一些调节参数,虽然默认参数的效果通常就很好。你可以尝试修改perplexity 和 early_exaggeration,但作用一般很小。