锐客网

  1. 首页 > 睿知 > it技术 > >

Python实现12种降维算法的示例代码

目录

  • 为什么要进行数据降维
  • 数据降维原理
  • 主成分分析(PCA)降维算法
  • 其它降维算法及代码地址
    • 1.KPCA(kernel PCA)
    • 2.LDA(Linear Discriminant Analysis)
    • 3.MDS(multidimensional scaling)
    • 4.ISOMAP
    • 5.LLE(locally linear embedding)
    • 6.t-SNE
    • 7.LE(Laplacian Eigenmaps)
    • 8.LPP(Locality Preserving Projections)
网上关于各种降维算法的资料参差不齐，同时大部分不提供源代码。这里有个 GitHub 项目整理了使用 Python 实现了 11 种经典的数据抽取(数据降维)算法，包括：PCA、LDA、MDS、LLE、TSNE 等，并附有相关资料、展示效果; 非常适合机器学习初学者和刚刚入坑数据挖掘的小伙伴。

为什么要进行数据降维 所谓降维，即用一组个数为 d 的向量 Zi 来代表个数为 D 的向量 Xi 所包含的有用信息，其中 d
通常，我们会发现大部分数据集的维度都会高达成百乃至上千，而经典的 MNIST，其维度都是 64。
Python实现12种降维算法的示例代码
文章图片

MNIST 手写数字数据集
但在实际应用中，我们所用到的有用信息却并不需要那么高的维度，而且每增加一维所需的样本个数呈指数级增长，这可能会直接带来极大的「维数灾难」; 而数据降维就可以实现：
  • 使得数据集更易使用
  • 确保变量之间彼此独立
  • 降低算法计算运算成本
去除噪音一旦我们能够正确处理这些信息，正确有效地进行降维，这将大大有助于减少计算量，进而提高机器运作效率。而数据降维，也常应用于文本处理、人脸识别、图片识别、自然语言处理等领域。

数据降维原理 往往高维空间的数据会出现分布稀疏的情况，所以在降维处理的过程中，我们通常会做一些数据删减，这些数据包括了冗余的数据、无效信息、重复表达内容等。
例如：现有一张 1024*1024 的图，除去中心 50*50 的区域其它位置均为零值，这些为零的信息就可以归为无用信息; 而对于对称图形而言，对称部分的信息则可以归为重复信息。
因此，大部分经典降维技术也是基于这一内容而展开，其中降维方法又分为线性和非线性降维，非线性降维又分为基于核函数和基于特征值的方法。
  • 线性降维方法：PCA 、ICA LDA、LFA、LPP(LE 的线性表示)
  • 非线性降维方法：
基于核函数的非线性降维方法——KPCA 、KICA、KDA
基于特征值的非线性降维方法(流型学习)——ISOMAP、LLE、LE、LPP、LTSA、MVU
哈尔滨工业大学计算机技术专业的在读硕士生 Heucoder 则整理了 PCA、KPCA、LDA、MDS、ISOMAP、LLE、TSNE、AutoEncoder、FastICA、SVD、LE、LPP 共 12 种经典的降维算法，并提供了相关资料、代码以及展示，下面将主要以 PCA 算法为例介绍降维算法具体操作。

主成分分析(PCA)降维算法 PCA 是一种基于从高维空间映射到低维空间的映射方法，也是最基础的无监督降维算法，其目标是向数据变化最大的方向投影，或者说向重构误差最小化的方向投影。它由 Karl Pearson 在 1901 年提出，属于线性降维方法。与 PCA 相关的原理通常被称为最大方差理论或最小误差理论。这两者目标一致，但过程侧重点则不同。
Python实现12种降维算法的示例代码
文章图片

最大方差理论降维原理
将一组 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。
Python实现12种降维算法的示例代码
文章图片

最小误差理论降维原理
【Python实现12种降维算法的示例代码】而最小误差则是使得平均投影代价最小的线性投影，这一过程中，我们则需要找到的是平方错误评价函数 J0(x0) 等参数。
主成分分析(PCA)代码实现
Python实现12种降维算法的示例代码
文章图片

关于 PCA 算法的代码如下：
from __future__ import print_functionfrom sklearn import datasetsimport matplotlib.pyplot as pltimport matplotlib.cm as cmximport matplotlib.colors as colorsimport numpy as np%matplotlib inlinedef 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]# 正规化数据集 Xdef normalize(X, axis=-1, p=2):lp_norm = np.atleast_1d(np.linalg.norm(X, p, axis))lp_norm[lp_norm == 0] = 1return X / np.expand_dims(lp_norm, axis)# 标准化数据集 Xdef 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 = Xn_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 = Noneself.eigen_vectors = Noneself.k = 2def 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_transformeddef main():# Load the datasetdata = https://www.it610.com/article/datasets.load_iris()X = data.datay = 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 distributionsfor 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 legendplt.legend(class_distr, y, loc=1)# Axis labelsplt.xlabel('Principal Component 1')plt.ylabel('Principal Component 2')plt.show()if __name__ == "__main__":main()

最终，我们将得到降维结果如下。其中，如果得到当特征数 (D) 远大于样本数 (N) 时，可以使用一点小技巧实现 PCA 算法的复杂度转换。
Python实现12种降维算法的示例代码
文章图片

PCA 降维算法展示
当然，这一算法虽然经典且较为常用，其不足之处也非常明显。它可以很好的解除线性相关，但是面对高阶相关性时，效果则较差; 同时，PCA 实现的前提是假设数据各主特征是分布在正交方向上，因此对于在非正交方向上存在几个方差较大的方向，PCA 的效果也会大打折扣。

其它降维算法及代码地址
1.KPCA(kernel PCA)
KPCA 是核技术与 PCA 结合的产物，它与 PCA 主要差别在于计算协方差矩阵时使用了核函数，即是经过核函数映射之后的协方差矩阵。
引入核函数可以很好的解决非线性数据映射问题。kPCA 可以将非线性数据映射到高维空间，在高维空间下使用标准 PCA 将其映射到另一个低维空间。
Python实现12种降维算法的示例代码
文章图片

KPCA 降维算法展示
# coding:utf-8# 实现KPCAfrom sklearn.datasets import load_irisfrom sklearn.decomposition import KernelPCAimport numpy as npimport matplotlib.pyplot as pltfrom scipy.spatial.distance import pdist, squareform'''author: heucoderemail: 812860165@qq.comdate: 2019.6.13'''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 xdef 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)) / NK = 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_nif __name__ == "__main__":data = https://www.it610.com/article/load_iris().dataY = load_iris().targetdata_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()


2.LDA(Linear Discriminant Analysis)
LDA 是一种可作为特征抽取的技术，其目标是向最大化类间差异，最小化类内差异的方向投影，以利于分类等任务即将不同类的样本有效的分开。LDA 可以提高数据分析过程中的计算效率，对于未能正则化的模型，可以降低维度灾难带来的过拟合。
Python实现12种降维算法的示例代码
文章图片

LDA 降维算法展示
#coding:utf-8import numpy as npfrom sklearn.discriminant_analysis import LinearDiscriminantAnalysisfrom sklearn.datasets import load_irisimport matplotlib.pyplot as plt'''author: heucoderemail: 812860165@qq.comdate: 2019.6.13'''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 matrixSw = 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 matrixSB = 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 += SBiS = np.linalg.inv(Sw)*SBeigVals,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_ndimif __name__ == '__main__':iris = load_iris()X = iris.dataY = iris.targetdata_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()

代码地址

3.MDS(multidimensional scaling)
MDS 即多维标度分析，它是一种通过直观空间图表示研究对象的感知和偏好的传统降维方法。该方法会计算任意两个样本点之间的距离，使得投影到低维空间之后能够保持这种相对距离从而实现投影。
由于 sklearn 中 MDS 是采用迭代优化方式，下面实现了迭代和非迭代的两种。
Python实现12种降维算法的示例代码
文章图片

MDS 降维算法展示
代码地址

4.ISOMAP
Isomap 即等度量映射算法，该算法可以很好地解决 MDS 算法在非线性结构数据集上的弊端。
MDS 算法是保持降维后的样本间距离不变，Isomap 算法则引进了邻域图，样本只与其相邻的样本连接，计算出近邻点之间的距离，然后在此基础上进行降维保距。
Python实现12种降维算法的示例代码
文章图片

ISOMAP 降维算法展示
# coding:utf-8import numpy as npfrom sklearn.datasets import make_s_curveimport matplotlib.pyplot as pltfrom sklearn.manifold import Isomapfrom mpl_toolkits.mplot3d import Axes3Ddef floyd(D,n_neighbors=15):Max = np.max(D)*1000n1,n2 = D.shapek = n_neighborsD1 = np.ones((n1,n1))*MaxD_arg = np.argsort(D,axis=1)for i in range(n1):D1[i,D_arg[i,0:k+1]] = D[i,D_arg[i,0:k+1]]for k in range(n1):for i in range(n1):for j in range(n1):if D1[i,k]+D1[k,j]
代码地址

5.LLE(locally linear embedding)
LLE 即局部线性嵌入算法，它是一种非线性降维算法。该算法核心思想为每个点可以由与它相邻的多个点的线性组合而近似重构，然后将高维数据投影到低维空间中，使其保持数据点之间的局部线性重构关系，即有相同的重构系数。在处理所谓的流形降维的时候，效果比 PCA 要好很多。
Python实现12种降维算法的示例代码
文章图片

LLE 降维算法展示
# coding:utf-8import numpy as npfrom sklearn.datasets import make_s_curveimport matplotlib.pyplot as pltfrom sklearn.manifold import LocallyLinearEmbeddingfrom mpl_toolkits.mplot3d import Axes3D'''author: heucoderemail: 812860165@qq.comdate: 2019.6.13'''def make_swiss_roll(n_samples=100, noise=0.0, random_state=None):#Generate a swiss roll dataset.t = 1.5 * np.pi * (1 + 2 * np.random.rand(1, n_samples))x = t * np.cos(t)y = 83 * np.random.rand(1, n_samples)z = t * np.sin(t)X = np.concatenate((x, y, z))X += noise * np.random.randn(3, n_samples)X = X.Tt = np.squeeze(t)return X, tdef 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 distdef get_n_neighbors(data, n_neighbors = 10):''':param data: (n_samples, n_features):param n_neighbors: n nearest neighbors:return: neighbors indexs'''dist = cal_pairwise_dist(data)dist[dist < 0] = 0dist = dist**0.5n = dist.shape[0]N = np.zeros((n, n_neighbors))for i in range(n):index_ = np.argsort(dist[i])[1:n_neighbors+1]N[i] = N[i] + index_return N.astype(np.int32)def lle(data, n_dims = 2, n_neighbors = 10):''':param data:(n_samples, n_features):param n_dims: target n_dims:param n_neighbors: n nearest neighbors:return: (n_samples, n_dims)'''N = get_n_neighbors(data, n_neighbors)n, D = data.shape# prevent Si to smallif n_neighbors > D:tol = 1e-3else:tol = 0# calculate WW = np.zeros((n_neighbors, n))I = np.ones((n_neighbors, 1))for i in range(n):Xi = np.tile(data[i], (n_neighbors, 1)).TNi = data[N[i]].TSi = np.dot((Xi-Ni).T, (Xi-Ni))# magic and why????Si = Si+np.eye(n_neighbors)*tol*np.trace(Si)Si_inv = np.linalg.pinv(Si)wi = (np.dot(Si_inv, I))/(np.dot(np.dot(I.T, Si_inv), I)[0,0])W[:, i] = wi[:,0]print("Xi.shape", Xi.shape)print("Ni.shape", Ni.shape)print("Si.shape", Si.shape)W_y = np.zeros((n, n))for i in range(n):index = N[i]for j in range(n_neighbors):W_y[index[j],i] = W[j,i]I_y = np.eye(n)M = np.dot((I_y - W_y), (I_y - W_y).T)eig_val, eig_vector = np.linalg.eig(M)index_ = np.argsort(np.abs(eig_val))[1:n_dims+1]print("index_", index_)Y = eig_vector[:, index_]return Yif __name__ == '__main__':# X, Y = make_s_curve(n_samples = 500,#noise = 0.1,#random_state = 42)X, Y = make_swiss_roll(n_samples = 500, noise=0.1, random_state=42)data_1 =lle(X, n_neighbors = 30)print(data_1.shape)data_2 = LocallyLinearEmbedding(n_components=2, n_neighbors = 30).fit_transform(X)plt.figure(figsize=(8,4))plt.subplot(121)plt.title("my_LLE")plt.scatter(data_1[:, 0], data_1[:, 1], c = Y)plt.subplot(122)plt.title("sklearn_LLE")plt.scatter(data_2[:, 0], data_2[:, 1], c = Y)plt.savefig("LLE.png")plt.show()

代码地址

6.t-SNE
t-SNE 也是一种非线性降维算法，非常适用于高维数据降维到 2 维或者 3 维进行可视化。它是一种以数据原有的趋势为基础，重建其在低纬度(二维或三维)下数据趋势的无监督机器学习算法。
下面的结果展示参考了源代码，同时也可用 tensorflow 实现(无需手动更新参数)。
Python实现12种降维算法的示例代码
文章图片

t-SNE 降维算法展示
代码地址

7.LE(Laplacian Eigenmaps)
LE 即拉普拉斯特征映射，它与 LLE 算法有些相似，也是以局部的角度去构建数据之间的关系。它的直观思想是希望相互间有关系的点(在图中相连的点)在降维后的空间中尽可能的靠近; 以这种方式，可以得到一个能反映流形的几何结构的解。
Python实现12种降维算法的示例代码
文章图片

LE 降维算法展示
# coding:utf-8import numpy as npimport matplotlib.pyplot as pltfrom sklearn.datasets import load_digitsfrom mpl_toolkits.mplot3d import Axes3D'''author: heucoderemail: 812860165@qq.comdate: 2019.6.13'''def make_swiss_roll(n_samples=100, noise=0.0, random_state=None):#Generate a swiss roll dataset.t = 1.5 * np.pi * (1 + 2 * np.random.rand(1, n_samples))x = t * np.cos(t)y = 83 * np.random.rand(1, n_samples)z = t * np.sin(t)X = np.concatenate((x, y, z))X += noise * np.random.randn(3, n_samples)X = X.Tt = np.squeeze(t)return X, tdef rbf(dist, t = 1.0):'''rbf kernel function'''return np.exp(-(dist/t))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 distdef cal_rbf_dist(data, n_neighbors = 10, t = 1):dist = cal_pairwise_dist(data)dist[dist < 0] = 0n = dist.shape[0]rbf_dist = rbf(dist, t)W = np.zeros((n, n))for i in range(n):index_ = np.argsort(dist[i])[1:1+n_neighbors]W[i, index_] = rbf_dist[i, index_]W[index_, i] = rbf_dist[index_, i]return Wdef le(data,n_dims = 2,n_neighbors = 5, t = 1.0):''':param data: (n_samples, n_features):param n_dims: target dim:param n_neighbors: k nearest neighbors:param t: a param for rbf:return:'''N = data.shape[0]W = cal_rbf_dist(data, n_neighbors, t)D = np.zeros_like(W)for i in range(N):D[i,i] = np.sum(W[i])D_inv = np.linalg.inv(D)L = D - Weig_val, eig_vec = np.linalg.eig(np.dot(D_inv, L))sort_index_ = np.argsort(eig_val)eig_val = eig_val[sort_index_]print("eig_val[:10]: ", eig_val[:10])j = 0while eig_val[j] < 1e-6:j+=1print("j: ", j)sort_index_ = sort_index_[j:j+n_dims]eig_val_picked = eig_val[j:j+n_dims]print(eig_val_picked)eig_vec_picked = eig_vec[:, sort_index_]# print("L: ")# print(np.dot(np.dot(eig_vec_picked.T, L), eig_vec_picked))# print("D: ")# D not equal I ???print(np.dot(np.dot(eig_vec_picked.T, D), eig_vec_picked))X_ndim = eig_vec_pickedreturn X_ndimif __name__ == '__main__':# X, Y = make_swiss_roll(n_samples = 2000)# X_ndim = le(X, n_neighbors = 5, t = 20)## fig = plt.figure(figsize=(12,6))# ax1 = fig.add_subplot(121, projection='3d')# ax1.scatter(X[:, 0], X[:, 1], X[:, 2], c = Y)## ax2 = fig.add_subplot(122)# ax2.scatter(X_ndim[:, 0], X_ndim[:, 1], c = Y)# plt.show()X = load_digits().datay = load_digits().targetdist = cal_pairwise_dist(X)max_dist = np.max(dist)print("max_dist", max_dist)X_ndim = le(X, n_neighbors = 20, t = max_dist*0.1)plt.scatter(X_ndim[:, 0], X_ndim[:, 1], c = y)plt.savefig("LE2.png")plt.show()

代码地址

8.LPP(Locality Preserving Projections)
LPP 即局部保留投影算法，其思路和拉普拉斯特征映射类似，核心思想为通过最好的保持一个数据集的邻居结构信息来构造投影映射，但 LPP 不同于 LE 的直接得到投影结果，它需要求解投影矩阵。
Python实现12种降维算法的示例代码
文章图片

LPP 降维算法展示
# coding:utf-8import numpy as npimport matplotlib.pyplot as pltfrom sklearn.decomposition import PCAfrom mpl_toolkits.mplot3d import Axes3Dfrom sklearn.datasets import load_digits, load_iris'''author: heucoderemail: 812860165@qq.comdate: 2019.6.13'''def make_swiss_roll(n_samples=100, noise=0.0, random_state=None):#Generate a swiss roll dataset.t = 1.5 * np.pi * (1 + 2 * np.random.rand(1, n_samples))x = t * np.cos(t)y = 83 * np.random.rand(1, n_samples)z = t * np.sin(t)X = np.concatenate((x, y, z))X += noise * np.random.randn(3, n_samples)X = X.Tt = np.squeeze(t)return X, tdef rbf(dist, t = 1.0):'''rbf kernel function'''return np.exp(-(dist/t))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 distdef cal_rbf_dist(data, n_neighbors = 10, t = 1):dist = cal_pairwise_dist(data)dist[dist < 0] = 0n = dist.shape[0]rbf_dist = rbf(dist, t)W = np.zeros((n, n))for i in range(n):index_ = np.argsort(dist[i])[1:1 + n_neighbors]W[i, index_] = rbf_dist[i, index_]W[index_, i] = rbf_dist[index_, i]return Wdef lpp(data,n_dims = 2,n_neighbors = 30, t = 1.0):''':param data: (n_samples, n_features):param n_dims: target dim:param n_neighbors: k nearest neighbors:param t: a param for rbf:return:'''N = data.shape[0]W = cal_rbf_dist(data, n_neighbors, t)D = np.zeros_like(W)for i in range(N):D[i,i] = np.sum(W[i])L = D - WXDXT = np.dot(np.dot(data.T, D), data)XLXT = np.dot(np.dot(data.T, L), data)eig_val, eig_vec = np.linalg.eig(np.dot(np.linalg.pinv(XDXT), XLXT))sort_index_ = np.argsort(np.abs(eig_val))eig_val = eig_val[sort_index_]print("eig_val[:10]", eig_val[:10])j = 0while eig_val[j] < 1e-6:j+=1print("j: ", j)sort_index_ = sort_index_[j:j+n_dims]# print(sort_index_)eig_val_picked = eig_val[j:j+n_dims]print(eig_val_picked)eig_vec_picked = eig_vec[:, sort_index_]data_ndim = np.dot(data, eig_vec_picked)return data_ndimif __name__ == '__main__':X = load_digits().datay = load_digits().target# X, y = make_swiss_roll(n_samples = 1000)dist = cal_pairwise_dist(X)max_dist = np.max(dist)print("max_dist", max_dist)data_2d = lpp(X, n_neighbors = 5, t = 0.01*max_dist)data_2 = PCA(n_components=2).fit_transform(X)plt.figure(figsize=(12,6))plt.subplot(121)plt.title("LPP")plt.scatter(data_2d[:, 0], data_2d[:, 1], c = y)plt.subplot(122)plt.title("PCA")plt.scatter(data_2[:, 0], data_2[:, 1], c = y)plt.show()

代码地址
以上就是Python实现12种降维算法的示例代码的详细内容，更多关于Python降维算法的资料请关注脚本之家其它相关文章！

    推荐阅读


    上一篇：python类名和类方法cls修改类变量的值

    下一篇：深入理解C#之枚举