ST-CFSFDP：快速搜索密度峰值的时空聚类算法

王培晓,张恒才,王海波,吴 升

1. 福州大学数字中国研究院(福建)，福建 福州 350002； 2. 中国科学院地理科学与资源研究所资源与环境信息系统国家重点实验室，北京 100101； 3. 海西政务大数据应用协同创新中心，福建 福州 350002； 4. 湖北工业大学经济与管理学院，湖北 武汉 430068

近年来，室内外定位技术迅猛发展，如北斗、GPS、WiFi、蓝牙、UWB(ultra-wideband)等，海量移动终端不断普及，如PDA、智能手机、平板电脑等，网络技术不断进步，室内外位置服务应用不断增多，如在线导航、基于位置的社交网络、基于位置的广告推送、商业物流调度与管理等，时空轨迹数据爆发式增长[1-2]，为人类出行模式及人类移动性、城市计算、社会计算、交通管理、城市规划、人口流动监测、公共安全保障等研究提供了重要支撑[3-4]。

聚类算法是发现时空模式、探寻移动规律的关键技术之一，旨在将实体划分为一系列具有一定分布模式的簇集，同一簇集中的实体具有较大的相似度，不同簇集中的实体具有较大差别[5-8]，已广泛应用于犯罪热点分析[9]、地震空间分布模式挖掘[10]、制图自动综合[11]、遥感影像处理[12]、公共设施选址[13]、地价评估[14]、用户停留区域识别[15]等研究。

目前，聚类算法大致分为[6,16]：基于划分、基于层次、基于密度、基于格网的聚类算法。基于划分聚类算法使用聚类中心(位于簇集中心附近的一个对象)来表示每个簇集，如经典的k-means[17]算法和K-Medoid[18]算法，具有计算逻辑简单、运算效率高等优点；基于层次聚类算法使用逐步合并或分裂的策略来实现聚类，如改进的CURE[19]算法，不仅能发现球形的空间簇集[20-21]，也能够发现较为复杂结构的空间簇集，但过多的超参数增加了算法的不确定性；基于密度聚类算法，如DBSCAN[22]将簇定义为密度相连的点的最大集合，但由于采用固定阈值聚类，难以适应空间实体密度的变化[23]，OPTICS[22,24]算法为聚类分析生成一个增广簇排序，解决了DBSCAN算法全局密度阈值的缺点；基于格网的聚类算法[25]其聚类效果依赖于网格的大小，如果网格划分过细，则时间复杂度会显著增加，如果网格划分过粗，则聚类质量难以得到保证。

此外，许多研究学者在经典聚类算法基础之上，提出了时空聚类算法，文献[26—27]在基于密度聚类的基础上提出了时空密度聚类ST-DBSCAN和ST-OPTICS算法，在时空密度聚类中，使用时间距离将空间邻域扩展到时空邻域，从而寻找时空邻域下密度相连的区域。文献[10,28]采用时空距离(时间距离与空间距离的函数)度量任意两个时空事件的差异，后将时空距离应用到传统的聚类算法中挖掘时空事件的模式；文献[29—30]在空间扫描统计的基础上提出了时空扫描统计算法，通过滑动扫描窗口(空间半径和时间间隔定义的圆柱体)揭示时空事件随时间聚类模式；文献[31]通过窗口距离与时空最近邻的概念重新定义了时空密度，从而提出了STSNN算法；文献[32]从统计学的视角提出了WNN算法，将时空邻域中的实体分为特征与噪声两个时空泊松点过程，使用密度相连(density-connective)的概念将特征分为多个簇集；文献[33—34]分别在基于划分和基于层次聚类基础上提出了WKM和ST-AGNES算法，将其分别应用于人类行为和用户停留区域识别中，并取得了较好的聚类效果。

经典的CFSFDP[35](clustering by fast search and find of density peaks)算法结合了基于划分和基于密度聚类算法的优点，不仅可以发现多密度任意形状的簇集，还可以通过决策图辅助识别聚类中心的个数。但该算法无法很好地应用于时空数据聚类研究之中，究其原因在于CFSFDP算法未考虑时间约束，如图1所示，并没有正确识别簇集。如错误地将簇集A(t5,t6,t7,t8)、C(t15,t16,t17,t18)合并为一个簇集，且并不能发现簇集之间的顺序关系。此外，在单簇集中存在多个密度峰值点时，该算法将会产生错误的时空聚类结果。

图1 是否考虑时间约束的两种聚类结果Fig.1 Comparison of clustering results with and without time constraints

鉴于此，本文提出一种快速搜索密度峰值的时空聚类算法(spatial-temporal clustering by fast search and find of density peaks，ST-CFSFDP)，在CFSFDP算法的基础上引入了时间约束，修改了样本属性值的计算策略。采用模拟数据和真实数据案例进行算法有效性的验证。模拟数据试验结果表明，与CFSFDP算法相比，ST-CFSFDP算法不仅可以克服单簇集中可能存在多密度峰值的不足，且可以区分并识别相同位置不同时间的簇集。以移动对象轨迹中停留点识别为例，ST-CFSFDP算法较经典的ST-DBCSAN、ST-OPTICS和ST-AGNES算法识别正确率分别可提高5.2%、4.2%和7.6%。

1 ST-CFSFDP算法

1.1 CFSFDP算法

ρi=∑jχ(dij-dc)

(1)

(2)

式中，dij为样本点i和样本点j之间的空间距离；dc是人为指定的超参数，称为截断距离；χ(x)表示0～1函数，当x<0时，χ(x)=1,当x≥0时，χ(x)=0，即局部密度ρi表示落在以样本点xi为圆心，截断距离dc为半径的圆形区域中样本点个数；δi表示局部密度大于xi的所有样本点中，距离xi最近的样本点与xi之间的距离,由于式(2)无法计算密度最大的样本点，因此密度最大的样本点δi=max (dij)。

CFSFDP算法通过各样本点的局部密度ρi和距离δi确定聚类中心(ρi和δi均相对较大的样本点)。如图2(a)所示，对于数据集X中每一个样本点xi，计算其局部密度ρi和距离δi，生成相应的三元组(ρi,δi,λi)，λi表示δi对应的节点编号，用于非聚类中心样本点的类别分配，使用ρi和δi构造数据集X的决策图(decision graph)，如图2(b)所示。容易发现，标号为①和⑩的样本点同时具有较大的ρ值和δ值，因此被确定为数据集X的两个聚类中心。为了更加方便地确定数据集X中的聚类中心，文献[35]提出了一个综合考虑ρ和δ的量γ。

γi=ρiδi

(3)

(4)

1.2 基于时空约束的ST-CFSFDP聚类算法

ST-CFSFDP算法首先在CFSFDP算法的基础上引入了第2个超参数tc，针对时间序列数据集X，修改了ρ与δ的计算策略。其中时间约束下的局部密度ρ计算方法如式(5)所示

ρi=∑jχ(dij-dc,tij-tc)

(5)

(6)

式中，tqiqj表示样本点xqi与样本点xqj的时间距离，即修改后的δ根据样本点的局部密度降序依次计算，即使单簇集中存在多个密度峰值，依旧仅有一个样本点被选做聚类中心。然而当单簇集的时间跨度过长时，仅使用时间距离计算δ同样会将单簇集分裂成多个小簇集。如图4(b)所示，假设簇集C内部的样本点4、k具有较高的局部密度ρ，由于两样本点的时间跨度较大，两点都将具有较大的δ，依据本文的计算策略γi=ρiδi，簇集C将产生两个聚类中心。为解决这个问题，本文添加了聚类中心合并的过程，将识别的聚类中心按时间顺序链接，分别计算相邻两聚类中心的距离，若小于距离邻域dc，将相邻两个聚类中心合并，经聚类中心合并后，聚类中心的个数将与簇集的个数一一对应。

在ST-CFSFDP算法中，剩余样本点的分配沿时间轴进行。如图5所示，针对时间顺序分布的数据集X，首先识别数据集中的聚类中心(红色样本点为合并后的聚类中心)，然后每个聚类中心沿时间轴向前、后搜索本簇集中的剩余样本点，以聚类中心CP2向前搜索为例，寻找CP2前一时刻的样本点t15，计算样本点与CP1和CP2的距离，如果距离CP2近，则将样本点t15归属于CP2，否则聚类中心CP2向前搜索完毕，以同样的策略向后搜索最终确定该簇集包含的所有样本。所有聚类中心依次搜索完毕后，未被分配类别的样本点为离群点，以图5黄色样本点为例，聚类中心CP2延时间轴向后搜索时，由于黄色样本点OP1距离聚类中心CP3较近，因此聚类中心CP2向后搜索结束，当聚类中心CP3延时间轴向前搜索时，黄色样本点OP2距离聚类中心CP2较近，因此聚类中心CP3向前搜索结束，此时黄色样本点OP1与OP2之间的样本点为离群点。ST-CFSFDP算法整体实现流程如下。

算法 快速搜索密度峰值的时空聚类算法

输入：序列数据sequence=={ti,pi},距离阈值disthr,时间阈值tthr

输出：每一个样本点的聚类类别labels={ci}

functionST_CFSFDP (sequence, disthr,tthr)

∥计算样本点的空间距离矩阵

stDisMat=computeSpatialDisMat(sequence)

∥计算样本点的时间距离矩阵

timeDisMat=computeTimeDisMat(sequence)

∥计算样本点在时空约束下的局部密度

densityArr=computeDensity(stDisMat,disthr,

timeDisMat,tthr)

∥获得局部密度降序排列的下标序列

densitySortArr=argsort(densityArr)

∥初始化数组，用于存放每个轨迹点的属性δ

closestDis=[]

fori=0→len(densitySortArr)do

∥获得当前样本点的索引

node=densitySortArr[i]

∥获得比当前样本点局部密度更大的索引集合

nodeIdArr=densitySortArr[i+1;]

∥计算每一个样本点的属性δ

closestDis[node]=compute(node, nodeIdArr,

stDisMat,timeDisMat)

endfor

∥计算每一个样本点的属性γ

gamma=closestDis*densityArr

∥通过决策图得到聚类的数目

classNum=getNumFromDecisionGraph(gamma)

∥获取样本点的聚类类别

labels=extract_clustering(gamma,classNum,

stDisMat,timeDisMat)

returnlabels

endfunction

1.3 算法时间复杂度分析

ST-CFSFDP算法的时间复杂度定性描述了该算法的运行时间。ST-CFSFDP算法的时间复杂度主要由两部分组成：计算样本点的局部密度ρ和距离δ，其中ρ和δ的计算都涉及样本点之间的距离计算。在数据集规模为n的情况下，①计算任意两个样本点之间的距离，时间复杂度为O(n2)；②计算n个样本点的局部密度δ，时间复杂度为O(n2)；③计算n个样本点的距离δ，时间复杂度同样为O(n2)。因此，ST-CFSFDP算法的时间复杂度为O(n2)。

图2 数据集X及其决策图Fig.2 Dataset X and its decision graph

图3 递减顺序排列的γFig.3 The value of γ in decreasing order

图4 CFSFDP算法改进Fig.4 Improvement of CFSFDP algorithm

图5 ST-CFSFDP算法剩余点分配Fig.5 Distribution of non-clustered center points of ST-CFSFDP algorithm

2 试验结果与讨论

2.1 试验数据

为分析ST-CFSFDP算法的性能，分别采用模拟数据集和真实用户轨迹进行试验。采用模拟数据集进行试验的主要目的是以二维图形的方式直观地展示聚类结果，采用真实用户轨迹进行试验的主要目的是使用相关指标定量地评价聚类算法的性能。

本文共模拟4组数据量不同的时序数据集X。每组模拟数据集的生成过程如下：首先生成指定的具有时间约束的几何形状，然后根据指定的几何形状等时间间隔采样而成，其中一组数据如图6(a)、(b)所示。

真实数据来源于济南市某广场移动用户的WiFi定位数据。用户移动定位数据是一种典型的时空数据，在一定程度反映了用户的个人生活习惯和日常行为模式。时空聚类算法常用于从用户移动定位数据识别用户的停留区域，从而进一步挖掘用户的兴趣爱好，因此本文将ST-CFSFDP算法应用于停留区域识别并定量的评价该算法的性能。移动用户定位数据覆盖广场5个楼层，平均采样为1～10 s不等，定位精度约为3 m。由于定位数据难以获取用户停留区域的标注数据，因此本文采用爬虫技术从百度地图上获取了该广场的商铺数据，借助商铺数据标注用户轨迹中的停留信息。标注过程如下：首先对商铺数据与蓝牙定位数据求交，然后统计某用户在商铺内部的持续停留时间，若用户在商铺内部的持续停留时间超过一定阈值，则将相应的轨迹点标注为停留点。依据上述规则，本文共标注472条用户轨迹(当用户轨迹中停留区域过少时，本文认为该轨迹价值不高，删除该轨迹)，轨迹点总记录量为935 639个。经标注后的某用户轨迹数据如表1所示，每条记录包括用户唯一标识ID、记录时间、用户位置(X、Y坐标及所在楼层ID)、停留的商铺ID(-1代表用户未停留)。

图6 模拟数据集Fig.6 Simulated data sets

表1 单用户轨迹数据实例

2.2 算法评价指标及对比试验选择

本文以准确率(accuracy)和召回率(recall)作为停留区域识别的定量评价指标。用户停留区域可理解为时间约束下用户轨迹中的簇集，进一步可抽象为[36]S=(userID,Lng,Lat,arrT,levT)，(Lng,Lat)表示停留区域内的某个点坐标，其应最大概率出现在停留区域内部，一般为区域内所有轨迹点的均值坐标，本文为聚类中心坐标；arrT表示用户抵达停留区域的时间；levT表示用户离开停留区域的时间。停留区域识别结果是否正确主要依赖3个方面：①识别的停留区域和标注的停留区域数量是否一致；②识别的停留区域(Lng,Lat)坐标是否处于标注商铺内部；③识别的停留区域起止时间(arrT,levT)是否与标注数据一致。结合上述3个方面，本文accuracy和recall计算方法如式(7)、式(8)所示。

(7)

(8)

式中，SClabel表示标注的停留区域个数；SCalgorithm表示算法识别的停留区域个数；SCcorrect表示算法判断正确的停留区域个数，本文采用SO和TO[15]两个函数共同判断某个停留区域是否识别正确(SO函数判断识别结果和标注数据在空间上是否邻近，TO函数判断识别结果和标注数据在时间上是否重合)。本文将ST-DBSCAN、ST-OPTICS、ST-AGNES、ST-CFSFDP时空聚类算法进行对比，探讨4种时空聚类算法在不同阈值下准确率和召回率的变化情况。

2.3 模拟数据集上的测试结果分析

CFSFDP与ST-CFSFDP算法在模拟数据集上的聚类结果如图7(e)、7(f)所示：①对于区域A，由于存在两个聚类中心，使用CFSFDP算法，分裂成两个小簇集，产生了错误的聚类结果，而使用ST-CFSFDP算法，依据改进的属性计算策略，只聚类成一个簇集，从而得到更合理的聚类结果；②对于区域B，由于CFSFDP算法未考虑数据的时间约束，因此只聚类出一个簇集，而ST-CFSFDP算法考虑数据的时间约束，所以可以聚类出相同位置但不同时间的两个簇集；③对于区域C，使用ST-CFSFDP算法识别出两个聚类中心，但仅产生一个簇集，原因是其中一个小簇集时间跨度过长，在该簇集中识别出两个聚类中心，但ST-CFSFDP算法最终会将这两个聚类中心合并(相邻聚类中心的距离小于距离阈值dc)，因此最终只产生一个簇集。综上所述，与CFSFDP算法相比，ST-CFSFDP算法不仅克服了单簇集中可能存在多密度峰值的不足，且实现了考虑时间约束的时空聚类。ST-CFSFDP算法与ST-DBSCAN、ST-OPTICS、ST-AGNES算法的聚类结果如图8所示，进一步说明ST-CFSFDP算法具有发现时空簇集的能力。

由上文分析可知，ST-CFSFDP算法具有良好的聚类性能，不仅解决了单簇集存在多密度峰值的问题，还能正确发现时空数据集的簇分布特征。进一步分析算法的运行时间，表2为两种聚类算法在4组模拟数据集上运行时间的对比。试验结果可以看出，ST-CFSFDP算法的时间开销要略大于CFSFDP算法，这主要是因为ST-CFSFDP算法在计算样本点局部密度时增加了时间开销。但是这两种算法的时间复杂度的均为O(n2)。与ST-DBSCAN、ST-OPTICS、ST-AGNES算法的运行时间相比，ST-CFSFDP的运行时间较小，能够较快地得到聚类的运行结果。

表2 时空聚类算法运行时间在模拟数据集对比分析

Tab.2 The running time of spatio-temporal clustering on simulated data sets

数据集X记录数CFSFDPST-CFSFDPST-DBSCANST-OPTICSST-AGNES3140.00450.00480.0052 0.03120.003334540.21590.29820.80581.97230.232463630.68390.84203.122110.6710.7918156234.00945.879115.627626.5154.8917

2.4 真实数据集上的测试结果分析

为评价ST-CFSFDP算法性能，本文参考室内房间宽度将初始距离阈值设置为1 m，并以0.5 m步长递增；初始时间阈值设置为50 s，以20 s步长递增；以启发式的方法[22]设置ST-DBSCAN、ST-OPTICS算法中minPts参数，minPts=ln(N)，N为某用户轨迹点数量。算法识别结果的准确率对比如图9所示。4种算法的识别准确率变化趋势在整体上较为相似，在时间阈值固定的情况下，准确率均随着距离阈值的增加呈现先增加后下降的趋势。综合4种算法在不同阈值下的表现可以发现：①ST-CFSFDP算法对参数敏感度低，ST-DBSCAN与ST-OPTICS算法对参数敏感度较高(准确率随着超参数的改变波动较大)；②相较ST-DBSCAN、ST-OPTICS、ST-AGNES算法，ST-CFSFDP算法准确率较高，在时间阈值90 s时最为明显，ST-CFSFDP算法准确率最高可达82.4%，与现有算法相比，最高可提升7.6%；③ST-DBSCAN、ST-OPTICS算法本质上是一种算法，所以两种算法的准确率高度一致；④ST-AGNES算法隐含时间约束，仅需要距离阈值即可完成聚类，因此算法准确率不受时间阈值的影响。

算法识别结果的召回率对比如图10所示。4种算法的召回率均随距离阈值增加呈现先增加后下降的趋势，与算法准确率一致。当时间阈值为90 s时，4种算法的召回率较高，原因在于人工标注时用户的停留时间阈值与90 s较接近。在此阈值下4种算法将会得到较精确的识别结果。与此同时ST-CFSFDP算法的召回率略高于其余3种算法且准确率与召回率波动最小，原因为ST-CFSFDP算法采用聚类中心作为用户最可能的停留位置，而聚类中心作为局部范围内密度最大的点，受阈值影响较小。本文将4种算法最高准确率相应的运行时间进行对比，结果如表3所示。ST-CFSFDP算法与ST-AGNES算法相比，ST-CFSFDP算法的运行时略高，但算法的识别正确率却可提升7.6%；与ST-DBSCAN算法相比，ST-CFSFDP算法的运行时间略有降低且识别准确率提升了5.2%；与ST-OPTICS算法相比，ST-CFSFDP算法的运行时间得到了大幅度提升且识别准确率也有一定程度的提高。

表3 4种算法运行时间对比分析

3 结论与展望

CFSFDP算法是一种新颖的空间聚类算法，通过计算样本点的属性值γ，能够快速发现数据集中的密度峰值点。然而当数据集的某簇集存在多密度峰值点时，该算法易产生错误的聚类结果，且CFSFDP算法无法实现考虑时间约束的时空聚类。针对以上两点不足，本文提出了时空聚类算法ST-CFSFDP。ST-CFSFDP算法在CFSFDP算法基础上加入时间约束，并修改了样本属性值的计算策略。为验证算法的有效性，首先采用模拟数据进行定性试验。结果表明，与原有的CFSFDP算法相比，ST-CFSFDP算法不仅可以克服单簇集中可能存在多密度峰值的不足，且可以区分并识别相同位置不同时间的簇集。最后本文将ST-CFSFDP算法应用于用户的停留区域识别，结果表明，ST-CFSFDP算法在时间阈值90 s、距离阈值5 m的识别正确率高达82.4%，较经典的ST-DBCSAN、ST-OPTICS和ST-AGNES算法分别提高了5.2%、4.2%和7.6%。

本文尚在以下方面存在不足，需在后续工作中进一步研究：受试验数据采样间隔的限制，现有算法仅采用秒级时间粒度的数据进行验证,当定位数据的采样间隔增大时，算法识别准确率及可靠性需要进一步探究。

致谢：感谢上海图聚智能科技股份有限公司为本研究提供室内定位轨迹试验数据!