【机器学习】3：Density Peaks聚类算法实现（局部密度聚类算法）

聚类算法主要包括哪些算法？

主要包括：K-means，DBSCAN，Density Peaks聚类（局部密度聚类），层次聚类，谱聚类。

若按照聚类的方式可划分成三类：

1. 类似于K-means，DBSCAN，Density Peaks聚类（局部密度聚类）的依据密度的聚类方式;
2. 类似于层次聚类的依据树状结构的聚类方式;
3. 类似于谱聚类的依据图谱结构的聚类方式。

什么是无监督学习？

• 无监督学习也是相对于有监督学习来说的，因为现实中遇到的大部分数据都是未标记的样本，要想通过有监督的学习就需要事先人为标注好样本标签，这个成本消耗、过程用时都很巨大，所以无监督学习就是使用无标签的样本找寻数据规律的一种方法
• 聚类算法就归属于机器学习领域下的无监督学习方法。

目的是：可以从庞大的样本集合中选出一些具有代表性的样本子集加以标注，再用于有监督学习-可以从无类别信息情况下，寻找表达样本集具有的特征

分类和聚类的区别是什么呢？

• 对于分类来说，在给定一个数据集，我们是事先已知这个数据集是有多少个种类的。比如一个班级要进行性别分类，我们就下意识清楚分为“男生”、“女生”两个类；该班又转入一个同学A，“男ta”就被分入“男生”类；
• 而对于聚类来说，给定一个数据集，我们初始并不知道这个数据集包含多少类，我们需要做的就是将该数据集依照某个“指标”，将相似指标的数据归纳在一起，形成不同的类；
• 分类是一个后续的过程，已知标签数据，再将测试样本分入同标签数据集中；聚类是不知道标签，将“相似指标”的数据强行“撸”在一起，形成各个类。

一、基于局部密度聚类算法——Density Peaks

Density Peaks聚类算法是在2014年 6 月份，由Alex Rodriguez 和 Alessandro Laio 在 Science 上发表了一篇名为《Clustering by fast search and find of density peaks》的文章，这为聚类算法的设计提供了一种新的思路。

1.1、定义局部密度大小——${\rho }_i$

Density Peaks聚类算法要是用文字描述是有一些费解的，我尽量用图解释一下：
下图是一个样本空间点的分布图，一共分布着28个点：

• 事先给定一个邻域半径$d_c$，定义任意i、j两点的距离用$d_{ij}$表示；
• ${\rho }_i$表示i点的密度大小，i点的密度大小是如何确定的呢——以i点为圆心，包含在半径大小为$d_c$的圆内点的个数即为i点的密度大小（与DBSCAN密度确定方法相似）；
• 数学公式如下右所示，即比较两点距离$d_{ij}$与领域半径$d_c$的大小关系，小于表示在圆内，计数1，大于表示圆外，计数0，最后求和。
  

可得的结论：${\rho }_i$越大表示点i的局部密度越大，越有可能成为聚类中心。

1.2、定义聚类中心距离——${\delta }_i$

密度峰聚类算法的巧妙之处：就是在于聚类中心距离 ${\delta }_i$的选定。
根据局部密度的定义，我们可以计算出上图中每个点的密度，依照密度确定聚类中心距离 ${\delta }_i$

1.首先将每个点的密度从大到小排列： ${\rho }_i$ > ${\rho }_j$ > ${\rho }_k$ > …；密度最大的点的聚类中心距离与其他点的聚类中心距离的确定方法是不一样的;
2.先确定密度最大的点的聚类中心距离–i点是密度最大的点，它的聚类中心距离${\delta }_i$等于与i点最远的那个点n到点i的直线距离 $d_{in}$；
3. 再确定其他点的聚类中心距离——其他点的聚类中心距离是等于在密度大于该点的点集合中，与该点距离最小的的那个距离。例如i、j、k的密度都比n点的密度大，且j点离n点最近，则n点的聚类中心距离等于$d_{jn}$
4. 依次确定所有的聚类中心距离$\delta$

聚类中心距离$\delta$的数学式如下：（虽然我觉得这个数学式表达的不是很贴切）

1.3、决策图确定聚类簇核心、簇边缘

---

Density Peaks聚类算法就是依据每个点的局部密度大小${\rho }_i$、聚类中心距离${\delta }_i$的数值，组合(${\rho }_i$，${\delta }_i$)投射到二维坐标系中。先上决策图：

（这里聚类中心距离经过归一化处理，将①号点的聚类中心距离${\delta }_{1,27}$定为1）

从B图中可以清楚的看出来：

分布在右上角区域的是聚类的核心点：周围密度很大，且没有其他核心点；
分布在靠近ρ轴的值是属于正常值：密度虽然大，但是周围有比它更合适作为核心点的点；
分布在靠近δ轴的值是属于噪声点：周围密度小，而且离其他点的距离还远。

现实意义就是：北京联合天津、廊坊等地构成帝都经济群，上海联合无锡、常州、苏州构成长江三角洲经济群，广州深圳形成珠三角经济群。北京、上海、广州深圳相当于聚类核心点，引领发展；天津虽然也是现代化城市，但是由于它离北京很近，更多的资源会流向聚类核心点北京，所以天津只能作为正常点。

1.4、Density Peaks聚类算法的意义

聚类算法中最困惑的地方就是选定K值等于多少才算合适，Density Peaks聚类算法给出了一种比较好的确定K值的方式：定义γi=ρi*δi，得到的乘积比较选取较大的K个点作为聚类中心。

1.5、代码

注意：请修改第10、11行的设定参数

# -*- coding：utf-8 -*-
# -*- author：zzZ_CMing  优快云 address:https://blog.youkuaiyun.com/zzZ_CMing
# -*- 2018/08/22；16:11
# -*- python3.5
import numpy as np
import matplotlib.pyplot as plt
import sklearn.datasets as ds
import matplotlib.colors

min_distance = 4.6           # 邻域半径
points_number = 5            # 随机点个数

# 计算各点间距离、各点点密度(局部密度)大小
def get_point_density(datas,labers,min_distance,points_number):
    # 将numpy.ndarray格式转为list格式，并定义元组大小
    data = datas.tolist()
    laber = labers.tolist()
    distance_all = np.random.rand(points_number,points_number)
    point_density = np.random.rand(points_number)

    # 计算得到各点间距离
    for i in range(points_number):
        for n in range(points_number):
            distance_all[i][n] = np.sqrt(np.square(data[i][0]-data[n][0])+np.square(data[i][1]-data[n][1]))
    print('距离数组:\n',distance_all,'\n')

    # 计算得到各点的点密度
    for i in range(points_number):
        x = 0
        for n in range(points_number):
            if distance_all[i][n] > 0 and distance_all[i][n]< min_distance:
                x = x+1
            point_density[i] = x
    print('点密度数组:', point_density, '\n')
    return distance_all, point_density


# 计算点密度最大的点的聚类中心距离
def get_max_distance(distance_all,point_density,laber):
    point_density = point_density.tolist()
    a = int(max(point_density))
    # print('最大点密度',a,type(a))

    b = laber[point_density.index(a)]
    # print("最大点密度对应的索引：",b,type(b))

    c = max(distance_all[b])
    # print("最大点密度对应的聚类中心距离",c,type(c))

    return c


# 计算得到各点的聚类中心距离
def get_each_distance(distance_all,point_density,data,laber):
    nn = []
    for i in range(len(point_density)):
        aa = []
        for n in range(len(point_density)):
            if point_density[i] < point_density[n]:
                aa.append(n)
        # print("大于自身点密度的索引",aa,type(aa))
        ll = get_min_distance(aa,i,distance_all, point_density,data,laber)
        nn.append(ll)
    return nn


# 获得：到点密度大于自身的最近点的距离
def get_min_distance(aa,i,distance_all, point_density,data,laber):
    min_distance = []
    """
    如果传入的aa为空，说明该点是点密度最大的点，该点的聚类中心距离计算方法与其他不同
    """
    if aa != []:
        for k in aa:
            min_distance.append(distance_all[i][k])
        # print('与上各点距离',min_distance,type(nn))
        # print("最小距离：",min(min_distance),type(min(min_distance)),'\n')
        return min(min_distance)
    else:
        max_distance = get_max_distance(distance_all, point_density, laber)
        return max_distance


def get_picture(data,laber,points_number,point_density,nn):
    # 创建Figure
    fig = plt.figure()
    # 用来正常显示中文标签
    matplotlib.rcParams['font.sans-serif'] = [u'SimHei']
    # 用来正常显示负号
    matplotlib.rcParams['axes.unicode_minus'] = False

    # 原始点的分布
    ax1 = fig.add_subplot(211)
    plt.scatter(data[:,0],data[:,1],c=laber)
    plt.title(u'原始数据分布')
    plt.sca(ax1)
    for i in range(points_number):
        plt.text(data[:,0][i],data[:,1][i],laber[i])

    # 聚类后分布
    ax2 = fig.add_subplot(212)
    plt.scatter(point_density.tolist(),nn,c=laber)
    plt.title(u'聚类后数据分布')
    plt.sca(ax2)
    for i in range(points_number):
        plt.text(point_density[i],nn[i],laber[i])

    plt.show()


def main():
    # 随机生成点坐标
    data, laber = ds.make_blobs(points_number, centers=points_number, random_state=0)
    print('各点坐标：\n', data)
    print('各点索引：', laber, '\n')

    # 计算各点间距离、各点点密度(局部密度)大小
    distance_all, point_density = get_point_density(data, laber, min_distance, points_number)
    # 得到各点的聚类中心距离
    nn = get_each_distance(distance_all, point_density, data, laber)
    print('最后的各点点密度：', point_density.tolist())
    print('最后的各点中心距离：', nn)

    # 画图
    get_picture(data, laber, points_number, point_density, nn)
    """
    距离归一化：就把上面的nn改为：nn/max(nn)
    """


if __name__ == '__main__':
    main()

结果展示：

1、由于随机生成的数据大多是分散的，所以得到的实验效果不是很理想；

2、有这方面需求的伙伴们可以自己做拓展——将随机数据生成改为引入自己在txt或excel中的数据做测试；

3、高维空间也是可以的，但时间复杂度和内存消耗我并没有估计，改进也留给你们了。

---

二、2024年内容：

import numpy as np
from sklearn.metrics import pairwise_distances
import matplotlib.pyplot as plt

# 生成距离矩阵
def compute_distance_matrix(X):
    return pairwise_distances(X, X)

# 计算每个点的局部密度
def compute_local_density(dist_matrix, dc):
    rho = np.sum(np.exp(-(dist_matrix / dc) ** 2), axis=1) - 1  # 减1是为了排除自身的密度
    return rho

# 计算每个点的delta值，找到其距离最近的密度更高的点
def compute_delta(dist_matrix, rho):
    delta = np.full_like(rho, np.inf)
    nearest_neighbor = np.full_like(rho, -1, dtype=int)

    sorted_indices = np.argsort(-rho)  # 按密度降序排列
    for i in range(1, len(rho)):
        for j in range(i):
            if dist_matrix[sorted_indices[i], sorted_indices[j]] < delta[sorted_indices[i]]:
                delta[sorted_indices[i]] = dist_matrix[sorted_indices[i], sorted_indices[j]]
                nearest_neighbor[sorted_indices[i]] = sorted_indices[j]
    
    delta[sorted_indices[0]] = np.max(delta)  # 密度最高的点的delta设为最大值
    return delta, nearest_neighbor

# 聚类函数
def density_peaks_clustering(X, dc, num_clusters):
    dist_matrix = compute_distance_matrix(X)
    rho = compute_local_density(dist_matrix, dc)
    delta, nearest_neighbor = compute_delta(dist_matrix, rho)

    # 选择聚类中心
    cluster_centers = np.argsort(-rho * delta)[:num_clusters]
    labels = -np.ones(len(X), dtype=int)

    # 将每个点分配到距离最近的密度更高的点所在的簇
    for center_idx, center in enumerate(cluster_centers):
        labels[center] = center_idx
    for i in range(len(X)):
        if labels[i] == -1:
            labels[i] = labels[nearest_neighbor[i]]
    
    return labels, cluster_centers, rho, delta

# 可视化结果
def plot_density_peaks(X, labels, rho, delta, cluster_centers):
    plt.figure(figsize=(10, 6))
    
    # 绘制数据点
    scatter = plt.scatter(X[:, 0], X[:, 1], c=labels, cmap='viridis', s=30)
    plt.colorbar(scatter)
    
    # 绘制聚类中心
    plt.scatter(X[cluster_centers, 0], X[cluster_centers, 1], s=200, c='red', marker='*', label='Cluster Centers')
    
    plt.xlabel('Feature 1')
    plt.ylabel('Feature 2')
    plt.title('Density Peaks Clustering')
    plt.legend()
    plt.show()

# 生成示例数据并执行算法
if __name__ == "__main__":
    from sklearn.datasets import make_blobs

    # 生成随机数据
    X, _ = make_blobs(n_samples=300, centers=4, cluster_std=0.5, random_state=0)
    
    # 设置密度截断距离和聚类数量
    dc = 0.5  # 选择合适的dc参数
    num_clusters = 4  # 设置预期聚类数量

    # 执行聚类算法
    labels, cluster_centers, rho, delta = density_peaks_clustering(X, dc, num_clusters)
    
    # 绘制聚类结果
    plot_density_peaks(X, labels, rho, delta, cluster_centers)

代码解释

1. 计算距离矩阵：compute_distance_matrix函数使用欧氏距离计算点之间的距离矩阵。
2. 计算局部密度：compute_local_density函数基于密度截断距离dc，计算每个点的局部密度rho，即每个点与其他点的距离之和，距离越小，密度越大。
3. 计算距离和密度峰值：compute_delta函数计算每个点到密度更高点的最小距离delta。其中密度最高的点的delta设为最大值。
4. 选择聚类中心和分配标签：density_peaks_clustering函数根据密度和距离选择前num_clusters个点作为聚类中心，并将每个点分配到其最邻近的密度更高的点所在的簇。
5. 可视化：plot_density_peaks函数显示聚类结果，用颜色区分不同的簇，红色星号标记聚类中心。

参数说明：

• dc：密度截断距离，用于确定局部密度。需要根据数据规模和分布来调节，通常通过计算数据集的距离分布或经验方法选取。
• num_clusters：期望的簇数量，影响算法如何选择聚类中心。

Density Peaks Clustering 在具有明显密度分布的聚类任务中效果较好，适合用于检测数据集中具有较高密度的中心点。