python 数组降维,python自编码器降维
数据降维算法是机器学习算法的一大家族，不同于分类、回归、聚类等算法。它的目标是将向量投影到低维空间中，以达到某种目的，如可视化或分类。本文将使用Python实现12种降维算法，有需要的可以参考。
00-1010为什么要降低数据的维度？主成分分析(PCA)，其他降维算法和代码地址1。KPCA(内核主成分分析)2。线性判别分析。MDS(多维标度)4。Isomap 5。LLE(局部线性嵌入)6。T-SNE7。乐(拉普拉斯八映射)8。LPP(保局投影)网上关于各种降维算法的信息参差不齐，大部分都不提供源代码。这里有一个GitHub项目，用Python整理并实现了11个经典的数据抽取(数据降维)算法，包括PCA、LDA、MDS、LLE、TSNE等。附相关材料和示范效果；非常适合机器学习的初学者和刚入坑开始数据挖掘的朋友。
目录
所谓降维，就是用一组编号为D的向量Zi代表编号为D的向量Xi所包含的有用信息，其中D
通常，我们会发现大多数数据集有数百甚至数千维，而经典的MNIST有64维。
MNIST手写数字数据集
但在实际应用中，我们使用的有用信息并不需要那么高的维度，每增加一个维度所需的样本数呈指数级增长，这可能直接导致极大的“维度灾难”。可以实现数据降维：
使数据集更容易使用，确保变量相互独立，降低算法的计算成本，去除噪声。一旦我们能够正确处理这些信息，正确有效地降维，将极大地有助于减少计算量，从而提高机器的效率。数据降维常用于文本处理、人脸识别、图片识别、自然语言处理等领域。
为什么要进行数据降维
往往高维空间中的数据会出现稀疏分布，所以在降维的过程中，我们通常会做一些数据删除，包括冗余数据、无效信息、重复表达式等等。
比如有一个现有的1024*1024的图，除了以50*50为中心的区域，其他位置都是零值，这些零值信息可以归为无用信息；对于对称图形，对称部分的信息可以归为重复信息。
因此，经典的降维技术大多基于这一内容，其中降维方法又分为线性和非线性降维方法，非线性降维方法又分为基于核函数的方法和基于特征值的方法。
线性降维方法：PCA、ICA LDA、LFA、LPP(LE的线性表示)非线性降维方法：基于核的非线性降维方法——KPCA、KICA、KDA
基于特征值的非线性降维方法(流模式学习)——ISOMAP，LLE，乐，MVU
哈工大计算机技术专业硕士生Heucoder整理了PCA、KPCA、LDA、MDS、ISOMAP、LLE、TSNE、AutoEncoder、FastICA、SVD、LE、LPP等12种经典降维算法，并提供了相关素材、代码和展示。下面将以PCA算法为例，主要介绍降维算法的具体操作。
数据降维原理
PCA是一种基于高维空间到低维空间映射的映射方法，也是最基本的无监督降维算法。它的目标是投影到数据变化最大的方向，或者最小化重建误差。它是由卡尔皮尔逊在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。
最小误差理论降维原理
而最小误差则是使得平均投影代价最小的线性投影，这一过程中，我们则需要找到的是平方错误评价函数 J0(x0) 等参数。
主成分分析(PCA)代码实现
关于 PCA 算法的代码如下：
from __future__ import print_functionfrom 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()
最终，我们将得到降维结果如下。其中，如果得到当特征数 (D) 远大于样本数 (N) 时，可以使用一点小技巧实现 PCA 算法的复杂度转换。
PCA 降维算法展示
当然，这一算法虽然经典且较为常用，其不足之处也非常明显。它可以很好的解除线性相关，但是面对高阶相关性时，效果则较差;同时，PCA 实现的前提是假设数据各主特征是分布在正交方向上，因此对于在非正交方向上存在几个方差较大的方向，PCA 的效果也会大打折扣。
其它降维算法及代码地址
1.KPCA(kernel PCA)
KPCA 是核技术与 PCA 结合的产物，它与 PCA 主要差别在于计算协方差矩阵时使用了核函数，即是经过核函数映射之后的协方差矩阵。
引入核函数可以很好的解决非线性数据映射问题。kPCA 可以将非线性数据映射到高维空间，在高维空间下使用标准 PCA 将其映射到另一个低维空间。
KPCA 降维算法展示
# 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
author: heucoder
email: 812860165@qq.com
date: 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 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()
2.LDA(Linear Discriminant Analysis)
LDA 是一种可作为特征抽取的技术，其目标是向最大化类间差异，最小化类内差异的方向投影，以利于分类等任务即将不同类的样本有效的分开。LDA 可以提高数据分析过程中的计算效率，对于未能正则化的模型，可以降低维度灾难带来的过拟合。
LDA 降维算法展示
#coding:utf-8import numpy as np
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis
from sklearn.datasets import load_iris
import matplotlib.pyplot as plt
author: heucoder
email: 812860165@qq.com
date: 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 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()
代码地址
3.MDS(multidimensional scaling)
MDS 即多维标度分析，它是一种通过直观空间图表示研究对象的感知和偏好的传统降维方法。该方法会计算任意两个样本点之间的距离，使得投影到低维空间之后能够保持这种相对距离从而实现投影。
由于 sklearn 中 MDS 是采用迭代优化方式，下面实现了迭代和非迭代的两种。
MDS 降维算法展示
代码地址
4.ISOMAP
Isomap 即等度量映射算法，该算法可以很好地解决 MDS 算法在非线性结构数据集上的弊端。
MDS 算法是保持降维后的样本间距离不变，Isomap 算法则引进了邻域图，样本只与其相邻的样本连接，计算出近邻点之间的距离，然后在此基础上进行降维保距。
ISOMAP 降维算法展示
# coding:utf-8import numpy as np
from sklearn.datasets import make_s_curve
import matplotlib.pyplot as plt
from sklearn.manifold import Isomap
from mpl_toolkits.mplot3d import Axes3D
def floyd(D,n_neighbors=15):
Max = np.max(D)*1000
n1,n2 = D.shape
k = n_neighbors
D1 = np.ones((n1,n1))*Max
D_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]<D1[i,j]:
D1[i,j] = D1[i,k]+D1[k,j]
return D1
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(dist, n_dims):
# dist (n_samples, n_samples)
dist = dist**2
n = dist.shape[0]
T1 = np.ones((n,n))*np.sum(dist)/n**2
T2 = np.sum(dist, axis = 1)/n
T3 = np.sum(dist, axis = 0)/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_]
return picked_eig_vector*picked_eig_val**(0.5)
def my_Isomap(data,n=2,n_neighbors=30):
D = cal_pairwise_dist(data)
D[D < 0] = 0
D = D**0.5
D_floyd=floyd(D, n_neighbors)
data_n = my_mds(D_floyd, n_dims=n)
return data_n
def scatter_3d(X, y):
fig = plt.figure(figsize=(6, 5))
ax = fig.add_subplot(111, projection=3d)
ax.scatter(X[:, 0], X[:, 1], X[:, 2], c=y, cmap=plt.cm.hot)
ax.view_init(10, -70)
ax.set_xlabel("$x_1$", fontsize=18)
ax.set_ylabel("$x_2$", fontsize=18)
ax.set_zlabel("$x_3$", fontsize=18)
plt.show()
if __name__ == __main__:
X, Y = make_s_curve(n_samples = 500,
noise = 0.1,
random_state = 42)
data_1 = my_Isomap(X, 2, 10)
data_2 = Isomap(n_neighbors = 10, n_components = 2).fit_transform(X)
plt.figure(figsize=(8,4))
plt.subplot(121)
plt.title("my_Isomap")
plt.scatter(data_1[:, 0], data_1[:, 1], c = Y)
plt.subplot(122)
plt.title("sklearn_Isomap")
plt.scatter(data_2[:, 0], data_2[:, 1], c = Y)
plt.savefig("Isomap1.png")
plt.show()
代码地址
5.LLE(locally linear embedding)
LLE 即局部线性嵌入算法，它是一种非线性降维算法。该算法核心思想为每个点可以由与它相邻的多个点的线性组合而近似重构，然后将高维数据投影到低维空间中，使其保持数据点之间的局部线性重构关系，即有相同的重构系数。在处理所谓的流形降维的时候，效果比 PCA 要好很多。
LLE 降维算法展示
# coding:utf-8import numpy as np
from sklearn.datasets import make_s_curve
import matplotlib.pyplot as plt
from sklearn.manifold import LocallyLinearEmbedding
from mpl_toolkits.mplot3d import Axes3D
author: heucoder
email: 812860165@qq.com
date: 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.T
t = np.squeeze(t)
return X, 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 dist
def 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] = 0
dist = dist**0.5
n = 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 small
if n_neighbors > D:
tol = 1e-3
else:
tol = 0
# calculate W
W = np.zeros((n_neighbors, n))
I = np.ones((n_neighbors, 1))
for i in range(n):
Xi = np.tile(data[i], (n_neighbors, 1)).T
Ni = data[N[i]].T
Si = 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 Y
if __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 实现(无需手动更新参数)。
t-SNE 降维算法展示
代码地址
7.LE(Laplacian Eigenmaps)
LE 即拉普拉斯特征映射，它与 LLE 算法有些相似，也是以局部的角度去构建数据之间的关系。它的直观思想是希望相互间有关系的点(在图中相连的点)在降维后的空间中尽可能的靠近;以这种方式，可以得到一个能反映流形的几何结构的解。
LE 降维算法展示
# coding:utf-8import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import load_digits
from mpl_toolkits.mplot3d import Axes3D
author: heucoder
email: 812860165@qq.com
date: 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.T
t = np.squeeze(t)
return X, t
def 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 dist
def cal_rbf_dist(data, n_neighbors = 10, t = 1):
dist = cal_pairwise_dist(data)
dist[dist < 0] = 0
n = 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 W
def 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 - W
eig_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 = 0
while eig_val[j] < 1e-6:
j+=1
print("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_picked
return X_ndim
if __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().data
y = load_digits().target
dist = 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 的直接得到投影结果，它需要求解投影矩阵。
LPP 降维算法展示
# coding:utf-8import numpy as np
import matplotlib.pyplot as plt
from sklearn.decomposition import PCA
from mpl_toolkits.mplot3d import Axes3D
from sklearn.datasets import load_digits, load_iris
author: heucoder
email: 812860165@qq.com
date: 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.T
t = np.squeeze(t)
return X, t
def 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 dist
def cal_rbf_dist(data, n_neighbors = 10, t = 1):
dist = cal_pairwise_dist(data)
dist[dist < 0] = 0
n = 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 W
def 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 - W
XDXT = 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 = 0
while eig_val[j] < 1e-6:
j+=1
print("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_ndim
if __name__ == __main__:
X = load_digits().data
y = 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降维算法的资料请关注盛行IT软件开发工作室其它相关文章！
郑重声明：本文由网友发布，不代表盛行IT的观点，版权归原作者所有，仅为传播更多信息之目的，如有侵权请联系，我们将第一时间修改或删除，多谢。