经典数据降维算法

hyl_2022

已于 2024-03-05 15:21:17 修改

阅读量1k

点赞数 14

文章标签：算法机器学习人工智能

于 2024-02-22 20:15:24 首次发布

本文链接：https://blog.csdn.net/weixin_45871964/article/details/136242052

版权

本文介绍了数据降维的基本概念，包括其优点和原理，并详细讲解了主成分分析（PCA）、核主成分分析（KPCA）、线性判别分析（LDA）和多维标度分析（MDS）的原理及Python代码示例。这些技术用于减少高维数据的复杂性，提升数据处理效率和可视化效果。

摘要由CSDN通过智能技术生成

数据降维

定义

将高维数据下降为低维数据。

优点

使得数据集更易使用
确保变量之间彼此独立
降低算法计算运算成本
去除噪音

原理

删除冗余的数据、无效信息、重复表达内容等。例如，对于对称图形而言，对称部分的信息则可以归为重复信息。大部分经典降维技术也是基于这一内容而展开，其中降维方法又分为线性和非线性降维，非线性降维又分为基于核函数和基于特征值的方法。

线性降维方法：

PCA 、ICA 、LDA、LFA、LPP(LE 的线性表示)
非线性降维方法：

基于核函数的非线性降维方法——KPCA 、KICA、KDA

基于特征值的非线性降维方法（流型学习）——ISOMAP、LLE、LE、LPP、LTSA、MVU

主成分分析（PCA，principal components analysis）降维算法

原理

PCA是一种简单的无监督降维算法，旨在通过将数据投影到变化最大或重构误差最小的方向来减少维度。它是由Karl Pearson于1901年提出的线性降维方法。PCA相关的原理通常被称为最大方差理论或最小误差理论，目标相同但过程略有不同。

将一组 N 维向量降维到 K 维（K 大于 0，小于 N），目标是选择 K 个单位正交基，确保各字段两两间 COV(X,Y) 为 0，同时让字段的方差尽可能大。所以，最大方差就是让投影数据的方差最大化。在这个过程中，我们需要找到数据集 Xmxn 的最佳投影空间 Wnxk，协方差矩阵等。算法流程如下：

输入：数据集 Xmxn；
计算数据集 X 的均值 Xmean，并让 Xnew=X−Xmean；
求解矩阵 Xnew 的协方差矩阵 Cov；
计算协方差矩阵 Cov 的特征值和对应的特征向量；
将特征值从大到小排序，选择最大的 k 个，然后将对应的 k 个特征向量作为列向量组成特征向量矩阵 Wnxk；
计算 XnewW，即将数据集 Xnew 投影到选取的特征向量上，得到我们需要的降维数据集 XnewW。

代码


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
%matplotlib inline

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()

KPCA（kernel PCA）

原理

KPCA 是核技术与 PCA 结合的产物，它与 PCA 主要差别在于计算协方差矩阵时使用了核函数，即是经过核函数映射之后的协方差矩阵。

引入核函数可以很好的解决非线性数据映射问题。kPCA 可以将非线性数据映射到高维空间，在高维空间下使用标准 PCA 将其映射到另一个低维空间。

代码

# coding:utf-8
# 实现KPCA

from sklearn.datasets import load_iris
from sklearn.decomposition import KernelPCA
import numpy as np
import matplotlib.pyplot as plt
from scipy.spatial.distance import pdist, squareform

def sigmoid(x, coef = 0.25):
    x = np.dot(x, x.T)
    return np.tanh(coef*x+1)

def linear(x):
    x = np.dot(x, x.T)
    return x

def rbf(x, gamma = 15):
    sq_dists = pdist(x, 'sqeuclidean')
    mat_sq_dists = squareform(sq_dists)
    return np.exp(-gamma*mat_sq_dists)

def kpca(data, n_dims=2, kernel = rbf):
    '''

    :param data: (n_samples, n_features)
    :param n_dims: target n_dims
    :param kernel: kernel functions
    :return: (n_samples, n_dims)
    '''

    K = kernel(data)
    #
    N = K.shape[0]
    one_n = np.ones((N, N)) / N
    K = K - one_n.dot(K) - K.dot(one_n) + one_n.dot(K).dot(one_n)
    #
    eig_values, eig_vector = np.linalg.eig(K)
    idx = eig_values.argsort()[::-1]
    eigval = eig_values[idx][:n_dims]
    eigvector = eig_vector[:, idx][:, :n_dims]
    print(eigval)
    eigval = eigval**(1/2)
    vi = eigvector/eigval.reshape(-1,n_dims)
    data_n = np.dot(K, vi)
    return data_n


if __name__ == "__main__":
    data = load_iris().data
    Y = load_iris().target
    data_1 = kpca(data, kernel=rbf)


    sklearn_kpca = KernelPCA(n_components=2, kernel="rbf", gamma=15)
    data_2 = sklearn_kpca.fit_transform(data)

    plt.figure(figsize=(8,4))
    plt.subplot(121)
    plt.title("my_KPCA")
    plt.scatter(data_1[:, 0], data_1[:, 1], c = Y)

    plt.subplot(122)
    plt.title("sklearn_KPCA")
    plt.scatter(data_2[:, 0], data_2[:, 1], c = Y)
    plt.show()

LDA（Linear Discriminant Analysis）

原理

LDA 是一种可作为特征抽取的技术，其目标是向最大化类间差异，最小化类内差异的方向投影，以利于分类等任务即将不同类的样本有效的分开。LDA 可以提高数据分析过程中的计算效率，对于未能正则化的模型，可以降低维度灾难带来的过拟合。
最小误差理论降维原理

代码

#coding:utf-8
import numpy as np
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis
from sklearn.datasets import load_iris
import matplotlib.pyplot as plt



def lda(data, target, n_dim):
    '''
    :param data: (n_samples, n_features)
    :param target: data class
    :param n_dim: target dimension
    :return: (n_samples, n_dims)
    '''

    clusters = np.unique(target)

    if n_dim > len(clusters)-1:
        print("K is too much")
        print("please input again")
        exit(0)

    #within_class scatter matrix
    Sw = np.zeros((data.shape[1],data.shape[1]))
    for i in clusters:
        datai = data[target == i]
        datai = datai-datai.mean(0)
        Swi = np.mat(datai).T*np.mat(datai)
        Sw += Swi

    #between_class scatter matrix
    SB = np.zeros((data.shape[1],data.shape[1]))
    u = data.mean(0)  #所有样本的平均值
    for i in clusters:
        Ni = data[target == i].shape[0]
        ui = data[target == i].mean(0)  #某个类别的平均值
        SBi = Ni*np.mat(ui - u).T*np.mat(ui - u)
        SB += SBi
    S = np.linalg.inv(Sw)*SB
    eigVals,eigVects = np.linalg.eig(S)  #求特征值，特征向量
    eigValInd = np.argsort(eigVals)
    eigValInd = eigValInd[:(-n_dim-1):-1]
    w = eigVects[:,eigValInd]
    data_ndim = np.dot(data, w)

    return data_ndim

if __name__ == '__main__':
    iris = load_iris()
    X = iris.data
    Y = iris.target
    data_1 = lda(X, Y, 2)

    data_2 = LinearDiscriminantAnalysis(n_components=2).fit_transform(X, Y)


    plt.figure(figsize=(8,4))
    plt.subplot(121)
    plt.title("my_LDA")
    plt.scatter(data_1[:, 0], data_1[:, 1], c = Y)

    plt.subplot(122)
    plt.title("sklearn_LDA")
    plt.scatter(data_2[:, 0], data_2[:, 1], c = Y)
    plt.savefig("LDA.png")
    plt.show()

MDS（multidimensional scaling）

原理

MDS 即多维标度分析，它是一种通过直观空间图表示研究对象的感知和偏好的传统降维方法。该方法会计算任意两个样本点之间的距离，使得投影到低维空间之后能够保持这种相对距离从而实现投影。

由于 sklearn 中 MDS 是采用迭代优化方式，下面实现了迭代和非迭代的两种。

代码

# coding:utf-8
import numpy as np
from sklearn.datasets import load_iris
from sklearn.manifold import MDS
import matplotlib.pyplot as plt

def cal_pairwise_dist(x):
    '''计算pairwise 距离, x是matrix
    (a-b)^2 = a^2 + b^2 - 2*a*b
    '''
    sum_x = np.sum(np.square(x), 1)
    dist = np.add(np.add(-2 * np.dot(x, x.T), sum_x).T, sum_x)
    #返回任意两个点之间距离的平方
    return dist


def my_mds(data, n_dims):
    '''

    :param data: (n_samples, n_features)
    :param n_dims: target n_dims
    :return: (n_samples, n_dims)
    '''

    n, d = data.shape
    dist = cal_pairwise_dist(data)
    dist[dist < 0 ] = 0
    T1 = np.ones((n,n))*np.sum(dist)/n**2
    T2 = np.sum(dist, axis = 1, keepdims=True)/n
    T3 = np.sum(dist, axis = 0, keepdims=True)/n

    B = -(T1 - T2 - T3 + dist)/2

    eig_val, eig_vector = np.linalg.eig(B)
    index_ = np.argsort(-eig_val)[:n_dims]
    picked_eig_val = eig_val[index_].real
    picked_eig_vector = eig_vector[:, index_]
    # print(picked_eig_vector.shape, picked_eig_val.shape)
    return picked_eig_vector*picked_eig_val**(0.5)

if __name__ == '__main__':
    iris = load_iris()
    data = iris.data
    Y = iris.target
    data_1 = my_mds(data, 2)

    data_2 = MDS(n_components=2).fit_transform(data)

    plt.figure(figsize=(8, 4))
    plt.subplot(121)
    plt.title("my_MDS")
    plt.scatter(data_1[:, 0], data_1[:, 1], c=Y)

    plt.subplot(122)
    plt.title("sklearn_MDS")
    plt.scatter(data_2[:, 0], data_2[:, 1], c=Y)
    plt.savefig("MDS_1.png")
    plt.show()