前言
DESCAN是一种密度聚类算法,其主要思想是密度和传销。本文不做详细的术语讲解。而是通过案例进行讲解,笔者认为这种方式更能通俗易懂。如果想看严谨的原理解释,那本文可能并不适合你,请移步其他博文。
算法流程
①选定一个半径r和最小密度点数MinPoint
②随机选择数据中的一个点记作核心点。
③从该核心点计算到其他点的距离(比如欧氏距离),选出其中距离小于半径r的点。
④如果这些点的数量(包括核心点)大于等于我们设定的MinPoint,我们就说这些点是密度可达点。并将这些密度可达点记作与核心点同一个簇(类别)。
⑤将得到的密度可达点作为新的核心点。
⑥循环③、④、⑤步,直到不再能够找到密度可达点。
⑦第⑥步结束后开启新一轮的循环,即重复②-⑥步,直到所有的点都被设定为核心点,也就是所有的点都计算完。
实现案例
先来看一下图:
上面每个小圆代表数据。我们的目标就是要将这些数据聚成两类。
根据步骤,我们选定一个半径为r=1,Minpoint=3。
第二步随机选取一个点作为核心点,假设选中x1,那么计算距离
可以看到,计算距离之后,对于所有的点,距离小于r的有{x0、x1、x3、x5}。此时的距离小于r的数目为4,因为4>MinPoint。所以此时可以说{x0、x3、x5}(除去x1核心点)是密度可达点。并且{x0、x3、x5}和核心点x1归为一簇,即{x0、x1、x3、x5}为一簇。
接下来就将{x0、x3、x5}作为新的核心点,分别计算其是否密度可达点。最终算出来
x0∈{x0、x1、x3、x5},点数(4)大于MinPoint,所以可以将{x1、x3、x5}与x0归为一簇,因为x0是从x1的密度可达点,所以x0与x1必须是同一类。又因为{x0、x1、x3、x5}包含了x0里的{x0、x1、x3、x5},所以不做修改。
x5∈{x0、x1、x5},点数(3)等于Minpoint,所以可以将{x0、x1}与x5归为一簇,因为x5是从x1的密度可达点,所以x5与x1必须是同一类。又因为{x0、x1、x3、x5}包含了x5里的{x0、x1、x5},所以不做修改。
对于x3
x3∈{x0、x1、x3、x4},点数(4)大于MinPoint,所以可以将{x0、x1、x4}与x3归为一簇,因为x3是从x1的密度可达点,所以x3与x1必须是同一类。又因为{x0、x1、x3、x5}没有包含{x4},所以修改原始的簇为{x0、x1、x3、x4、x5}。
x3中包含有新的密度可达点x4。那么就将x4作为新的核心点。计算后结果也没有新的密度可达点。所以将{x0、x1、x3、x4、x5}归为一簇。
接下来再随机选择一个没有被遍历过的点,假设选中x11。
发现其仅仅包含一个数据点。数据点数目小于3,所以不做分类。
因此直接进行下一步。随机选中一个没有被迭代过的点作为核心点。假设选中x10
可以看到点数等于Minpoint,因此将其余作为密度可达点继续迭代。最终得到x{6、7,8,9,10}作为一簇。
最后完成迭代,结束算法。
代码实现
注意两点。
①因为要遍历所有的点,因此代码中用self.is_check代表是否某个点被遍历过。
②代码计算距离的公式用欧氏距离
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
import imageio
plt.ion()
a=0
fig,ax=plt.subplots()
image=[]
class DESCAN():
def __init__(self,x):
self.x=x
#假定数据对应的标签y为-1,训练过程会更改
self.y=np.ones(shape=(self.x.shape[0],1))*-1
#计算某个点到其他点的距离
self.distance=np.ones(shape=(self.x.shape[0],self.x.shape[0]))
for index,i in enumerate(self.x):
#比如第i行代表第i个数据到其他点的距离。
self.distance[index,:]=np.linalg.norm(self.x-i,ord=2,axis=1).reshape(-1)
def fit(self,r,numpoit):
self.r=r#半径
self.numpoint=numpoit#最小点数
# 生成一个矩阵用来记录该点是否被遍历过,事先是定全为1,被遍历过后设置为0.
self.is_check = np.ones(shape=self.x.shape[0])
class_num=0#标签值
#设置循环直到所有的数据都被遍历过
while True:
# 计算概率。因为被遍历过的数据点对应的is_check会被置为0,也就是概率为0,不会让下面随机选的时候选到它
p=self.is_check/sum(self.is_check)
# 随机选择一个点,返回对应索引,p是用来排除已经被计算过的点
index=np.random.choice(self.x.shape[0],1,p=p)
self.is_check[index]=0#令被选中的点更改为已被遍历过
self.cycle(index,class_num,True)#代入循环式
if (self.is_check == 0).all(): # 判断所有的数据是否都被遍历过,是则跳出循环
plt.ioff() #绘图函数,与模型无关
plt.show() #绘图函数,与模型无关
break
# 判断是否标签值y中是否有新的类别。意为如果self.y中只要包含有一个我们传过去的class_num,说明有满足条件的簇.
#因此我们就要寻找下一个簇。因此class_num+=1
if (self.y==class_num).any():
class_num+=1
# imageio.mimsave("result.gif",image) #保存图像函数,与模型无关
def cycle(self,index,class_num,First=False):
'''
:param index: 核心点的索引
:param class_num: 类别标签
:param First: 用于判断是从fit函数过来的,还是从cycle循环本身循环过来的
:return:
'''
#因此index在第二次循环的时候会有可能有多个index,因此作循环。
for i in index:#迭代选取一个核心点的索引,寻找密度可达点
# 取出一个核心点
kernel_point = self.x[i, :]
# 取出该核心点和其他点的距离。也就是第i行
distances=self.distance[i,:].squeeze()
# 找出距离小于r且没有被遍历过的数据点。如果没有没被遍历过的点那就没有意义了。
new_index = np.argwhere(((distances <= self.r) & (self.is_check == 1)))
plot_figure(self.x,self.y,kernel_point,self.r)#绘图函数,与模型无关
# 判断下一步距离小于r的数据点数能否达到我们预设的最少点数
if new_index.size != 0: # 判断是否有没有被遍历过的数据点,有的话则将该这些点作为新的核心点。
# 找出距离小于r的点的对应的索引,然后判断是否设定的最小点数。
point_index = np.argwhere((distances <= self.r))
if point_index.shape[0]>=self.numpoint:
############################
# 将这些点赋给对应的标签值,因此此处只更新新的没有遍历过的点。但是一开始从fit传过来的点没有赋予标签值。
# 因此用First判断,如果是从fit函数过来的,就要赋给标签值进行额外更新。
if First:
self.y[index, :] = class_num
self.y[new_index] = class_num
############################
self.is_check[new_index]=0#设定这些核心点为被计算过
self.cycle(new_index,class_num)#将这些核心点放到下一轮循环找出其他的直接密度可达点
def plot_figure(x,y,kernel_point,r):
global a
'''绘图函数,与模型无关'''
map_color={-1:"k",0:"r",1:"g",2:"b",3:"y",4:"c"}
color_initial=[map_color[i] for i in y.squeeze()]
a+=1
kernel_point = kernel_point.squeeze()
plt.cla()
plt.scatter(x[:,0],x[:,1],c=color_initial)
circle=plt.Circle((kernel_point[0],kernel_point[1]),radius=r,fill=False)
plt.axis("equal")
ax.add_artist(circle)
# plt.savefig("1.png") #用于保存图像
# image.append(imageio.v3.imread("1.png")) #用于保存图像
plt.pause(0.01)
def main():
#生成第一类数据g
x1=stats.norm.rvs(0,1,(100,2))
#生成第二类数据r
x2=stats.norm.rvs(10,1,(100,2))
#生成第三类数据b
x3=stats.norm.rvs(20,1,(100,2))
#合并数据
x=np.concatenate([x1,x2,x3],axis=0)
model=DESCAN(x)
model.fit(1,3)
if __name__ == '__main__':
main()
结束语
文章并不严谨。有任何错误的地方。还望指出。阿里嘎多。