1. 维数灾难

物体在高维空间表现的十分不同

• 在高维超正方体中，大多数点都分布在边界处：在二维平面的一个正方形单元（一个 1×1 的正方形）中随机取一个点，那么随机选的点离所有边界大于 0.001（靠近中间位置）的概率为 0.4%（ 1 - 0.998^2 ）。在一个 1,0000 维的单位超正方体（一个 1×1×...×1 的立方体，有 10,000 个 1），这种可能性超过了 99.999999%。
• 高维数据集有很大风险分布的非常稀疏：在一个平方单位中随机选取两个点，那么这两个点之间的距离平均约为0.52，三维空间中这个距离是0.66，在一个 1,000,000 维超立方体中随机抽取两点，他们的距离大概是408.25。大多数训练实例可能彼此远离，训练集的维度越高，过拟合的风险就越大。

2. 降维的主要方法

投影 Projection

在大多数现实问题中，训练实例并不是在所有维度上均匀分布的。许多特征几乎是常数（如之前讨论的 MNIST，最外面一圈的边界像素基本都是255），所有训练实例实际上位于（或接近）高维空间的低维子空间内。

看个例子

一个分布接近于2D子空间的3D数据集

2D投影结果

# 准备数据
np.random.seed(4)
m = 60
w1, w2 = 0.1, 0.3
noise = 0.1

angles = np.random.rand(m) * 3 * np.pi / 2 - 0.5
X = np.empty((m, 3))
X[:, 0] = np.cos(angles) + np.sin(angles)/2 + noise * np.random.randn(m) / 2
X[:, 1] = np.sin(angles) * 0.7 + noise * np.random.randn(m) / 2
X[:, 2] = X[:, 0] * w1 + X[:, 1] * w2 + noise * np.random.randn(m)

# sklearn PCA
from sklearn.decomposition import PCA
'''
pca=PCA(n_components=0.95)  
将 n_components 设置为 0.0 到 1.0 之间的浮点数，表明希望保留的方差比率
'''
pca = PCA(n_components = 2)
X2D = pca.fit_transform(X)

'''
PC的方向
'''
pca.components_
# array([[-0.93636116, -0.29854881, -0.18465208],
#       [ 0.34027485, -0.90119108, -0.2684542 ]])
'''
每个主成的方差比例，越大说明越能区分。所有维度的方差比例之和为1
'''
pca.explained_variance_ratio_
# array([0.84248607, 0.14631839])

但如果子空间是旋转扭曲的，那么投影就不是一个好的降维方法，比如瑞士卷：

from sklearn.datasets import make_swiss_roll
X, t = make_swiss_roll(n_samples=1000, noise=0.2, random_state=42)

from sklearn.datasets import make_swiss_roll
X, t = make_swiss_roll(n_samples=1000, noise=0.2, random_state=42)

axes = [-11.5, 14, -2, 23, -12, 15]

fig = plt.figure(figsize=(6, 5))
ax = fig.add_subplot(111, projection='3d')

ax.scatter(X[:, 0], X[:, 1], X[:, 2], c=t, cmap=plt.cm.hot)
ax.view_init(10, -70)
ax.set_xlabel("$x_1$", fontsize=18)
ax.set_ylabel("$x_2$", fontsize=18)
ax.set_zlabel("$x_3$", fontsize=18)
ax.set_xlim(axes[0:2])
ax.set_ylim(axes[2:4])
ax.set_zlim(axes[4:6])

plt.show()

瑞士卷

简单地将数据集投射到一个平面上（例如，直接丢弃 x3 ）会将瑞士卷的不同层叠在一起，如下图左侧所示。但是，我们真正想要的是展开瑞士卷所获取到的类似下图右侧的 2D 数据集。

瑞士卷的投影与展开

流形学习 Manifold Learning

瑞士卷一个是二维流形的例子。简而言之，二维流形是一种二维形状，它可以在更高维空间中弯曲或扭曲。更一般地，一个 d 维流形是类似于 d 维超平面的 n 维空间（其中 d < n ）的一部分。在我们瑞士卷这个例子中， d = 2 ， n = 3 ：它有些像 2D 平面，但是它实际上是在第三维中卷曲。许多降维算法通过对训练实例所在的流形进行建模从而达到降维目的；这叫做流形学习。它依赖于流形猜想（manifold assumption），也被称为流形假设（manifold hypothesis），它认为大多数现实世界的高维数据集大都靠近一个更低维的流形。这种假设经常在实践中被证实。

3. 主成分分析 PCA

找到接近数据集分布的超平面，然后将所有的数据都投影到这个超平面上。

保留最大方差 variance

选择保持最大方差的轴，因为它比其他投影损失更少的信息。

子空间选择

左上图一个简单的二维数据集，以及三个不同的轴（即一维超平面）。图右边是将数据集投影到每个轴上的结果。投影到上保留了最大方差。

主成分

奇异值分解（SVD），可以将训练集矩阵 X 分解为三个矩阵 U·Σ·V^T 的点积，其中 V^T 包含我们想要的所有主成分

X_centered = X - X.mean(axis=0)
U, s, Vt = np.linalg.svd(X_centered)
c1 = Vt.T[:, 0]
c2 = Vt.T[:, 1]

投影到d维空间

过将数据集投影到由前 d 个主成分构成的超平面上。

W2 = Vt.T[:, :2]
X2D = X_centered.dot(W2)

SKlearn中的PCA见上文的例子

PCA压缩

在降维之后，训练集占用的空间要少得多。
将 PCA 应用于 MNIST 数据集，同时保留 95% 的方差。你应该发现每个实例只有 154 多个特征，而不是原来的 784 个特征。压缩比率不到20%。

pca = PCA(n_components=0.95)
X_reduced = pca.fit_transform(X_train)

pca.n_components_
# 154

pca = PCA(n_components = 154)
X_reduced = pca.fit_transform(X_train)
X_recovered = pca.inverse_transform(X_reduced)

def plot_digits(instances, images_per_row=5, **options):
    size = 28
    images_per_row = min(len(instances), images_per_row)
    images = [instance.reshape(size,size) for instance in instances]
    n_rows = (len(instances) - 1) // images_per_row + 1
    row_images = []
    n_empty = n_rows * images_per_row - len(instances)
    images.append(np.zeros((size, size * n_empty)))
    for row in range(n_rows):
        rimages = images[row * images_per_row : (row + 1) * images_per_row]
        row_images.append(np.concatenate(rimages, axis=1))
    image = np.concatenate(row_images, axis=0)
    plt.imshow(image, cmap = matplotlib.cm.binary, **options)
    plt.axis("off")

plt.figure(figsize=(7, 4))
plt.subplot(121)
plot_digits(X_train[::2100])
plt.title("Original", fontsize=16)
plt.subplot(122)
plot_digits(X_recovered[::2100])
plt.title("Compressed", fontsize=16)

PCA压缩

增量PCA

增量 PCA（IPCA）算法：您可以将训练集分批，并一次只对一个批量使用 IPCA 算法。这对大型训练集非常有用，并且可以在线应用 PCA（即在新实例到达时即时运行）。

from sklearn.decomposition import IncrementalPCA

n_batches = 100
inc_pca = IncrementalPCA(n_components=154)
for X_batch in np.array_split(X_train, n_batches):
    print(".", end="") # not shown in the book
    inc_pca.partial_fit(X_batch)

X_reduced = inc_pca.transform(X_train)

随机PCA

随机 PCA。这是一种随机算法，可以快速找到前d个主成分的近似值。它的计算复杂度是 O(m × d^2) + O(d^3) ，而不是 O(m × n^2) + O(n^3) ，所以当 d 远小于 n 时，它比之前的算法快得多。

在MNIST数据集上做测试

import time

for n_components in (2, 10, 154):
    print("n_components =", n_components)
    regular_pca = PCA(n_components=n_components)
    inc_pca = IncrementalPCA(n_components=n_components, batch_size=500)
    rnd_pca = PCA(n_components=n_components, random_state=42, svd_solver="randomized")

    for pca in (regular_pca, inc_pca, rnd_pca):
        t1 = time.time()
        pca.fit(X_train)
        t2 = time.time()
        print("    {}: {:.1f} seconds".format(pca.__class__.__name__, t2 - t1))

'''
n_components = 2
    PCA: 2.5 seconds
    IncrementalPCA: 57.0 seconds
    PCA: 2.3 seconds
n_components = 10
    PCA: 3.3 seconds
    IncrementalPCA: 58.8 seconds
    PCA: 3.1 seconds
n_components = 154
    PCA: 8.5 seconds
    IncrementalPCA: 87.4 seconds
    PCA: 9.6 seconds
'''

核PCA （Kernal PCA）

核技巧，一种将实例隐式映射到非常高维空间（称为特征空间）的数学技术，从而可以执行复杂的非线性投影来降低维度。

X, t = make_swiss_roll(n_samples=1000, noise=0.2, random_state=42)

from sklearn.decomposition import KernelPCA

from sklearn.decomposition import KernelPCA

lin_pca = KernelPCA(n_components = 2, kernel="linear", fit_inverse_transform=True)
rbf_pca = KernelPCA(n_components = 2, kernel="rbf", gamma=0.0433, fit_inverse_transform=True)
sig_pca = KernelPCA(n_components = 2, kernel="sigmoid", gamma=0.001, coef0=1, fit_inverse_transform=True)

y = t > 6.9

plt.figure(figsize=(11, 4))
for subplot, pca, title in ((131, lin_pca, "Linear kernel"), (132, rbf_pca, "RBF kernel, $\gamma=0.04$"), (133, sig_pca, "Sigmoid kernel, $\gamma=10^{-3}, r=1$")):
    X_reduced = pca.fit_transform(X)
    if subplot == 132:
        X_reduced_rbf = X_reduced
    
    plt.subplot(subplot)
    #plt.plot(X_reduced[y, 0], X_reduced[y, 1], "gs")
    #plt.plot(X_reduced[~y, 0], X_reduced[~y, 1], "y^")
    plt.title(title, fontsize=14)
    plt.scatter(X_reduced[:, 0], X_reduced[:, 1], c=t, cmap=plt.cm.hot)
    plt.xlabel("$z_1$", fontsize=18)
    if subplot == 131:
        plt.ylabel("$z_2$", fontsize=18, rotation=0)
    plt.grid(True)

plt.show()

不同核的效果

选择一种核并调整参数

from sklearn.model_selection import GridSearchCV
from sklearn.linear_model import LogisticRegression
from sklearn.pipeline import Pipeline

clf = 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)

print(grid_search.best_params_)
# {'kpca__gamma': 0.043333333333333335, 'kpca__kernel': 'rbf'}

4. LLE (Locally Linear Embedding) 局部线性嵌入

LLE是另一种非常有效的非线性降维方法。这是一种流形学习技术，不依赖于投影。
简而言之，LLE 首先测量每个训练实例与其最近邻之间的线性关系，然后寻找能最好地保留这些局部关系的训练集的低维表示 。这使得它特别擅长展开扭曲的流形，尤其是在没有太多噪音的情况下。

X, t = make_swiss_roll(n_samples=1000, noise=0.2, random_state=41)

from sklearn.manifold import LocallyLinearEmbedding

lle = LocallyLinearEmbedding(n_components=2, n_neighbors=10, random_state=42)
X_reduced = lle.fit_transform(X)

plt.title("Unrolled swiss roll using LLE", fontsize=14)
plt.scatter(X_reduced[:, 0], X_reduced[:, 1], c=t, cmap=plt.cm.hot)
plt.xlabel("$z_1$", fontsize=18)
plt.ylabel("$z_2$", fontsize=18)
plt.axis([-0.065, 0.055, -0.1, 0.12])
plt.grid(True)

plt.show()

LLE

5. 其他降维方法

• 多维缩放 MDS： 在尝试保持实例之间距离的同时降低了维度
• Isomap：通过将每个实例连接到最近的邻居来创建图形，然后在尝试保持实例之间的测地距离时降低维度。
• t-SNE（t-Distributed Stochastic Neighbor Embedding）：同时试图保持相似的实例临近并将不相似的实例分开。它主要用于可视化。

from sklearn.manifold import MDS

mds = MDS(n_components=2, random_state=42)
X_reduced_mds = mds.fit_transform(X)

from sklearn.manifold import Isomap

isomap = Isomap(n_components=2)
X_reduced_isomap = isomap.fit_transform(X)

from sklearn.manifold import TSNE

tsne = TSNE(n_components=2, random_state=42)
X_reduced_tsne = tsne.fit_transform(X)

titles = ["MDS", "Isomap", "t-SNE"]

plt.figure(figsize=(11,4))

for subplot, title, X_reduced in zip((131, 132, 133), titles,
                                     (X_reduced_mds, X_reduced_isomap, X_reduced_tsne)):
    plt.subplot(subplot)
    plt.title(title, fontsize=14)
    plt.scatter(X_reduced[:, 0], X_reduced[:, 1], c=t, cmap=plt.cm.hot)
    plt.xlabel("$z_1$", fontsize=18)
    if subplot == 131:
        plt.ylabel("$z_2$", fontsize=18, rotation=0)
    plt.grid(True)
plt.show()

其他降维方法