面向地理大数据的新型密度山峰聚类算法（一）——算法描述

最新推荐文章于 2024-05-30 20:45:00 发布

万万VV

最新推荐文章于 2024-05-30 20:45:00 发布

阅读量1.2k

点赞数 1

文章标签：聚类机器学习 Python 大数据

本文链接：https://blog.csdn.net/weixin_39626452/article/details/101199909

版权

摘要

通过实现和改进在《Science》杂志上发表的一种新型快速密度山峰聚类算法[1],与现有的常用聚类方法进行分析比较，发现其优点和不足的地方，并实际应用于地理定位大数据的分析中，以快速发现位置数据中任意形状的聚类簇模式和噪声。首先，基于类簇中心被具有较低局部密度的邻居点包围，且与具有更高密度的任何点有相对较大的距离这样的假设提出基本算法。其次，通过Python实现该算法并采用基本聚类测试数据集进行效果测试和算法评估分析。最后，用在基准测试数据上验证了所提算法的聚类效果,以及在位置大数据上的实验结果统计和原先普遍所采用的的方法，如基于密度的方法DBSCAN[2]、基于划分的方法K-means[3]、基于网格的方法（如GDCA[4]和GG[5]）基于路径的谱聚类算法[6]等进行比较。

原论文阅读链接：Clustering by fast search and find of density peaks

一、引言

聚类分析旨在依赖所要聚类元素之间的相似度以期将他们分成不同的类。该方法的应用领域范围包括航空航天、生物化学、文献计量学、模式识别等。基于聚类中心具有比相邻点密度更高以及它们之间具有相当远的距离的特点，作者依此提出了一种新的聚类方法。该种方法最大的突破就是跳出了K-means[3]方法的桎梏，采用了一种较新的思想，且群的数目直观可见、噪声点自动标出并排除在分析之外、群的区分不依赖于它们的形状和所处的空间维度。

本文针对位置大数据的高效密度聚类问题，从密度聚类算法内在特性优化视角，实现和改进了一种新型的快速查找密度山峰的聚类算法，可在一般计算平台实现堆位置数据的密度蕨类任务。最后,在多个基准测试数据和大规模位置大数据上进行了综合实验评估,验证了所提算法的聚类效果和性能,并与常用聚类方法进行了对比分析。

二、算法描述

1、基本思想

该算法基于以下假设：簇的中心被局部密度更小的点包围而且这些中心距离局部密度比它们大的点相当远。即聚类中心同时具有以下两个特点：

本身的密度大，即它被密度均不超过它的邻居包围；
与其它密度更大的数据点之间的“距离”相对更大；

对于每一个数据点，我们计算两个数值，：局部密度Pi 以及它与更高密度点的距离δi .这两个值只依赖与数据点之间的距离，我们假定该距离满足三角不等式。我们将采用Cut-off kernel和Gaussian kernel两种计算方式数据点的局部密度定义为：

其中，dc 为截断距离, ρi为到点i 的距离比dc 小的点的数目。该算法支持在不同的点出Pi 的相对大小敏感。这就意味着对于大型数据集，该分析的结果对于dc 的选择鲁棒性很强。

δi 是点i 到与比它密度高的点的最小距离：

对于密度最高的点，我们一般记为δi=maxj(dij) 。我们注意到局部或者全局密度最高的点的δi 会比周围的邻居点大。因此，那些δi 的值异常答的点被视作簇的中心。

2、聚类方法

2.1 找出聚类中心

如图1[11]所示，所有点的密度值按照由高到低排列，“1”表示密度最高的点。B图中画图每个点和的函数关系，从中可以看出“9”和“10”号点拥有相近的密度值但是其δi 不同，这里“9”属于“1”号类别。“26”，“27”和“28”号点有一个相对较大的δi ，但是其ρi 太小，这主要是因为它们是孤立点。

通过生成的Decision-Graph(delta-rho)图，根据前面的结论，δi 和Pi 值均比其它大的点就是我们要找的聚类中心，然后根据此决策图给定δmin 和pmin ，筛选出同时满足()和()条件的点作为聚类中心点。

确定聚类中心，并初始化数据点归类属性标记,具体为：

通过生成的Decision-Graph(delta-rho)图，根据前面的结论，δi 和Pi 值均比其它大的点就是我们要找的聚类中心，然后根据此决策图给定δmin和pmin，筛选出同时满足(ρi>ρmin )和(δi>δmin )条件的点作为聚类中心点。

确定聚类中心，并初始化数据点归类属性标记,具体为：

对于那些在决策图中无法用肉眼判断出聚类中心的情形，作者在文中给出了一种确定聚类中心个数的提示：计算一个将ρ 值和δ 值综合考虑的量

显然，值越大，越有可能是聚类中心。因此，只需对进行降序排列，然后从前往后截取若干个数据点作为聚类中心。

将排序后的在坐标平面画出来，如图2所示。由此可见L非聚类中心的值比较平滑，而从非聚类中心过渡到聚类中心时，值有一个明显的跳跃、这个跳跃用肉眼或数值检测都可以判断出来。作者还提到，对于人工随机生成的数据集，的分布还满足幂次定律，即近似呈直线，且斜率依赖于数据维度。

2.2 剩余点的类别指派

当聚类中心确定之后，剩下的非聚类中心点的类别标签指定按照以下原则：

当前点的类别标签等于高于当前点密度的最近的点的标签一致。从而对所有点的类别进行了指定。如图3[11]所示，编号表示密度高低，“1”表示密度最高，以此类推。“1”和“2”均为聚类中心,"3"号点的类别标签应该为与距离其最近的密度高于其的点一致，因此“3”号点属于聚类中心1，由于“4”号点最近的密度比其高的点为“3”号点，因此其类别标签与”3“号相同，也为聚类中心1。

3．去除噪声点

2.3 类别间边界确定

一个簇中的数据点可以分为簇核心部分和簇的光环（作为噪声点）。在对每一个点指派所属类别之后，这里文章没有人为直接用噪音信号截断的方法去除噪音点，而是先算出类别之间的边界，边界区域由这样的数据点构成：它们本身属于该簇，但在与其距离不超过dc的范围内，存在属于其它簇的数据点。

2.4 划分出簇的光环

利用边界区域，这个簇就可以计算出一个平均局部密度ρb，密度值小于该值的点则被分为光环部分，即噪声点。如图4所示，橙色圈内的为簇核心部分，外的数据点为光环。

4.算法改进和优化

4.1局部密度定义

在原文中局部密度Pi 使用的是Cur-off kernel的计算方式，在实际代码实现中我们采用的是Gaussian kernel的计算方式

对比定义（1.1）和（1.2）易知，Cut-off kernel为离散值，Gaussian kernel为连续值，因此，相对来说。后者产生冲突（即不同的数据点具有相同的局部密度值）的概率更小。此外，对于（1.2）仍满足原定义。

4.2 聚类中心个数的确定

利用ρ 和δ 定义更合适的函数，鉴于这两个值可能处于不同的数量级。因此，首先对做一次归一化，都归一到[0,1]区间，然后定义 =，画出曲线图，非聚类中心和聚类中心的点会产生一个突然的跳跃，这样就可自动判断出聚类中心，测试结果如图5所示：

4.3 噪声点的去除
光晕（噪声点）即指的是较大而较小的点，通过曲线图即可以判断出噪声点，效果如图6所示：

4.4 距离计算方式

在作者提供的matlab源代码实现过程中计算距离采用的是欧氏距离的计算方式，由于所采用的的训练测试数据集为二维、三维或者通过降维处理之后的数据，因此采用欧氏距离没有很大的问题。但是如果将该算法应用到实际中，例如图像、音视频识别等具有更高维的数据时欧氏距离将不再使用。

算法Python复现：

参考了一个博主的代码并进行了一定改动。这篇论文的原作者个人主页上有matlab实现的代码。

# -*- coding: utf-8 -*-
"""
Created on Mon Oct 20 13:38:18 2014

@author: Shaobo
The best of all
"""

import numpy as np
from numpy import *
import matplotlib.pyplot as plt
import random
import math



MAX = 1000000

def nearestNeighbor(index):
    dd = MAX
    neighbor = -1
    for i in list(range(length)):
        if dist[index, i] < dd and rho[index] < rho[i]:
            dd = dist[index, i]
            neighbor = i
    if result[neighbor] == -1:
        result[neighbor] = nearestNeighbor(neighbor)
    return result[neighbor]

#Read data
fileName = input("Enter the file's name:")
location = []
label = []
for line in open(fileName, "r"):
    items = line.strip("\n").split(",")#按“，”来分割
    label.append(int(items.pop()))
    tmp = []
    for item in items:
        tmp.append(float(item))
    location.append(tmp)
location = np.array(location)
label = np.array(label)
length = len(location)

#Caculate distance计算距离
dist = np.zeros((length, length))
ll = []
begin = 0
while begin < length-1:#循环数据
    end = begin + 1
    while end < length:
        dd = np.linalg.norm(location[begin]-location[end])
        dist[begin][end] = dd
        dist[end][begin] = dd
        ll.append(dd)
        end = end + 1
    begin = begin + 1
ll = np.array(ll)

# Algorithm        
#percent = float(input("Enter the average percentage of neighbours: "))
percent = 2.0
position = int(len(ll) * percent / 100)
sortedll = np.sort(ll)
dc = sortedll[position] #阈值

#求点的局部密度(local density)
rho = np.zeros((length, 1))
begin = 0
while begin < length-1:
    end = begin + 1
    while end < length:
        rho[begin] = rho[begin] + math.exp(-(dist[begin][end]/dc) ** 2)
        rho[end] = rho[end] + math.exp(-(dist[begin][end]/dc) ** 2)
        #if dist[begin][end] < dc:
        #    rho[begin] = rho[begin] + 1
        #    rho[end] = rho[end] + 1
        end = end + 1
    begin = begin + 1

#求比点的局部密度大的点到该点的最小距离
delta = np.ones((length, 1)) * MAX
maxDensity = np.max(rho)
begin = 0
while begin < length:
    if rho[begin] < maxDensity:
        end = 0
        while end < length:
            if rho[end] > rho[begin] and dist[begin][end] < delta[begin]:
                delta[begin] = dist[begin][end]
            end = end + 1
    else:
        delta[begin] = 0.0
        end = 0
        while end < length:
            if dist[begin][end] > delta[begin]:
                delta[begin] = dist[begin][end]
            end = end + 1
    begin = begin + 1

rate1 = 0.1
#Aggregation Spiral 0.6
#Jain Flame 0.8
#D31 0.75
#R15 0.6
#Compound 0.5
#Pathbased 0.2
#T 0.008 0.1
thRho = rate1 * (np.max(rho) - np.min(rho)) + np.min(rho)

rate2 = 0.016
#Aggregation Spiral 0.2
#Jain Flame 0.2
#D31 0.05
#R15 0.1
#Compound 0.08
#Pathbased 0.4
#T 0.028 0.016
thDel = rate2 * (np.max(delta) - np.min(delta)) + np.min(delta)

#确定聚类中心
result = np.ones(length, dtype=np.int) * (-1)
center = 0
#items = range(length)
#random.shuffle(items)
for i in list(range(length)): #items:
    if rho[i] > thRho and delta[i] > thDel:
        result[i] = center
        center = center + 1
#赋予每个点聚类类标
for i in list(range(length)):
    dist[i][i] = MAX

for i in list(range(length)):
    if result[i] == -1:
        result[i] = nearestNeighbor(i)
    else:
        continue

plt.plot(rho, delta, '.')
plt.xlabel('rho'), plt.ylabel('delta')
plt.show()

R = list(range(256))
random.shuffle(R)
R = np.array(R)/255.0
G = list(range(256))
random.shuffle(G)
G = np.array(G)/255.0
B = list(range(256))
random.shuffle(B)
B = np.array(B)/255.0
colors = []
for i in list(range(256)):#256
    colors.append((R[i], G[i], B[i]))#绘制点颜色

#——————————————————————
fig=plt.figure()
rect=[0.1,0.1,0.8,0.8]
axprops=dict(xticks=[],yticks=[])
ax0=fig.add_axes(rect,label='ax0',**axprops)
imgp=plt.imread('bg.png')
ax0.imshow(imgp)

ax1=fig.add_axes(rect,label='ax1',frameon=False)
#------------------------------------------

for i in list(range(400)):#256
    index = result[i]
    ax1.scatter(location[i][0], location[i][1], color = colors[index], marker = '.')
plt.xlabel('x'), plt.ylabel('y')
plt.show()

#——————————————————————

fig=plt.figure()
rect=[0.1,0.1,0.8,0.8]
axprops=dict(xticks=[],yticks=[])
ax0=fig.add_axes(rect,label='ax0',**axprops)
imgp=plt.imread('bg.png')
ax0.imshow(imgp)
#---------------------------------------

ax1=fig.add_axes(rect,label='ax1',frameon=False)

for i in list(range(length)):#256
    index = label[i]
    ax1.scatter(location[i][0], location[i][1], color = colors[index], marker = '.')
plt.xlabel('x'), plt.ylabel('y')
plt.show()

万万VV

关注

1
点赞
踩
10

收藏

觉得还不错? 一键收藏
0
评论
面向地理大数据的新型密度山峰聚类算法（一）——算法描述

摘要通过实现和改进在《Science》杂志上发表的一种新型快速密度山峰聚类算法[1],与现有的常用聚类方法进行分析比较，发现其优点和不足的地方，并实际应用于地理定位大数据的分析中，以快速发现位置数据中任意形状的聚类簇模式和噪声。首先，基于类簇中心被具有较低局部密度的邻居点包围，且与具有更高密度的任何点有相对较大的距离这样的假设提出基本算法。其次，通过Python实现该算法并采用基...
复制链接

扫一扫