第五篇 降维:深入主成分分析

降维:深入主成分分析

这里我们将探索主成分分析(PCA), 这是非常有用的线性维度下降技术。

我们先导入一些标准模块:

1
2
3
4
5
6
7
8
9
from __future__ import print_function, division
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from scipy import stats
# use seaborn plotting style defaults
import seaborn as sns; sns.set()

介绍主成分分析(Principal components analysis)

主成分分析是一种非常强大的用来对数据降维的非监督方法。最简单的是看下面一组二维数据的可视化过程:

1
2
3
4
np.random.seed(1)
X = np.dot(np.random.random(size=(2, 2)), np.random.normal(size=(2, 200))).T
plt.plot(X[:, 0], X[:, 1], 'o')
plt.axis('equal');

png

我们能看到在数据里面包含了确定的趋势。PCA能在数据中发现关键的成分,并能够描述那些成分在数据分布中的重要性。

1
2
3
4
5
6
7
from sklearn.decomposition import PCA
# 特征维度数目
pca = PCA(n_components=2)
pca.fit(X)
# 降维后的各主成分的方差值
print(pca.explained_variance_)
print(pca.components_)
[ 0.75871884  0.01838551]
[[-0.94446029 -0.32862557]
 [-0.32862557  0.94446029]]

看下这些值是什么意思,让我们把它们作为向量绘制出来:

1
2
3
4
5
plt.plot(X[:, 0], X[:, 1], 'o', alpha=0.5)
for length, vector in zip(pca.explained_variance_, pca.components_):
v = vector * 3 * np.sqrt(length)
plt.plot([0, v[0]], [0, v[1]], '-k', lw=3)
plt.axis('equal');

png

注意到一个向量比另外一个向量要长。直觉告诉我们这个方向比另外一个方向更加重要。各个主成分的方差(explainedvariance)量化了这些方向的重要程度。

我们可以认为第二重要的成分可以完全忽略也不会损失什么信息。我们看一下如果有95%的差异性被被保留的情况。

1
2
3
4
clf = PCA(0.95) # keep 95% of variance
X_trans = clf.fit_transform(X)
print(X.shape)
print(X_trans.shape)
(200, 2)
(200, 1)

通过丢弃5%的差异性,数据被压缩了50%。让我们来看一下压缩后的数据:

1
2
3
4
X_new = clf.inverse_transform(X_trans)
plt.plot(X[:, 0], X[:, 1], 'o', alpha=0.2)
plt.plot(X_new[:, 0], X_new[:, 1], 'ob', alpha=0.8)
plt.axis('equal');

png

浅色的点是原始数据,而深色的点是投影后的数据。我们可以看到通过丢弃数据集中5%的差异,然后重新投影,最重要的特征能够保留,并且数据也被压缩了50%!

应用PCA到数字图形上

刚刚列举的基于二维数据的降维看上去有点抽象,但是在高维数据的可视化中投影和降维是非常有用的。让我们看下如何把PCA应用在我们之前见过的数字的例子上:

1
2
3
4
5
from sklearn.datasets import load_digits
digits = load_digits()
X = digits.data
y = digits.target
print(X)
[[  0.   0.   5. ...,   0.   0.   0.]
 [  0.   0.   0. ...,  10.   0.   0.]
 [  0.   0.   0. ...,  16.   9.   0.]
 ..., 
 [  0.   0.   1. ...,   6.   0.   0.]
 [  0.   0.   2. ...,  12.   0.   0.]
 [  0.   0.  10. ...,  12.   1.   0.]]
1
2
3
4
pca = PCA(2) # project from 64 to 2 dimensions
Xproj = pca.fit_transform(X)
print(X.shape)
print(Xproj.shape)
(1797, 64)
(1797, 2)
1
2
3
plt.scatter(Xproj[:, 0], Xproj[:, 1], c=y, edgecolor='none', alpha=0.5,
cmap=plt.cm.get_cmap('nipy_spectral', 10))
plt.colorbar();

png

我们能从中发现数字之间的关系。本质上,我们把这些数字图形的64维空间摊平了。

成分(Components)是什么意思?

PCA是非常有用的降维算法,因为它是特征向量的直观表示。输入数据可以表示为一个向量,以下面的数字为例:

$$
x = [x_1, x_2, x_3 \cdots]
$$

它的真正含义是:

$$
image(x) = x_1 \cdot{\rm (pixel~1)} + x_2 \cdot{\rm (pixel~2)} + x_3 \cdot{\rm (pixel~3)} \cdots
$$

如果我们把像素空间的维度减少到6,我们只能复原部分的图像:

1
2
3
4
from fig_code.figures import plot_image_components
sns.set_style('white')
plot_image_components(digits.data[0])

png

像素级别的表示不是唯一的选择。我们也使用其他的basis函数,像下面的函数:

$$
image(x) = {\rm mean} + x_1 \cdot{\rm (basis~1)} + x_2 \cdot{\rm (basis~2)} + x_3 \cdot{\rm (basis~3)} \cdots
$$

PCA做的事情是对basis函数进行优化,以便用少的特征得到近似的估计。

1
2
from fig_code.figures import plot_pca_interactive
plot_pca_interactive(digits.data)

png

这里我们看到我们用了6个成分对输入数字做了近似的估计!

因此我们从两方面来看待PCA。一方面把它看做维度约减,另一方面它能被作为一种数据压缩,去掉没用的噪音。通过这种方式,PCA也能被当做filtering

选择成分的数量

我们到底要丢弃多少信息呢?我们通过分析成分方差来找出答案:

1
2
3
4
5
sns.set()
pca = PCA().fit(X)
plt.plot(np.cumsum(pca.explained_variance_ratio_))
plt.xlabel('number of components')
plt.ylabel('cumulative explained variance');

png

这里我们二维的投影损失了非常多的信息,我们需要20种成分来保留90%的差异性。

用PCA做数据压缩

我们刚才提到,PCA能作为数据压缩的一种方式。使用小的n_components允许你用主要的向量来表示高维的点。

这里通过改变成分的数量来观察数字0。

1
2
3
4
5
6
7
8
9
10
11
12
fig, axes = plt.subplots(8, 8, figsize=(8, 8))
fig.subplots_adjust(hspace=0.1, wspace=0.1)
for i, ax in enumerate(axes.flat):
pca = PCA(i + 1).fit(X)
im = pca.inverse_transform(pca.transform(X[20:21]))
ax.imshow(im.reshape((8, 8)), cmap='binary')
ax.text(0.95, 0.05, 'n = {0}'.format(i + 1), ha='right',
transform=ax.transAxes, color='green')
ax.set_xticks([])
ax.set_yticks([])

png

让我们通过IPython的interact功能来看下当成分数量(n_components)变化时数字的效果。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
from IPython.html.widgets import interact
def plot_digits(n_components):
fig = plt.figure(figsize=(8, 8))
plt.subplot(1, 1, 1, frameon=False, xticks=[], yticks=[])
nside = 10
pca = PCA(n_components).fit(X)
Xproj = pca.inverse_transform(pca.transform(X[:nside ** 2]))
Xproj = np.reshape(Xproj, (nside, nside, 8, 8))
total_var = pca.explained_variance_ratio_.sum()
im = np.vstack([np.hstack([Xproj[i, j] for j in range(nside)])
for i in range(nside)])
plt.imshow(im)
plt.grid(False)
plt.title("n = {0}, variance = {1:.2f}".format(n_components, total_var),
size=18)
plt.clim(0, 16)
interact(plot_digits, n_components=[1, 64], nside=[1, 8]);

png


Jupyter实现


知识共享许可协议
本作品采用知识共享署名-非商业性使用-禁止演绎 3.0 未本地化版本许可协议进行许可。