说明：本文依据《Sklearn 与 TensorFlow 机器学习实用指南》完成，所有版权和解释权均归作者和翻译成员所有，我只是搬运和做注解。

第八章降维

到达第一部分机器学习的最终章，降维，最早也是在这里开始应用的，当时是使用Sklearn中LDA模型完成主题抽取，现在又回到这里，另外可能看帖的童鞋也发现了，最近的状态有问题，一直在往前推，但是对代码的分析变得很少。我只是想尽快搞定，然后转到NLP上去。

很多机器学习的问题都会涉及到有着几千甚至数百万维的特征的训练实例。这不仅让训练过程变得非常缓慢，同时还很难找到一个很好的解，这种问题通常被称为维数灾难（curse of dimentionality）。

降维会让项目更复杂因而更难维护。所有应该先尝试使用原始的数据训练，如果训练速度太慢的话再考虑使用降维。在某些情况下，降低训练集数据的维度可能会筛选掉一些噪音和不必要的细节，这可能会让你的结果比降维之前更好。

源代码已经同步在github中

https://github.com/jwc19890114/-02-learning-file-100days

2.PCA主成分分析

主成分分析是目前最为流行的一种降维方法，先找到数据集分布的超平面，然后投影在这个平面上。

保留最大方差

在将训练集投影到较低维超平面之前，需要选择正确的超平面。下图左侧是一个简单的二维数据集，以及三个不同的轴（即一维超平面）。右侧图是将数据集投影到每个轴上的结果。正如你所看到的，投影到实线上保留了最大方差，而在点线上的投影只保留了非常小的方差，投影到虚线上保留的方差则处于上述两者之间。

准备工作

import numpy as np
import os
np.random.seed(42)
# 绘图设置
%matplotlib inline
import matplotlib as mpl
import matplotlib.pyplot as plt
mpl.rc('axes', labelsize=14)
mpl.rc('xtick', labelsize=12)
mpl.rc('ytick', labelsize=12)
# 保存图片路径
PROJECT_ROOT_DIR = "."
CHAPTER_ID = "unsupervised_learning"
def save_fig(fig_id, tight_layout=True):
 path = os.path.join(PROJECT_ROOT_DIR, "images", CHAPTER_ID, fig_id + ".png")
 print("Saving figure", fig_id)
 if tight_layout:
 plt.tight_layout()
 plt.savefig(path, format='png', dpi=300)
import warnings
warnings.filterwarnings(action="ignore", message="^internal gelsd")

主成分（Principle Componets）

PCA寻找训练集中可获得最大方差的轴。在2D的实例中，找到的就是点线。但如果在一个更高维的数据集中，PCA可以找到与前两个轴正交的第三个轴，以及与数据集中维数相同的第四个轴，第五个轴等。定义第i个轴的单位矢量被称为第i个主成分（PC）。

PCA使用SVD进行分解

注意，svd方法返回的U，s和Vt，其中Vt就是矩阵V的转置。

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

在这段代码中使用的是numpy提供的svd生成主成分，所以需要做一次数据集中心化处理。因为PCA假定数据集以原点为中心。之后使用的Sklearn的PCA类负责为您的数据集中心化处理。如果自己手撸PCA或者使用其他库，不要忘记首先要先对数据做中心化处理。

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)
#数据中心化处理
X_centered=X-X.mean(axis=0)
#使用SVD（奇异值分解）作为标准矩阵分解技术，将训练集矩阵X分解为三个矩阵U·Σ·V^T的点积，其中V^T包含我们想要的所有主成分
U,s,Vt=np.linalg.svd(X_centered)
c1=Vt.T[:,0]
c2=Vt.T[:,1]
m,n=X.shape
S=np.zeros(X_centered.shape)
S[:n,:n]=np.diag(s)

投影到D维空间

一旦确定了所有的主成分，就可以通过将数据集投影到由前d个主成分构成的超平面上，从而将数据集的维数降至d维。选择这个超平面可以确保投影将保留尽可能多的方差。计算投影的话就直接计算矩阵X和Wd的点积，其中Wd就是包含前面d个主成分矩阵，这里就是Vt的前两个元素Vt.T[:,2]。

到这里基本上就完成核心算法了，接下来是可视化操作，然而，代码还是有点难啊。

#将训练集投影到由前两个主成分定义的超平面上
W2=Vt.T[:,2]
X2D=X_centered.dot(W2)
X2D_using_svd=X2D
angle=np.pi/5
streth=5
m=200
np.random.seed(42)
X = np.random.randn(m, 2) / 10
X = X.dot(np.array([[stretch, 0],[0, 1]])) # stretch
X = X.dot([[np.cos(angle), np.sin(angle)], [-np.sin(angle), np.cos(angle)]]) # 旋转X
u1 = np.array([np.cos(angle), np.sin(angle)])
u2 = np.array([np.cos(angle - 2 * np.pi/6), np.sin(angle - 2 * np.pi/6)])
u3 = np.array([np.cos(angle - np.pi/2), np.sin(angle - np.pi/2)])
#绘制图像
X_proj1 = X.dot(u1.reshape(-1, 1))
X_proj2 = X.dot(u2.reshape(-1, 1))
X_proj3 = X.dot(u3.reshape(-1, 1))
plt.figure(figsize=(8,4))
plt.subplot2grid((3,2), (0, 0), rowspan=3)
plt.plot([-1.4, 1.4], [-1.4*u1[1]/u1[0], 1.4*u1[1]/u1[0]], "k-", linewidth=1)
plt.plot([-1.4, 1.4], [-1.4*u2[1]/u2[0], 1.4*u2[1]/u2[0]], "k--", linewidth=1)
plt.plot([-1.4, 1.4], [-1.4*u3[1]/u3[0], 1.4*u3[1]/u3[0]], "k:", linewidth=2)
plt.plot(X[:, 0], X[:, 1], "bo", alpha=0.5)
plt.axis([-1.4, 1.4, -1.4, 1.4])
plt.arrow(0, 0, u1[0], u1[1], head_width=0.1, linewidth=5, length_includes_head=True, head_length=0.1, fc='k', ec='k')
plt.arrow(0, 0, u3[0], u3[1], head_width=0.1, linewidth=5, length_includes_head=True, head_length=0.1, fc='k', ec='k')
plt.text(u1[0] + 0.1, u1[1] - 0.05, r"$\mathbf{c_1}$", fontsize=22)
plt.text(u3[0] + 0.1, u3[1], r"$\mathbf{c_2}$", fontsize=22)
plt.xlabel("$x_1$", fontsize=18)
plt.ylabel("$x_2$", fontsize=18, rotation=0)
plt.grid(True)
plt.subplot2grid((3,2), (0, 1))
plt.plot([-2, 2], [0, 0], "k-", linewidth=1)
plt.plot(X_proj1[:, 0], np.zeros(m), "bo", alpha=0.3)
plt.gca().get_yaxis().set_ticks([])
plt.gca().get_xaxis().set_ticklabels([])
plt.axis([-2, 2, -1, 1])
plt.grid(True)
plt.subplot2grid((3,2), (1, 1))
plt.plot([-2, 2], [0, 0], "k--", linewidth=1)
plt.plot(X_proj2[:, 0], np.zeros(m), "bo", alpha=0.3)
plt.gca().get_yaxis().set_ticks([])
plt.gca().get_xaxis().set_ticklabels([])
plt.axis([-2, 2, -1, 1])
plt.grid(True)
plt.subplot2grid((3,2), (2, 1))
plt.plot([-2, 2], [0, 0], "k:", linewidth=2)
plt.plot(X_proj3[:, 0], np.zeros(m), "bo", alpha=0.3)
plt.gca().get_yaxis().set_ticks([])
plt.axis([-2, 2, -1, 1])
plt.xlabel("$z_1$", fontsize=18)
plt.grid(True)
plt.show()

从图中可以看到，选择保持最大方差的轴看起来是最合理的，因为其投影损失更少的信息。证明这种选择的另一种方法是，选择这个轴使得将原始数据集投影到该轴上的均方距离最小。这是就 PCA 背后的思想，相当简单。

使用Sklearn完成PCA降维

引入mnist数据，这里要对mnist数据进行降维处理。

from sklearn.decomposition import PCA
from sklearn.datasets import fetch_mldata
from sklearn.model_selection import train_test_split
mnist = fetch_mldata('MNIST original',data_home="./MNIST_data")
X=mnist["data"]
y=mnist["target"]
X_train,X_test,y_train,y_test=train_test_split(X,y)
pca=PCA()
pca.fit(X_train,y_train)
#cumsum(n)：实现n轴上的累加：以最外面的数组元素为单位，以[[1,2,3],[8,9,12]]为开始实现后面元素的对应累加
#pca.explained_variance_ratio_：返回的是所保留的n个成分各自的方差百分比,这里可以理解为单个变量方差贡献率
cumsum=np.cumsum(pca.explained_variance_ratio_)
print(pca.explained_variance_ratio_)

选择正确的维度

教程中提到倾向于选择加起来到方差解释率能够达到足够占比（95%）的维度的数量，而不是任意选择要降低到的维度数量。

下面的代码在不降维的情况下进行 PCA，然后计算出保留训练集方差 95% 所需的最小维数。

d=np.argmax(cumsum>=0.95)+1
print(d) #=>154

可以将这个d=154作为要降的维度数再次带入PCA模型中进行训练，同时教程中提到另一个方法，不指定你想要保留的主成分个数，而是将n_components设置为 0.0 到 1.0 之间的浮点数，表明希望保留的方差比率，这里使用n_components=0.95。可以发现现在的理想降维度数已经由28*28=784降为154。

pca = PCA(n_components=0.95)
X_reduced = pca.fit_transform(X_train)
print(pca.n_components_) #=>154

PCA压缩

np.sum(pca.explained_variance_ratio_)
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 = mpl.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)