PCA（主成分分析）是什么

降维算法：保留原始数据的关键信息
注意：在降维过程中，维度通常是指特征的数量，而不是空间的维度。

主成分分析（PCA）是一个常见的特征提取方法，它通过找到数据中的主成分（方差最大的方向）来实现降维。
通俗理解，就是寻找新的坐标系，使得数据尽可能分布在一个或几个坐标轴上，同时尽可能保留原先数据分布的主要信息，使原先高维度的信息，在转换后能用低维度的信息来保存。而新坐标系的坐标轴，称为主成分(Principal components, PC)，这也就是PCA的名称来源。
为什么找方差最大的方向？
在信号处理领域中认为信号具有较大的方差，噪声有较小的方差；信噪比就是信号与噪声的方差比，越大意味着数据的质量越好。

PCA（主成分分析）用来干什么的，好处

通过保留主要特征来捕捉数据的重要结构，使得模型更快收敛。
展开来说，PCA去除噪声和不重要的特征，将多个指标转换为少数几个主成分，这些主成分是原始变量的线性组合，且彼此之间互不相关，其能反映出原始数据的大部分信息，而且可以提升数据处理的速度。

PCA（主成分分析）怎么实现的（三步）

PCA方法的总结
标准化（去中心化）；
从协方差矩阵或相关矩阵中获取特征向量和特征值，或执行奇异值分解SVD；
按降序对特征值进行排序，并选择$k$个与最大特征值相对应的特征向量，其中$k$是新特征子空间的维数（$k\le p$）；
将选择的$k$个特征向量构造为投影矩阵$W$；
通过$W$对原始数据集$X$进行变换，得到$k$维特征子空间$Y$。

准备数据集

关于Iris数据集
scikit-learn中的鸢尾花数据集（Iris dataset），一共包含150张图片。
3种鸢尾花：山鸢尾（50张）、杂色鸢尾（50张）、维吉尼亚鸢尾（50张）。
4个特征：萼片长度（厘米）、萼片宽度（厘米）、花瓣长度（厘米）、花瓣宽度（厘米）。

加载数据集
加载数据集

import pandas as pd

# load-dataset
df = pd.read_csv(
    filepath_or_buffer='https://archive.ics.uci.edu/ml/machine-learning-databases/iris/iris.data',
    header=None,
    sep=',')
df.columns = ['sepal_len', 'sepal_wid', 'petal_len', 'petal_wid', 'class']
df.dropna(how="all", inplace=True)  # drops the empty line at file-end

df.tail()

运行结果

划分成 数据$X$、标签$y$

# split data table into data X and class labels y
X = df.iloc[:, 0:4].values
y = df.iloc[:, 4].values

数据集为$150\times 4$的矩阵，不同列代表不同特征，每行为1个样本，每个样本行x可以被描述为一个4维向量。

可视化数据分布

from matplotlib import pyplot as plt
import seaborn as sns

# 设置 seaborn 的样式为 whitegrid
sns.set_style("whitegrid")

label_dict = {
    1: 'Iris-Setosa',
              2: 'Iris-Versicolor',
              3: 'Iris-Virgnica'}

feature_dict = {
    0: 'sepal length [cm]',
                1: 'sepal width [cm]',
                2: 'petal length [cm]',
                3: 'petal width [cm]'}

plt.figure(figsize=(8, 6))
for cnt in range(4):
    plt.subplot(2, 2, cnt + 1)
    for lab in ('Iris-setosa', 'Iris-versicolor', 'Iris-virginica'):
        plt.hist(X[y == lab, cnt],
                 label=lab,
                 bins=10,
                 alpha=0.3, )
    plt.xlabel(feature_dict[cnt])
plt.legend(loc='upper right', fancybox=True, fontsize=8)

plt.tight_layout()
plt.show()

运行结果

标准化

归一化：将数据映射到[0,1]或[-1,1]区间范围内，不同特征的量纲不同，值范围大小不同，存在奇异值，对训练有影响。

标准化：将数据映射到满足标准正态分布的范围内，使数据满足均值为0，标准差为1。

去中心化：使数据满足均值为0，但对标准差没有要求。

每种方法对应的使用场景
1.若对数据的范围没有限定要求，则选择标准化进行数据预处理
2.若要求数据在某个范围内取值，则采用归一化。
3若数据不存在极端的极大极小值时，采用归一化。
4.若数据存在较多的异常值和噪音，采用标准化。

标准化过程如下：
1、计算每个特征（每一列，共p个特征）的均值$\overline{x_i}$j和标准差$S_j$，公式如下：$$\overline{x_i}=\frac{1}{n}\sum _{i=1}^n{x}_{ij}$$ $$S_j=\sqrt{\frac{{\sum }_{i=1}^n(x_{ij}-\overline{x_j}{)}^2}{n-1}}$$
​
2、将每个样本的每个特征进行标准化处理，得到标准化特征矩阵$X_{stand}$：
$$X_{ij}=\frac{x_{ij}-\overline{x_j}}{S_j}$$

$$X_{stand}=\left[\begin{array}{cccc}{X}_{11}& X_{12}& \cdots & X_{1p}\\ X_{21}& X_{22}& \cdots & X_{2p}\\ ⋮& ⋮& \ddots & ⋮\\ X_{n1}& X_{n2}& \cdots & X_{np}\end{array}\right]=\left[{\mathbf{X}}_1,{\mathbf{X}}_2,\cdots ,{\mathbf{X}}_p\right]$$

标准化是沿着每个特征（即每列）独立进行的。换句话说，对于 X 中的每一列（即每一个特征），StandardScaler 会计算该列的均值和标准差，然后将该列中的每个值减去均值并除以标准差。

将数据单位厘米，转成单位尺度，且满足均值为0，标准差为1。

# 标准化(均值为0，标准差为1)
from sklearn.preprocessing import StandardScaler
X_std = StandardScaler().fit_transform(X)

第一步：特征分解——计算特征向量和特征值

第一种：特征分解

协方差矩阵（或相关矩阵）的特征向量和特征值代表PCA的“核心”：特征向量描述了数据的主要变化方向，而特征值则表示了在这些方向上的方差大小。

协方差矩阵

共p个特征，则计算得到的协方差矩阵为$p\times p$的方阵。
协方差矩阵计算公式
协方差矩阵$R$：
$$R=\left[\begin{array}{cccc}{r}_{11}& r_{12}& \cdots & r_{1p}\\ r_{21}& r_{22}& \cdots & r_{2p}\\ ⋮& ⋮& \ddots & ⋮\\ r_{p1}& r_{p2}& \cdots & r_{pp}\end{array}\right]$$
两个特征 j、k 之间的协方差计算公式：
$$r_{jk}=\frac{1}{n-1}\sum _{i=1}^n\left(x_{ij}-{\bar{x}}_j\right)\left(x_{ik}-{\bar{x}}_k\right).$$
$$=\frac{1}{n-1}\left((\mathbf{X}-\bar{\mathbf{x}}{)}^T(\mathbf{X}-\bar{\mathbf{x}})\right)$$
其中，$\bar{\mathbf{x}}$ 是均值向量
$$\bar{\mathbf{x}}=\frac{1}{n}\sum _{i=1}^n{x}_i.$$
均值向量是⼀个p 维向量，其中该向量中的每个值代表数据集中特征列的样本均值。

计算协方差矩阵

# 计算协方差矩阵
import numpy as np
mean_vec = np.mean(X_std, axis=0)
cov_mat = (X_std - mean_vec).T.dot((X_std - mean_vec)) / (X_std.shape[0]-1)
print('Covariance matrix \n%s' %cov_mat)

等价于np.cov

print('NumPy covariance matrix: \n%s' %np.cov(X_std.T))

对协方差矩阵进行特征分解——计算特征向量、特征值

# 对协方差矩阵进行特征分解
cov_mat = np.cov(X_std.T)  # 在标准化后的数据上，计算协方差矩阵

eig_vals, eig_vecs = np.linalg.eig(cov_mat)

print('Eigenvectors \n%s' %eig_vecs)  # 特征向量
print('\nEigenvalues \n%s' %eig_vals)  # 特征值

相关矩阵

特别是在“金融”领域，通常使用相关矩阵来代替协方差矩阵。但是，协方差矩阵的特征分解（如果输入数据已标准化）会产生与相关矩阵的特征分解相同的结果，因为相关矩阵可以理解为归一化的协方差矩阵。

对相关矩阵进行特征分解——计算特征向量、特征值
对标准化后的数据

# 对标准化后的数据计算相关矩阵并进行特征分解
cor_mat1 = np.corrcoef(X_std.T)  # 对标准化后的数据计算相关矩阵

eig_vals, eig_vecs = np.linalg.eig(cor_mat1)

print('Eigenvectors \n%s' %eig_vecs)  # 特征向量
print('\nEigenvalues \n%s' %eig_vals)  # 特征值

对原始数据

# 对原始数据计算相关矩阵并进行特征分解
cor_mat2 = np.corrcoef(X.T)  # 对原始数据计算相关矩阵

eig_vals, eig_vecs = np.linalg.eig(cor_mat2)

print('Eigenvectors \n%s' %eig_vecs)  # 特征向量
print('\nEigenvalues \n%s' %eig_vals)  # 特征值

总结

以上三种方法的特征分解结果相同，即计算出的特征向量和特征值都一样：
数据标准化后协方差矩阵的特征分解；
数据标准化后相关矩阵的特征分解；
原始数据相关矩阵的特征分解。

第二种：奇异值分解SVD

PCA通过特征分解找到使方差最大化的正交基，实现数据的降维和特征提取。
然而，特征分解只适用于方阵（协方差矩阵为方阵），而现实中的数据矩阵往往不是方阵。这时，奇异值分解（SVD）就派上了用场。SVD是一种能适用于任意矩阵的分解方法，它并不要求分解的矩阵为方阵。

得到的左奇异矩阵$u$是特征分解的特征向量；得到的奇异值$s$与特征分解的特征值类似表示方差的大小。

# 奇异值分解SVD
u,s,v = np.linalg.svd(X_std.T)
u,s

第二步：选择主成分

对特征进行排序

按降序对特征值进行排序，并选择$k$个与最大特征值相对应的特征向量，其中$k$是新特征子空间的维数（$k\le p$）

# Make a list of (eigenvalue, eigenvector) tuples
eig_pairs = [(np.abs(eig_vals[i]), eig_vecs[:,i]) for i in range(len(eig_vals))]

# Sort the (eigenvalue, eigenvector) tuples from high to low
eig_pairs.sort(key=lambda x: x[0], reverse=True)

# Visually confirm that the list is correctly sorted by decreasing eigenvalues
print('Eigenvalues in descending order:')
for i in eig_pairs:
    print(i[0])

为新特征子空间选择多少个主成分（k值取多少）合适？
计算主成分贡献率及累计贡献率
第$i$个主成分的贡献率，${\lambda }_i$代表第$i$个特征向量的特征值：
$$\frac{{\lambda }_i}{{\sum }_{k=1}^p{\lambda }_k}$$
前$i$个主成分的累计贡献率：
$$\frac{{\sum }_{j=1}^i{\lambda }_j}{{\sum }_{k=1}^p{\lambda }_k}$$
一般取累计贡献率超过80%的特征值所对应的第一、第二、…、第m（m ≤ p）个主成分

tot = sum(eig_vals)
var_exp = [(i / tot)*100 for i in sorted(eig_vals, reverse=True)]
cum_var_exp = np.cumsum(var_exp)

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

plt.bar(range(4), var_exp, alpha=0.5, align='center',
        label='individual explained variance')
plt.step(range(4), cum_var_exp, where='mid',
            label='cumulative explained variance')
plt.ylabel('Explained variance ratio')
plt.xlabel('Principal components')
plt.legend(loc='best')
plt.tight_layout()

上图清楚地表明，大部分方差（准确地说是方差的72.77%）可以单独由第一个主成分来解释。第二个主成分仍然包含一些信息（23.03%），而第三和第四主成分可以安全地删除而不会丢失太多信息。前两个主成分合计包含95.8%的信息。

构建投影矩阵

将选择的$k$个特征向量连接为投影矩阵$W$，该矩阵用于将数据映射到新的特征子空间。
在这个例子里，选择具有最高特征值的前2个特征向量来构建我们的特征向量，将 4 维特征空间简化为 2 维特征子空间 $p\times k$维投影矩阵$W$。

matrix_w = np.hstack((eig_pairs[0][1].reshape(4,1),
                      eig_pairs[1][1].reshape(4,1)))

print('Matrix W:\n', matrix_w)

第三步：投影到新特征空间

通过投影矩阵$W$对原始数据集$X$进行变换，得到$k$维特征子空间$Y$。
在最后⼀步中，我们将使用4 × 2维投影矩阵$W$通过方程$Y=X\times W$将样本转换到新的子空间，其中$Y$是$150\times 2$的矩阵。

Y = X_std.dot(matrix_w)

plt.figure(figsize=(6, 4))
for lab, col in zip(('Iris-setosa', 'Iris-versicolor', 'Iris-virginica'),
                    ('blue', 'red', 'green')):
    plt.scatter(Y[y==lab, 0],
                Y[y==lab, 1],
                label=lab,
                c=col)
plt.xlabel('Principal Component 1')
plt.ylabel('Principal Component 2')
plt.legend(loc='lower center')
plt.tight_layout()
plt.show()

scikit-learn 库中的 PCA

from sklearn.decomposition import PCA as sklearnPCA
sklearn_pca = sklearnPCA(n_components=2)
Y_sklearn = sklearn_pca.fit_transform(X_std)

plt.figure(figsize=(6, 4))
for lab, col in zip(('Iris-setosa', 'Iris-versicolor', 'Iris-virginica'),
                    ('blue', 'red', 'green')):
    plt.scatter(Y_sklearn[y==lab, 0],
                Y_sklearn[y==lab, 1],
                label=lab,
                c=col)
plt.xlabel('Principal Component 1')
plt.ylabel('Principal Component 2')
plt.legend(loc='lower center')
plt.tight_layout()
plt.show()