python_spectral_clustering.zip

""" 谱聚类 在所有样本中 找出不包含 xi 的 k 个最邻近的样本点 用 wij 构造图 """ import numpy as np import copy import matplotlib.pyplot as plt from K_Means import kmeans from sklearn.cluster import KMeans f1 = open("data_moon.csv", "rb") case_train = np.loadtxt(f1, delimiter=',', skiprows=0) f1.close() case_train = np.array(case_train) # 谱聚类 class spectral: def __init__(self, sigma, knn_num): self.sigma = sigma self.knn_num = knn_num # wij矩阵元素 def wij(self, xi, xj, sigma): out = np.exp(- (pow((xi - xj)[0], 2) + pow((xi - xj)[1], 2)) / (2 * sigma * sigma)) return out # Calculate 2_Dim Distance # 每个点以 1x2 or 2x1 维度输入 # o-距离 def calculate_2d_dist(self, point1, point2): dist = pow((point1 - point2)[0], 2) + pow((point1 - point2)[1], 2) dist = np.sqrt(dist) return dist # knn 找到k近邻的组合 # k : k个近邻 def knn(self, all_sample): sample_len = len(all_sample) knn_output = [] index = [] for n in range(sample_len): dist_list = [] knn_group = [] knn_group_index = [] for i in range(sample_len): dist = self.calculate_2d_dist(all_sample[n], all_sample[i]) # 第n个样本到全部样本的距离 dist_list.append(dist) dist_list = np.array(dist_list) min_index = dist_list.argsort() # min 自小至大 # [-3:][::-1] # max for j in range(self.knn_num + 1): # 把k个元素提取到 knn——group if min_index[j] != n: # 如果是他自己就不记录 knn_group.append(all_sample[min_index[j]]) # 只有k个了 knn_group_index.append(min_index[j]) knn_output.append(knn_group) # 此时输出有 sample_len 个元素，美元素包涵k个最近点 index.append(knn_group_index) # 每元素包涵k个近邻点的索引 return knn_output, index # knn——group i组内的i＊k个样本 / knn——group对应index 第i组第k个样本的索引 # 构建图 矩阵 # 输入 k_index : def get_wd(self, all_sample, k_group, k_index): data_len = len(all_sample) k_len = len(k_group[0]) # 实际上是knn中的k个近邻 w = np.zeros((data_len, data_len)) # 初始化 权重 nxn d = np.zeros((data_len, data_len)) # 得到 d矩阵 所有与该顶点相连接的边的权重之和 for i in range(data_len): # 计算当前结点到每一个邻近的wij temp = 0 for k in range(k_len): # 计算到k个近邻的相似度 w[i][k_index[i][k]] = self.wij(all_sample[i], k_group[i][k], self.sigma) temp = temp + w[i][k_index[i][k]] # 到 i 的权重和 to d d[i][i] = temp # 保证W矩阵的对称性 w = np.array(w) w = (w.transpose() + w) / 2 return w, d """ 无向图切图 构造拉普拉斯矩阵 通过k个邻近值，计算出亲和度矩阵w和d矩阵 再构建拉普拉斯矩阵 输入 k_num 切成k类 w_matrix d_matrix return 排序好的特征向量 矩阵 """ def cut_graph(self, k_num, w_matrix, d_matrix): w_shape = w_matrix.shape min_vectors = [] # 存贮特征值最小的 对应k列 # RatioCut """ l, w = np.linalg.eig(w_matrix) min_index = l.argsort() """ # NCut # Lsym = I − D−1/2 * W * D−1/2 d_12 = np.diag(np.power(d_matrix + 1e-5, -0.5)) # 防止除以0溢出 转为对角是因为别的地方是0 D = np.diag(d_12) lsym = np.eye(w_shape[0]) - np.dot(np.dot(D, w_matrix), D) #print(lsym) l, w = np.linalg.eig(lsym) l = np.array(l) min_index = abs(l).argsort() # 所有的特征值排序 从小到大 应该按绝对值排序了 # print(l) # 发现特征值竟然有负的！？？？ for j in range(k_num): min_vectors.append(w[min_index[j]]) # 前k个 vectors = np.array(min_vectors).transpose() return vectors # 按行做标准化 def normalize(self, vectors): vec_len = len(vectors) dim_len = len(vectors[0]) new_vectors = [] for i in range(vec_len): # 遍历维度 temp = 0 # 行元素平方和 temp_vec = [] for j in range(dim_len): temp = temp + pow(vectors[i][j], 2) temp = np.sqrt(temp) if temp != 0: # 如果平方和等于0 那就不算了 for m in range(dim_len): elem = vectors[i][m] / temp temp_vec.append(elem) new_vectors.append(temp_vec) else: new_vectors.append(vectors[i]) new_vectors = np.array(new_vectors) #new_vectors = abs(new_vectors) # 我取了绝对值？？？ return new_vectors # 类数／原始索引 def test_acc(self, data_index): last_num = 0 class_num = len(data_index) right_num = 0 all_num = 0 for i in range(class_num): this_num = len(data_index[i]) for j in range(this_num): if data_index[i][j] == (last_num * i + j): right_num = right_num + 1 all_num = all_num + 1 last_num = this_num acc = right_num / all_num return acc def change_k(): data = case_train acc_list = [] # 不同的k for k_val in range(4, 10): # sigma ／ k (4-9) sp = spectral(1, k_val) # 输入 all_sample / return : k 近邻组 k近邻对应索引 k_group, k_index = sp.knn(data) # w d 矩阵 w_matrix, d_matrix = sp.get_wd(data, k_group, k_index) # in ： k_num 类 , w_matrix, d_matrix vectors = sp.cut_graph(2, w_matrix, d_matrix) # 切图 返回 nxk 特征向量集合 vectors = sp.normalize(vectors) # input : data_set, k_num , step_num / return : class_list, final_center, point_cloud k0 = kmeans(vectors, k_num=2, step_num=10) # 应该改为返回 分类后的原始数据 class_list, trace, final_center, index_list = k0.kmeans_loop() acc = sp.test_acc(index_list) acc_list.append(acc) return acc_list def change_sigma(): data = case_train acc_list = [] # 不同的sigma for n_val in np.arange(0.8, 1.3, 0.1): # sigma(0.8, 1.2) ／ k (4-9) sp = spectral(n_val, 5) # 输入 all_sample / return : k 近邻组 k近邻对应索引 k_group, k_index = sp.knn(data) # w d 矩阵 w_matrix, d_matrix = sp.get_wd(data, k_group, k_index) # in ： k_num 类 , w_matrix, d_matrix vectors = sp.cut_graph(2, w_matrix, d_matrix) # 切图 返回 nxk 特征向量集合 vectors = sp.normalize(vectors) # input : data_set, k_num , step_num / return : class_list, final_center, point_cloud k0 = kmeans(vectors, k_num=2, step_num=10) # 应该改为返回 分类后的原始数据 class_list, trace, final_center, index_list = k0.kmeans_loop() acc = sp.test_acc(index_list) acc_list.append(acc) return acc_list def Run(): acc_list_k = change_k() acc_list_sig = change_sigma() print(acc_list_k) print(acc_list_sig) fig, axs = plt.subplots(1, 2, figsize=(8, 4)) # , sharey=True) # 原始 ／ 聚类后 axs[0].plot(list(range(4,10)), acc_list_k, "-") axs[0].set_title('k') axs[1].plot(np.arange(0.8, 1.3, 0.1), acc_list_sig, "--") axs[1].set_title('Sigma') plt.show() if __name__ == '__main__': Run()