【代码复现】SMDP

1、算法

关于密度峰值聚类的半监督多示例学习（semi-supervised multi-instance learning,SMDP）的具体介绍在这一篇中已经有了详细的介绍。代码的关键部分也是与论文中一一对应的。

2、代码

2.1 util.py

这一部分主要是一些基本操作，其中包括了：加载数据集、获取训练集与测试集的索引。

import numpy as np
from scipy.io import loadmat


def load_data(path=''):#加载数据集
    data = loadmat(path)['data']#加载数据集，并命名为data
    bags, labels = [], []#创建包与标签数组
    for i in range(len(data)):#将包与标签装入各自数组中
        bags.append(data[i, 0][:, :-1])
        labels.append(data[i, 1][0, 0])
    labels = np.array(labels)#将标签
    return bags, labels


def get_index(num_bags=92, para_k=10, seed=None):#获取索引
    if seed is not None:#若没有指定的随机数
        np.random.seed(seed)#生成指定随机数
    temp_rand_idx=np.random.permutation(num_bags)#随机生成索引

    temp_fold = int(np.ceil(num_bags/para_k))#找到折叠的位置，将数据集分为训练集与测试集两部分
    ret_tr_idx={}#要返回的训练集索引
    ret_te_idx={}#要返回的测试集索引
    for i in range(para_k):
        temp_tr_idx=temp_rand_idx[0:i*temp_fold].tolist()#将0~折叠位置的索引作为训练集索引
        temp_tr_idx.extend(temp_rand_idx[(i+1)*temp_fold:])
        ret_tr_idx[i]=temp_tr_idx
        ret_te_idx[i]=temp_rand_idx[i*temp_fold: (i + 1) * temp_fold].tolist()#剩余索引作为测试集
    return ret_tr_idx,ret_te_idx#返回训练集索引与测试集索引


if __name__=='__main__':
    load_data('')

2.2 DP_TWT.py

这一部分主要是进行Density peaks(DP) clustering algorithm方法的设置，具体原理见本文。这一段代码的关键就在于找到实例密度$\rho$以及实例到其master的距离$\delta$。最后将以上两个变量相乘得到$\lambda$用于选出聚类中心。

import numpy as np
import math
import time


class DP:
    def __init__(self):
        self.dis_matrix = np.array([])#距离矩阵
        self.num_bags = 0#包数量
        self.d_c , self.n_c = 0,0#
        self.local_density_list = []#实例密度
        self.dis_to_master_list = []#到上司的距离
        self.lambda_list = []#密度与距离之积
        self.center_list = []#聚类中心列表

    def train(self, matrix, r=0.2, u=0.4, kernel='gaussian'):
        '''

        :param matrix: 距离矩阵
        :param r: cut-off distance阶段距离
        :param u: cut-off radio划分比率
        :param kernel: 本次DP聚类所用核函数，为高斯函数
        :return:
        '''
        self.dis_matrix = np.array(matrix)
        self.num_bags = len(self.dis_matrix)
        self.d_c = r * np.max(matrix)
        self.n_c = int(self.num_bags * u)
        self.local_density_list = self.get_local_density_list(kernel)
        self.dis_to_master_list = self.get_dis_to_master_list()

        self.lambda_list = self.get_lambda_list()

        self.center_list = self.get_center_list()

    def get_local_density_list(self, kernel):#用于计算密度
        if kernel == 'gaussian':#如果使用高斯核函数
            temp_local_density_list = np.sum(1 / (np.exp(np.power(self.dis_matrix / self.d_c, 2))), axis=1)
            return temp_local_density_list
        elif kernel == 'cutoff':#如果使用截断核函数
            temp_local_density_list = []
            for i in range(0, self.num_bags):
                temp_local_density_list.append(self.cutoff_kernel(i))
            return temp_local_density_list

    def get_dis_to_master_list(self):
        temp_dis_to_master_list = []#初始化数组
        for i in range(0, self.num_bags):
            density_list = []  # 记录所有局部密度比第i个包更大的包的序号
            i_density = self.local_density_list[i]#计算第i个包的密度
            for j in range(0, self.num_bags):#找其他的包
                if self.local_density_list[j] > i_density:#如果找到一个包密度比第i个包还大的
                    density_list.append(j)#将其记录到数组中
            dis_list = []  # 记录所有局部密度大于第i个包的包与i包的距离
            for k in density_list:#计算该包到密度较大的包之间的距离
                dis_list.append(self.dis_matrix[i, k])#将距离记录到数组中
            if dis_list:
                dis_list.sort()#对距离进行排序
                temp_dis_to_master_list.append(dis_list[0])  # 返回i包与其master的距离
            else:
                temp_dis_to_master_list.append(float('inf'))  # 表示此包为局部密度最大的包
        return np.array(temp_dis_to_master_list)#返回到其上司距离的数组

    def get_lambda_list(self):
        lambda_list = np.multiply(self.local_density_list, self.dis_to_master_list)
        #将距离与密度相乘，得到lambda
        return lambda_list

    def get_center_list(self):
        temp_center_list = np.argsort(self.lambda_list)[-self.n_c:]
        #通过进行排序，并输出对应索引数组
        # [10, 50, 30]
        # np.argsort([10, 30, 50]) -> 0, 2, 1
        return np.array(temp_center_list)

2.3 SMDP.py

这一部分是关键，通过豪斯多夫距离来测量实例之间的距离。通过交叉验证的方式来检测算法的准确度。

import numpy as np
import os
from utils import load_data, get_index
from sklearn.metrics import euclidean_distances as eucl
from tqdm import tqdm
from DP_TWT import DP
from sklearn.svm import SVC
from sklearn.neighbors import KNeighborsClassifier


class SMDP:
    def __init__(self,path,measure,kernel,r,u,k,n):
        '''

        :param path: 数据路径
        :param measure: 距离测量方式
        :param kernel: 聚类核函数
        :param r:cut-off distance阶段距离
        :param u:cut-off radio划分比率
        :param k:划分训练集与测试集
        :param n:交叉验证次数
        '''
        self.path=path
        self.bags, self.labels = load_data(path)
        self.num_bags=len(self.labels)
        self.measure=measure
        self.kernel=kernel
        self.r,self.u=r,u
        self.k,self.n=k,n
        self.dis_matrix=self.get_dis_matrix()#计算包间距离，并存入距离矩阵中
        dp=DP()#初始化DP聚类
        dp.train(self.dis_matrix, r=self.r, u=self.u, kernel=kernel)#选定核函数，对包进行聚类
        center_idx=dp.center_list#将聚类中心存入索引中
        self.vectors=self.dis_matrix[:,center_idx]#将索引中找到的距离中心在距离矩阵中找出来，并存为向量


    def one_cv(self):#交叉验证操作
            tr_idx_list,te_idx_list=get_index(self.num_bags,k)
            acc_list=[]#准确度存储矩阵
            for i in range(self.k):
                acc=self.run(tr_idx_list[i],te_idx_list[i])#计算本次验证准确度
                acc_list.append(acc)#将结果加入矩阵中
            return np.mean(acc_list)#计算准确度均值

    def n_cv(self):#进行多次交叉验证操作
            acc_list=[]#准确度存储矩阵
            for i in range(self.n):#进行n次操作
                acc = self.one_cv()#进行本次交叉验证
                acc_list.append(acc)#将结果加入矩阵中
            return float(np.mean(acc_list)),float(np.std(acc_list,ddof=1))

    def run(self, tr_idx, te_idx):
            tr_vec = self.vectors[tr_idx]#创建训练数据向量矩阵
            te_vec = self.vectors[te_idx]#创建测试数据向量矩阵
            classfier=KNeighborsClassifier(n_neighbors=3)#使用KNN算法进行预测
            classfier.fit(tr_vec,self.label[tr_idx])
            pred_labels = classfier.predict(te_vec)#得到预测标签
            acc=np.sum(pred_labels==self.label[te_idx]/len(te_idx))#计算预测准确率
            return acc*100#输出结果

    def ave_hausdorff(self, i, j):#平均豪斯多夫距离，计算实例i与实例j之间的距离
            if i==j:
                return 0
            temp_1 = np.sum(np.min(eucl(self.bags[i], self.bags[j]), axis=1))
            temp_2 = np.sum(np.min(eucl(self.bags[j], self.bags[i]), axis=1))
            result = (temp_1 + temp_2) / (self.bags[i].shape[0] + self.bags[j].shape[0])
            return result

    def get_dis_matrix(self):
        dis_path = '../distance/' + self.path.split('/')[-1].split('.')[0] + '_' + self.measure + '.npy'
        if os.path.exists(dis_path):  # 若距离数据路径存在
            print("load computed distance matrix")
            dis_matrix = np.load(dis_path)  # 读取数据
            return dis_matrix
        # 若路径不存在，则必须重新计算
        print('computing distance matrix')
        dis_matrix = np.zeros((self.num_bags, self.num_bags))  # 创建包大小的零矩阵作为距离矩阵
        if self.measure == 'ave_h':  # 若距离度量方式为平均豪斯夫距离
            for i in tqdm(range(self.num_bags), desc='computing...'):  # 循环的同时也显示一个进度条
                for j in range(i, self.num_bags):
                    dis_matrix[i, j] = self.ave_hausdorff(i, j)  # 计算实例间的距离
                    dis_matrix[j, i] = self.ave_hausdorff(i, j)
        np.save(dis_path, dis_matrix)  # 存储到具体路径中
        return dis_matrix  # 返回距离矩阵


if __name__ == '__main__':
    path = ''#数据文件路径
    measure = 'ave_h'#距离度量方式为平均豪斯夫
    kernel = 'gaussian'#密度核函数为高斯核
    r, u = 0.2, 0.6
    k, n = 10, 10  # n次k vc
    smdp = SMDP(path, measure, kernel, r, u, k, n)
    acc, std = smdp.n_cv()
    print(path.split('/')[-1])
    print('$%.1f_{%.1f}$' % (acc, std))