机器学习----降低维度（Dimensionality Reduction）算法原理及python实现

通常情况下，在收集数据集时会有很多的特征，这代表着数据是高冗余的表示，但是对于某个工程来说其实可能并不需要那么多的特征。所以就需要给数据进行降维（Dimensionality Reduction）。

• 降维可以简化数据，使数据集更易使用，降低时间开销，而且能减少一部分噪音的影响，使最后的效果变好。

目前，主要降维方法有：主成分分析法(Principal Component Analysis, PCA)、因子分析法(Factor Analysis)、独立成分分析(Independent Component Analysis, ICA)

在这里我们主要实现主成分分析法(Principal Component Analysis, PCA)：

主成分分析是最常用的一种降维方法，不同于MDS采用距离保持的方法，主成分分析（PCA）直接通过一个线性变换，将原始空间中的样本投影到新的低维空间中。

• 简单来理解这一过程便是：PCA采用一组新的基来表示样本点，其中每一个基向量都是原来基向量的线性组合，通过使用尽可能少的新基向量来表出样本，从而达到降维的目的。

假设使用d’个新基向量来表示原来样本，实质上是将样本投影到一个由d’个基向量确定的一个超平面上（即舍弃了一些维度），要用一个超平面对空间中所有高维样本进行恰当的表达，最理想的情形是：若这些样本点都能在超平面上表出且这些表出在超平面上都能够很好地分散开来。

要求这个超平面大概具有如下两个性质：
最近重构性：样本点到超平面的距离足够近，即尽可能在超平面附近；
最大可分性：样本点在超平面上的投影尽可能地分散开来，即投影后的坐标具有区分性。

• PCA算法过程：
  输入为样本集D={x1,x2,…,xm}以及低维空间维数d’；输出为投影矩阵W*=（w1,w2,…,wd’）。
  1.对所有样本进行中心化：找到样本集的中点，将样本集整体移动直到中点移到原点位置；
  2.计算样本的协方差矩阵XXT；
  3.对协方差矩阵XXT做特征值分解；
  4.取最大的d’个特征值所对应的特征向量w1,w2,…,wd’。

PCA 仅需保留W*与样本的均值向量即可通过简单的向量减法和矩阵-向量乘法将新样本投影至低维空间中。显然，低维空间与原始高维空间必有不同， 因为对应于最小的d-d’个特征值的特征向量被舍弃了，这是降维导致的结果. 但舍弃这部分信息往往是必要的：
1.一方面舍弃这部分信息之后能使样本的采样密度增大，这正是降维的重要动机；
2.另一方面，当数据受到噪声影响时，最小的特征值所对应的特征向量往往与噪声有关，将它们舍弃能在一定程度上起到去噪的效果。

PCA的python实现代码为（参考大佬）：

from __future__ import print_function
from sklearn import datasets
import matplotlib.pyplot as plt
import matplotlib.cm as cmx
import matplotlib.colors as colors
import numpy as np




def shuffle_data(X, y, seed=None):
    if seed:
        np.random.seed(seed)

    idx = np.arange(X.shape[0])
    np.random.shuffle(idx)

    return X[idx], y[idx]



# 正规化数据集 X
def normalize(X, axis=-1, p=2):
    lp_norm = np.atleast_1d(np.linalg.norm(X, p, axis))
    lp_norm[lp_norm == 0] = 1
    return X / np.expand_dims(lp_norm, axis)


# 标准化数据集 X
def standardize(X):
    X_std = np.zeros(X.shape)
    mean = X.mean(axis=0)
    std = X.std(axis=0)

    # 做除法运算时请永远记住分母不能等于0的情形
    # X_std = (X - X.mean(axis=0)) / X.std(axis=0)
    for col in range(np.shape(X)[1]):
        if std[col]:
            X_std[:, col] = (X_std[:, col] - mean[col]) / std[col]

    return X_std


# 划分数据集为训练集和测试集
def train_test_split(X, y, test_size=0.2, shuffle=True, seed=None):
    if shuffle:
        X, y = shuffle_data(X, y, seed)

    n_train_samples = int(X.shape[0] * (1-test_size))
    x_train, x_test = X[:n_train_samples], X[n_train_samples:]
    y_train, y_test = y[:n_train_samples], y[n_train_samples:]

    return x_train, x_test, y_train, y_test



# 计算矩阵X的协方差矩阵
def calculate_covariance_matrix(X, Y=np.empty((0,0))):
    if not Y.any():
        Y = X
    n_samples = np.shape(X)[0]
    covariance_matrix = (1 / (n_samples-1)) * (X - X.mean(axis=0)).T.dot(Y - Y.mean(axis=0))

    return np.array(covariance_matrix, dtype=float)


# 计算数据集X每列的方差
def calculate_variance(X):
    n_samples = np.shape(X)[0]
    variance = (1 / n_samples) * np.diag((X - X.mean(axis=0)).T.dot(X - X.mean(axis=0)))
    return variance


# 计算数据集X每列的标准差
def calculate_std_dev(X):
    std_dev = np.sqrt(calculate_variance(X))
    return std_dev


# 计算相关系数矩阵
def calculate_correlation_matrix(X, Y=np.empty([0])):
    # 先计算协方差矩阵
    covariance_matrix = calculate_covariance_matrix(X, Y)
    # 计算X, Y的标准差
    std_dev_X = np.expand_dims(calculate_std_dev(X), 1)
    std_dev_y = np.expand_dims(calculate_std_dev(Y), 1)
    correlation_matrix = np.divide(covariance_matrix, std_dev_X.dot(std_dev_y.T))

    return np.array(correlation_matrix, dtype=float)



class PCA():
    \"\"\"
    主成份分析算法PCA，非监督学习算法.
    \"\"\"
    def __init__(self):
        self.eigen_values = None
        self.eigen_vectors = None
        self.k = 2

    def transform(self, X):
        \"\"\"
        将原始数据集X通过PCA进行降维
        \"\"\"
        covariance = calculate_covariance_matrix(X)

        # 求解特征值和特征向量
        self.eigen_values, self.eigen_vectors = np.linalg.eig(covariance)

        # 将特征值从大到小进行排序，注意特征向量是按列排的，即self.eigen_vectors第k列是self.eigen_values中第k个特征值对应的特征向量
        idx = self.eigen_values.argsort()[::-1]
        eigenvalues = self.eigen_values[idx][:self.k]
        eigenvectors = self.eigen_vectors[:, idx][:, :self.k]

        # 将原始数据集X映射到低维空间
        X_transformed = X.dot(eigenvectors)

        return X_transformed


def main():
    # Load the dataset
    data = datasets.load_iris()
    X = data.data
    y = data.target

    # 将数据集X映射到低维空间
    X_trans = PCA().transform(X)

    x1 = X_trans[:, 0]
    x2 = X_trans[:, 1]

    cmap = plt.get_cmap(\'viridis\')
    colors = [cmap(i) for i in np.linspace(0, 1, len(np.unique(y)))]

    class_distr = []
    # Plot the different class distributions
    for i, l in enumerate(np.unique(y)):
        _x1 = x1[y == l]
        _x2 = x2[y == l]
        _y = y[y == l]
        class_distr.append(plt.scatter(_x1, _x2, color=colors[i]))

    # Add a legend
    plt.legend(class_distr, y, loc=1)

    # Axis labels
    plt.xlabel(\'Principal Component 1\')
    plt.ylabel(\'Principal Component 2\')
    plt.show()


if __name__ == \"__main__\":
    main()

运行结果为：