一种基于DTW和距离衰减的多维时间序列聚类方法与流程

一种基于dtw和距离衰减的多维时间序列聚类方法
技术领域
1.本发明涉及城市发展研究技术领域，特别涉及一种基于dtw和距离衰减的多维时间序列聚类方法。


背景技术：

2.城市建设用地是城市的发展重点、发展特点、发展任务的空间载体，城市间不同类别建设用地其规模、占比及变化趋势的相似程度能够在相当程度上体现城市发展的相似性。城市建设用地的变化研究对于空间规划改革下的城市发展有着重要意义。在国土空间规划的新语境下，探寻如何高效利用城市开发边界内的土地是此一轮城市发展的重要任务。而在城市发展过程中，如果忽略城市建设用地变化特征，盲目模仿发达城市投放土地进行建设，将带来浪费土地及资金、甚至城市建设“开倒车”的局面。因此，将城市间建设用地变化相似性进行分类，对于城市发展模式借鉴、城市群内部类型结构划分、试点城市批次选择、以及用地变化轨迹分析等研究是一项基础而重要的工作。
3.目前，鲜有研究涉及描述我国城市间建设用地变化的相似性，即：一些城市在时间维度上其各类建设用地此消彼长所呈现的变化趋势的共性。关于城市建设用地变化进行的研究常聚焦于以下几个方面：(1)将城市建设用地包含在更广义的“土地利用和土地覆盖变化(land-use and land-cover change，lulcc)”中，针对其变化特征，研究对于气候、环境等的影响；(2)针对个别城市建设用地，对其规模、时空形态、景观学特征等方面的变化特征进行研究；(3)对以上两种研究对象使用统计学模型、元胞自动机(cellular automation)或机器学习算法进行模拟及预测；(4)基于某背景、某视角，对建设用地的控制、引导、供给、开发使用、实施保障等方面探讨具体的方法及措施。
4.将相似的序列分组到同一个集群中是时间序列聚类(time-series clustering)问题一个重要目标。目前，对于时间序列聚类的计算，通常认为对于真实数据，dtw(dynamic time warping，即动态时间规整)。效果较优(tavenard，2020)。dtw方法基于动态规划(dp)的思想，以两条一维时间序列的长度构造矩阵网格，通过基于时间粒度的样本距离计算，动态寻找最小化的样本之间距离总和，即两条时间序列的相似程度。对于多维时间序列，dtw可以使用两种不同方式进行处理：独立与非独立翘曲(dependent orindependentwarping)。前者将研究中最底层数据所属的各种范畴之间看做相互独立，将属于各自范畴的单条时间序列进行两两比较，然后再在上层范畴中将其加，作为两条多维时间序列的相似度结果，简称为dtwi；后者将研究中最底层数据所属的各个范畴间看做具有紧密联系、需要同时参与计算的整体进行dtw处理，简称为dtwd。虽然绝大多数研究默认这两种方式对结果并无影响，但roemer(2016)发现，两种处理方式孰优孰劣并非以数据集的不同而异，其效果分异于所研究问题的类别，甚至依不同个案而不同。在现实世界中，需要进行相似度计算的多维时间序列，其维度之间的大小与比例变化通常具有内在的逻辑。因此针对dtwd的研究与改进逐渐受到重视。
5.通常在应用dtw时，对于对应时间切片上的数据，很多研究者直接使用欧氏距离
(euclidean distance)表示其相似度却并未给出理由。这种处理方式是研究实验逻辑的重要缺失。
6.此外，为解决两时间序列形状相似但所处时间位置有差异会被认为有较大不同的问题，经宇毅(2019)提出了有效的改进算法，即通过计算两序列的最长公共子串引入距离衰减系数。但这一算法只能应对一维时间序列的相似度计算，若要应用至多维时间序列，使其能够处理属于多维时间序列的城市用地数据、并使城市间发展路径相同、但所处发展阶段不同这种现实情况得到适当的处理，需要进一步改进。
7.在分类过程中，不同的聚类算法参数、对应时间切片上数值距离的算法，以及不同的衰减系数算法等这些要素经过不同的组合，会得出不同的聚类结果及效果。而现有的对于时间序列进行聚类的工作流程往往会忽视这一点，常采用研究人员自身所熟悉的工作流或算法组合，从而得出相对偏颇的结果，影响分类的准确性。


技术实现要素：

8.本发明要解决的技术问题，在于提供一种基于dtw和距离衰减的多维时间序列聚类方法，改进了dtwd，将距离衰减系数方法改进并应用到多维时间序列的相似性计算中，并结合机器学习聚类算法对dtwd得出的两两城市间用地变化相似度矩阵进行聚类，最终得到轮廓系数较高的分类结果。
9.本发明是这样实现的：
10.一种基于dtw和距离衰减的多维时间序列聚类方法，包括：
11.步骤10、将城市用地数据进行预处理后得到每个城市各类用地面积数值在时间维度上的序列，即多维时间序列；每个维度代表城市用地的一种特定类型，每个时间切片代表城市用地在某一时间点上的各类用地面积；
12.步骤20、通过dtwd计算任意两城市多维时间序列的距离，并求出公共子串判定阈值，而后计算公共子串长度；根据公共子串长度求得两多维时间序列的距离衰减系数，并与dtwd计算得到的两多维时间序列的距离相乘，得到经衰减后的两多维时间序列相似度，即两城市间的用地变化相似度，完成所有序列的两两比较得到城市间用地变化相似度矩阵；
13.步骤30、通过层次聚类算法对多城市间用地相似度矩阵进行聚类，然后采用轮廓系数来评价聚类效果，选取最佳轮廓系数的分类方案作为城市间用地变化相似性分类结果。
14.进一步地，所述步骤20具体包括：先将每个多维时间序列的各自维度求标准差再取平均值，求得每个多维时间序列的平均波动，即城市各类用地的平均变动幅度；在进行dtwd计算、使用欧氏距离计算对应时间切片上的数值距离时，用参与计算的两多维时间序列中平均标准差较大的值作为两序列公共子串的判定阈值；若对应时间切片上的数值距离小于这一阈值，代表在本次单时间切片比较时两城市的用地构成相似，则判定为公共子串长度加1；其后，由公共子串长度计算出距离衰减系数，并与dtwd计算得出的两多维时间序列的距离相乘，得到经衰减后的两多维时间序列相似度，即两城市间的用地变化相似度，以达成衰减系数在多维时间序列dtwd计算中的使用，完成所有序列的两两比较得到城市间用地变化相似度矩阵。
15.进一步地，具体包括如下步骤：
16.设n为多维时间序列的维数，即城市用地种类数，设m为多维时间序列的时间切片数量，即城市用地数据的时间长度；q和c为分别为参与两两对比的城市；
17.定义其中，为城市q第i类用地面积的均值，i＝1,2，
…
，n，j＝1,2，
…
，m；
18.stdqi定义成其中stdqi为城市q第i类用地的标准差；
19.然后即城市q用地面积的平均波动；
20.同样的定义：
21.其中，为城市c第i类用地的均值；
22.stdci定义成其中stdci为城市c第i类用地的标准差；
23.然后即城市c用地面积的平均波动；
24.最后，将公共子串判定阈值ts定义成max{μ(q)，μ(c)}；
25.令路π＝[π0,
…
，π
p
]，其中每个πk＝(ak，bk)，那么令π
min
满足
[0026]
其中代表城市q在年份ak上各类用地的数值，代表城市c在年份bk上各类用地的数值，d1即经dtwd求出的两城市用地数据的最小距离，d
euc
是欧氏距离；
[0027]
记公共子串长度其中#是指集合中元素的个数，以及
[0028][0029]
其中len(q)＝len(c)＝m，α1为参加本次dtw计算的两城市用地数据的距离衰减系数；
[0030]
最后定义
[0031]
sim1＝d1·
α1[0032]
即经衰减后的两城市间用地变化相似度。
[0033]
进一步地，所述步骤20具体包括：首先求两多维时间序列每个对应维度上较大的标准差，即两城市相同种类用地的平均波动的较大值作为这一维度的独立阈值；当两城市相同时间切片上同种类用地的数值的欧式距离都小于其所属维度的阈值时，就判定为公共子串长度加1；然后，再在dtwd计算时使用欧式距离、余弦距离或平方欧式距离单独计算两城市对应时间切片上的数值距离；其后，由公共子串长度计算出距离衰减系数，并与dtwd计算得出的两多维时间序列的距离相乘，得到经衰减后的两多维时间序列相似度，即两城市
间的用地变化相似度，以达成衰减系数在多维时间序列dtwd计算中的使用；完成所有序列的两两比较得到城市间用地变化相似度矩阵。
[0034]
进一步地，包括如下步骤：
[0035]
设n为多维时间序列的维数，即城市用地种类数，设m为多维时间序列的时间切片数量，即城市用地数据的时间长度；q和c为分别为两两对比的城市；
[0036]
定义tsi＝max{stdqi，stdci}，stdqi为城市q第i类用地的标准差，stdci为城市c第i类用地的标准差，在计算公共子串长度时，考虑集合
[0037][0038]
于是定义公共子串长度l2＝#v
[0039]
同样令路π＝[π0，
…
，π
p
]，其中每个πk＝(ak，bk)，那么令π
min
满足
[0040][0041]
其中x取euc,cos，sqeuc中的任一种；
[0042]
计算
[0043][0044]
其中len(q)＝len(c)＝m
[0045]
最后定义
[0046]
sim2＝d2·
α2[0047]
即经衰减后的两城市间用地变化相似度。
[0048]
进一步地，所述步骤20具体包括：将每个多维时间序列的各自维度求标准差再取平均值，求得每个多维时间序列的平均波动，即城市各类用地的平均变动幅度；使用dtwd计算对应时间切片上数值的欧氏距离，若求得的距离小于两序列中较大的平均标准差，则公共子串长度加1，再由公共子串长度计算出距离衰减系数；其后，再次使用dtwd、并在其中使用欧式距离、余弦距离或平方欧式距离计算两对应时间切片上的用地面积数值距离，得出两多维时间序列的距离；最后，将距离衰减系数与两多维时间序列的距离相乘，得到经衰减后的任意两城市间的用地相似度，以达成衰减系数在多维时间序列dtwd计算中的使用；完成所有序列的两两比较得到城市间用地变化相似度矩阵。
[0049]
进一步地，包括如下步骤：
[0050]
设n为多维时间序列的维数，即城市用地种类数，设m为多维时间序列的时间切片数量，即城市用地数据的时间长度；q和c为分别为两两对比的城市；
[0051]
定义其中，为城市q第i类用地面积的均值，i＝1,2，
…
，n，j＝1,2，
…
，m；
[0052]
stdqi定义成其中stdqi为城市q第i类用地的标准差；
[0053]
然后即城市q用地面积的平均波动；
[0054]
同样的定义
[0055]
其中，为城市c第i类用地的均值；
[0056]
stdci定义成其中stdci为城市c第i类用地的标准差；
[0057]
然后即城市c用地面积的平均波动；
[0058]
最后将阈值ts定义成max{μ(q)，μ(c)}；
[0059]
接着定义q
·
，j
＝(q
1，j
，q
2,j
，
…
，q
n，j
)，以及c
·
，j
＝(c
1，j
，c
2,j
，
…
，c
n,j
)
[0060]
令路π＝[π0,
…
，π
p
]，其中每个πk＝(ak，bk)，
[0061]
并记公共子串长度其中#是指集合中元素的个数，d
euc
是欧氏距离；
[0062]
同时，令π
min
满足
[0063][0064]
其中x取euc，cos，sqeuc
[0065]
同样计算
[0066][0067]
其中len(q)＝len(c)＝m
[0068]
最后定义
[0069]
sim3＝d3·
α3[0070]
即经衰减后的两城市间用地变化相似度。
[0071]
进一步地，所述步骤30中，所述层次聚类算法通过枚举最小化成对聚类间最远样本距离complete、最小化成对聚类间平均样本距离值average与最小化成对聚类间最近样本距离值single三个参数，得出多种分类可能并进行对比；轮廓系数s的计算公式为：
[0072][0073]
其中，a为某个样本与其所在簇内其他样本的平均距离，b为某个样本与下一个最邻近簇样本的平均距离。
[0074]
本发明实施例的技术方案至少具有如下优点：
[0075]
1、通过将距离衰减系数方法改进并应用到多维时间序列的相似性计算中，解决了多维时间序列因形状相似但位置不同导致的相似度降低问题，使所处发展阶段不同的同类型城市间用地变化相似性得以体现；结合机器学习聚类算法探索城市间建设用地变化的类别划分，最终得到轮廓系数较高的聚类结果，即多维时间序列的聚类效果得到了整体性的
提高；
[0076]
2、通过使用层次聚类模块对已计算出的城市间建设用地变化相似度矩阵进行聚类得出最终分类结果，减少了计算复杂度，并采用多种簇间的连接方式得到更多的分类可能并进行对比，然后采用轮廓系数来评价聚类效果，不需要已知分类样例，直接以簇内差异尽可能小、簇间差异尽可能大为标准衡量聚类效果，得到整体最好的聚类效果。
附图说明
[0077]
下面参照附图结合实施例对本发明作进一步的说明。
[0078]
图1为本发明实施例方法的整体流程示意图；
[0079]
图2为本发明实施例进行预处理和清洗后的指定城市用地数据示意图。
具体实施方式
[0080]
本发明实施例提供一种基于dtw和距离衰减的多维时间序列聚类方法，改进了dtwd，将距离衰减系数方法改进并应用到多维时间序列的相似性计算中，并结合机器学习聚类算法对dtwd得出的两两城市间用地变化相似度矩阵进行聚类，最终得到轮廓系数较高的分类结果。
[0081]
本发明实施例的总体构思如下：
[0082]
对多个城市用地相似性的分类问题，可进一步理解为：每个城市所包含的多个种类建设用地面积数值的时间序列，其变化特征与其他城市的相似程度问题。
[0083]
在应用dtwd计算任意两城市间用地变化的相似性、即两条多维时间序列的相似度时，求两序列在某一时间粒度上的切片的距离这一重要步骤并不限于其他研究者常用的欧氏距离一种。根据wolfram(wolfram research,inc.:www.wolframalpha.com)的指南手册，基于数值数据的距离或相似度测度有12种：去除含有特殊含义的算法(例如街区距离manhattan distance、棋盘距离chessboard distance、二进制距离binary distance、物种组成差异braycurtis distance、加权曼哈顿距离canberra distance、相关系数距离correlation distance、时间翘曲warping distance与canonical warping distance等)，以及去除在两序列对应时间切片上若有相同值会导致计算结果出现nan的标准化欧氏距离(normalized squared euclidean distance)后，本技术采用剩下的3种算法计算两单个时间切片之间的距离，即：欧氏距离(euclidean distance)、平方欧氏距离(squared euclidean distance)、余弦距离。虽然调整余弦距离(adjusted cosine distance)比余弦距离更能反映尺度上的不同(badrul m.sarwar，2001)，但考虑到现实中有可能出现城市类型相同但规模不同的情况，即时间序列的形状相似但数值大小有整体性差异，因此并未采用。
[0084]
在进行多维时间序列的dtwd计算时，需要考虑城市所处发展阶段不同但类型相同这一常见现象。例如某新兴矿业城市和某矿产资源枯竭城市，前者的用地构成变化所显现的特征，可能与后者在十几年前所呈现出的特征相似。亦或a城市大力借鉴沿海b发达城市的建设经验，日后必反映在用地构成上的趋同，但在时间上有一定滞后。若不将这一常见现象进行适当处理，即将形状相近但位置不同的两序列判定为一定程度上的相似，则会得出和现实情况相反的计算结果。
[0085]
因此，需要对现有的仅在一维时间序列的dtw计算中可使用的距离衰减系数算法进行改进，使其计算理念能够应用于城市用地这种多维时间序列的相似性计算中。在两条一维时间序列的dtw计算时，已有的衰减系数法使用两时间序列中较大的那个标准差作为判定两序列公共子串的阈值。而在多维时间序列情况下，每个维度上的数据代表不同种类的用地面积，其各自代表的范畴不同，不能使用上述方法来计算公共子串的判定阈值，故本技术实施例提出了三种阈值计算及公共子串判定方式，实现对城市建设用地这类多维时间序列的距离衰减。
[0086]
其后，对多城市间的用地变化相似度矩阵进行聚类，可使用基于相似度矩阵的聚类算法——层次聚类得出聚类结果。同时，将不同参数枚举，寻找最佳分类效果方案。
[0087]
本实施例提供了一种基于dtw和距离衰减的多维时间序列聚类方法，如图1和图2所示，包括：
[0088]
步骤10、将城市用地数据进行预处理后得到每个城市各类用地面积数值在时间维度上的序列，即多维时间序列；每个维度代表城市用地的一种特定类型，每个时间切片代表城市用地在某一时间点上的各类用地面积；
[0089]
步骤20、通过dtwd计算任意两城市多维时间序列的距离，并求出公共子串判定阈值，而后计算公共子串长度；根据公共子串长度求得两多维时间序列的距离衰减系数，并与dtwd计算得到的两多维时间序列的距离相乘，得到经衰减后的两多维时间序列相似度，即两城市间的用地变化相似度，完成所有序列的两两比较得到城市间用地变化相似度矩阵；
[0090]
步骤30、通过层次聚类算法对多城市间用地相似度矩阵进行聚类，然后采用轮廓系数来评价聚类效果，选取最佳轮廓系数的分类方案作为城市间用地变化相似性分类结果。
[0091]
通过将距离衰减系数方法改进并应用到多维时间序列的相似性计算中，解决了多维时间序列因形状相似但位置不同导致的相似度降低问题，使所处发展阶段不同的同类型城市间用地变化相似性得以体现；结合机器学习聚类算法探索城市间建设用地变化的类别划分，最终得到轮廓系数较高的聚类结果，即多维时间序列的聚类效果得到了整体性的提高。
[0092]
步骤20包括了本发明提出的改进dtwd，通过将距离衰减系数计算中的公共子串判定阈值方法改进并应用到多维时间序列的相似性计算的方法，可以称为attenuating dynamic time warping with dependent warping(a-dtwd)，其有三种具体的实现方法，分别如下：
[0093]
第一种实现方法：先将每个多维时间序列的各自维度求标准差再取平均值，求得每个多维时间序列的平均波动，即城市各类用地的平均变动幅度；在进行dtwd计算、使用欧氏距离计算对应时间切片上的数值距离时，用参与计算的两多维时间序列中平均标准差较大的值作为两序列公共子串的判定阈值；若对应时间切片上的数值距离小于这一阈值，代表在本次单时间切片比较时两城市的用地面积构成相似，则判定为公共子串长度加1；其后，由公共子串长度计算出距离衰减系数，并与dtwd计算得出的两多维时间序列的距离(即两城市间用地数据距离)相乘，得到经衰减后的两多维时间序列相似度，即两城市间的用地变化相似度，以达成衰减系数在多维时间序列dtwd计算中的使用。完成所有序列的两两比较得到城市间用地变化相似度矩阵。这一方法只能在使用欧氏距离计算单时间切片上数值
距离时使用，具体可以包括如下步骤：
[0094]
设n为多维时间序列的维数，即城市用地种类数，设m为多维时间序列的时间切片数量，即城市用地数据的时间长度；q和c为分别为参与两两对比的城市；
[0095]
定义其中，为城市q第i类用地面积的均值，i＝1,2,
…
,n，j＝1,2,
…
,m；
[0096]
stdqi定义成其中stdqi为城市q第i类用地的标准差；
[0097]
然后即城市q用地面积的平均波动；
[0098]
同样的定义：
[0099]
其中，为城市c第i类用地的均值；
[0100]
stdci定义成其中stdci为城市c第i类用地的标准差；
[0101]
然后即城市c用地面积的平均波动；
[0102]
最后，将公共子串判定阈值ts定义成max{μ(q)，μ(c)}；
[0103]
令路π＝[π0,
…
，π
p
]，其中每个πk＝(ak，bk)，那么令π
min
满足
[0104]
其中代表城市q在年份ak上各类用地的数值，代表城市c在年份bk上各类用地的数值，d1即经dtwd求出的两城市用地数据的最小距离，d
euc
是欧氏距离；
[0105]
记公共子串长度其中#是指集合中元素的个数，以及
[0106][0107]
其中len(q)＝len(c)＝m，α1为参加本次dtw计算的两城市用地数据的距离衰减系数；
[0108]
最后定义
[0109]
sim1＝d1·
α1[0110]
即经衰减后的两城市间用地变化相似度。
[0111]
第二种实现方法：首先求两多维时间序列每个对应维度上较大的标准差，即两城市相同种类用地的平均波动的较大值作为这一维度的独立阈值。当两城市相同时间切片上同种类用地的数值的欧式距离都小于其所属维度的阈值时，就判定为公共子串长度加1；然后，再在dtwd计算时使用欧式距离、余弦距离或平方欧式距离单独计算两城市对应时间切片上的数值距离。其后，由公共子串长度计算出距离衰减系数，并与dtwd计算得出的两多维
时间序列的距离(即两城市间用地数据距离)相乘，得到经衰减后的两多维时间序列相似度，即两城市间的用地变化相似度，以达成衰减系数在多维时间序列dtwd计算中的使用。完成所有序列的两两比较得到城市间用地变化相似度矩阵。即从之前的“先算数值距离，再用阈值判定”变为“先用阈值判定，再算数值距离”，使阈值的计算不再受限于单时间切片上数值距离使用的算法。
[0112]
具体可以包括如下步骤：
[0113]
设n为多维时间序列的维数，即城市用地种类数，设m为多维时间序列的时间切片数量，即城市用地数据的时间长度；q和c为分别为两两对比的城市；
[0114]
定义tsi＝max{stdqi，stdci}，stdqi为城市q第i类用地的标准差，stdci为城市c第i类用地的标准差，在计算公共子串长度时，考虑集合
[0115][0116]
于是定义公共子串长度l2＝#v
[0117]
同样令路π＝[π0,
…
，π
p
]，其中每个πk＝(ak，bk)，那么令π
min
满足
[0118][0119]
其中x取euc,cos，sqeuc中的任一种；
[0120]
计算
[0121][0122]
其中len(q)＝len(c)＝m
[0123]
最后定义
[0124]
sim2＝d2·
α2[0125]
即经衰减后的两城市间用地变化相似度。
[0126]
第三种实现方法：综合了前两种方法的思想，亦是“先用阈值判定，再算数值距离”。只是这里的阈值不再是各自维度拥有的独立阈值，而是使用第一种改进方式中的两序列较大平均标准差作为阈值。这一阈值衡量多维时间序列的整体平均波动，相较于前者对各个维度进行独立限制，相对宽松。
[0127]
将每个多维时间序列的各自维度求标准差再取平均值，求得每个多维时间序列的平均波动，即城市各类用地的平均变动幅度；使用dtwd计算对应时间切片上数值的欧氏距离，若求得的距离小于两序列中较大的平均标准差，则公共子串长度加1，再由公共子串长度计算出距离衰减系数；其后，再次使用dtwd、并在其中使用欧式距离、余弦距离或平方欧式距离计算两对应时间切片上的用地面积数值距离，得出两多维时间序列的距离；最后，将距离衰减系数与两多维时间序列的距离相乘，得到经衰减后的任意两城市间的用地相似度，以达成衰减系数在多维时间序列dtwd计算中的使用。完成所有序列的两两比较得到城市间用地变化相似度矩阵。具体可以包括如下步骤：
[0128]
设n为多维时间序列的维数，即城市用地种类数，设m为多维时间序列的时间切片数量，即城市用地数据的时间长度；q和c为分别为两两对比的城市；
[0129]
定义其中，为城市q第i类用地面积的均值，i＝1,2，
…
，n，j＝1,2，
…
，m；
[0130]
stdqi定义成算中stdqi为城市q第i类用地的标准差；
[0131]
然后即城市q用地面积的平均波动；
[0132]
同样的定义
[0133]
算中，为城市c第i类用地的均值；
[0134]
stdci定义成其中stdci为城市c第i类用地的标准差；
[0135]
然后即城市c用地面积的平均波动；
[0136]
最后将阈值ts定义成max{μ(q)，μ(c)}；
[0137]
接着定义q
·
，j
＝(q
1，j
，q
2,j
，
…
，q
n，j
)，以及c
·
，j
＝(c
1，j
，c
2,j
，
…
，c
n,j
)
[0138]
令路π＝[π0,
…
，π
p
]，其中每个πk＝(ak，bk)，
[0139]
并记公共子串长度其中#是指集合中元素的个数，d
euc
是欧氏距离；
[0140]
同时，令π
min
满足
[0141][0142]
其中x取euc，cos，sqeuc
[0143]
同样计算
[0144][0145]
其中len(q)＝len(c)＝m
[0146]
最后定义
[0147]
sim3＝d3·
α3[0148]
即经衰减后的两城市间用地变化相似度。
[0149]
就这样将所有城市两两成对计算，可以得出多城市间的用地变化相似度矩阵，然后使用合适的聚类算法将其形成最终的分类结果。
[0150]
传统机器学习算法分为无监督与有监督两大类，区别在于是否先使用具有已知标签的样本去训练得到一个最优模型，再利用此模型将输入的未知标签的数据进行判断，从而实现分类。对于本技术方案而言，并没有可以事先参考的用地模式划分实例，即没有训练数据，因此只能使用无监督聚类算法。
[0151]
在目前流行的scikit-learn机器学习库的无监督聚类算法中，能够直接使用两两
相似矩阵作为输入的算法有谱聚类(spectral clustering)、ap聚类(affinity propagation)与层次聚类(hierarchical clustering)。与前两种聚类算法不同，层次聚类通过不断的合并或者分割前置聚类来构建最终聚类，其结果呈树状。基于这一特点，层次聚类的结果相较其他部分聚类(partial clustering)输出的非结构化扁平簇(flat cluster)结果能够提供更多的信息。这一树状结构为不同类别数量的聚类结果之间提供了连续性，意味着在划分城市用地模式类别时，不同类别数量下的结果不但不会自相矛盾，反而会由类别数量从少至多呈现出“更进一步的细分”的连贯性，为计算结果的解释与应用提供了便利。为减少计算复杂度，本实施例选择使用scikit-learn中的层次聚类模块中的agglomerativeclustering功能对已计算出的相似度矩阵进行计算得出最终分类结果。
[0152]
在使用agglomerativeclustering功能时，本实施例对簇间的连接方式(linkage，意即两个类别之间的距离)这一参数通过枚举complete(最小化成对聚类间最远样本距离)、average(最小化成对聚类间平均样本距离值)与single(最小化成对聚类间最近样本距离值)得出更多的分类方案并进行对比。同时，采用轮廓系数(silhouette score)(rousseeuw，1987)来评价聚类效果，其优点是不需要已知分类样例，直接以簇内差异尽可能小、簇间差异尽可能大为标准衡量聚类效果，公式如下：
[0153][0154]
其中，a为某个城市与其所在类别内其他城市的平均距离(相似度)，b为某个城市与下一个最相似类别中各城市的平均距离，s为各个城市的轮廓系数。本实施例取用所有城市的平均轮廓系数值作为以不同类别数量为目标的聚类效果的评价，范围为[-1,1]，越接近1则说明聚类效果越好，反之越差。计算从分类数量为2开始进行，直至每个城市自成一类结束，最终将不同类别数量下聚类的平均轮廓系数进行加和，作为不同公共子串判定阈值方法及类别间距离算法参数下的整体聚类效果评分。计算结果如表1：
[0155]
表1不同方法及参数下聚类轮廓系数对比
[0156]
[0157][0158]
由计算出的结果可知，轮廓系数随类别数量的增加而整体呈下降趋势，中间有或大或小的波动。随着对应时间切片上数值距离算法、簇间连接方式和公共子串判定阈值方法的变化，聚类效果是不尽相同的，有时会提高聚类效果，有时则会拉低总轮廓系数。
[0159]
通过使用层次聚类模块对已计算出的城市间建设用地变化相似度矩阵进行聚类得出最终分类结果，减少了计算复杂度，并采用多种簇间的连接方式得到更多的分类可能并进行对比，然后采用轮廓系数来评价聚类效果，不需要已知分类样例，直接以簇内差异尽可能小、簇间差异尽可能大为标准衡量聚类效果，得到整体最好的聚类效果。
[0160]
无论是采用何种两个类别之间的距离算法参数、以及两城市对应时间切片上两组数值距离算法，经过本发明提出的公共子串判定阈值方法进行距离衰减后，均得到了相较不经过距离衰减的聚类更高的轮廓系数分值。可以得出，经过本发明所提出的距离衰减方法处理后，多维时间序列的分类效果得到了整体性的提高。
[0161]
同时，相较传统研究中默认使用欧氏距离的做法，在本发明提出的方法及参数枚举方法下，两个类别之间的距离算法使用average参数、对应时间切片之间的距离计算采用平方欧氏距离、并采用第三种公共子串判定阈值方法进行距离衰减，得到了本次实验的最高分，即最佳的城市间建设用地变化相似性分类结果。但需要指出的是，这一结果并不代表此组算法及参数组合就是普适最优的。根据数据的所属范畴及应用场景不同，在计算中枚举不同的对应时间切片之间的距离算法、公共子串判定阈值方法和聚类参数是必要的。
[0162]
虽然以上描述了本发明的具体实施方式，但是熟悉本技术领域的技术人员应当理解，我们所描述的具体的实施例只是说明性的，而不是用于对本发明的范围的限定，熟悉本领域的技术人员在依照本发明的精神所作的等效的修饰以及变化，都应当涵盖在本发明的权利要求所保护的范围内。

技术特征：
1.一种基于dtw和距离衰减的多维时间序列聚类方法，其特征在于，包括：步骤10、将城市用地数据进行预处理后得到每个城市各类用地面积数值在时间维度上的序列，即多维时间序列；每个维度代表城市用地的一种特定类型，每个时间切片代表城市用地在某一时间点上的各类用地面积；步骤20、通过dtw
d
计算任意两城市多维时间序列的距离，并求出公共子串判定阈值，而后计算公共子串长度；根据公共子串长度求得两多维时间序列的距离衰减系数，并与dtw
d
计算得到的两多维时间序列的距离相乘，得到经衰减后的两多维时间序列相似度，即两城市间的用地变化相似度，完成所有序列的两两比较得到城市间用地变化相似度矩阵；步骤30、通过层次聚类算法对多城市间用地相似度矩阵进行聚类，然后采用轮廓系数来评价聚类效果，选取最佳轮廓系数的分类方案作为城市间用地变化相似性分类结果。2.根据权利要求1所述的方法，其特征在于：所述步骤20具体包括：先将每个多维时间序列的各自维度求标准差再取平均值，求得每个多维时间序列的平均波动，即城市各类用地的平均变动幅度；在进行dtw
d
计算、使用欧氏距离计算对应时间切片上的数值距离时，用参与计算的两多维时间序列中平均标准差较大的值作为两序列公共子串的判定阈值；若对应时间切片上的数值距离小于这一阈值，代表在本次单时间切片比较时两城市的用地构成相似，则判定为公共子串长度加1；其后，由公共子串长度计算出距离衰减系数，并与dtw
d
计算得出的两多维时间序列的距离相乘，得到经衰减后的两多维时间序列相似度，即两城市间的用地变化相似度，以达成衰减系数在多维时间序列dtw
d
计算中的使用，完成所有序列的两两比较得到城市间用地变化相似度矩阵。3.根据权利要求2所述的方法，其特征在于，具体包括如下步骤：设n为多维时间序列的维数，即城市用地种类数，设m为多维时间序列的时间切片数量，即城市用地数据的时间长度；q和c为分别为参与两两对比的城市；定义其中，为城市q第i类用地面积的均值，i＝1,2,
…
,n，j＝1,2,
…
,m；stdq
i
定义成其中stdq
i
为城市q第i类用地的标准差；然后即城市q用地面积的平均波动；同样的定义：,其中，为城市c第i类用地的均值；stdc
i
定义成其中stdc
i
为城市c第i类用地的标准差；然后即城市c用地面积的平均波动；最后，将公共子串判定阈值ts定义成max{μ(q),μ(c)}；令路π＝[π0,
…
,π
p
],其中每个π
k
＝(a
k
,b
k
),那么令π
min
满足
其中代表城市q在年份a
k
上各类用地的数值，代表城市c在年份b
k
上各类用地的数值，d1即经dtw
d
求出的两城市用地数据的最小距离,d
euc
是欧氏距离；记公共子串长度其中#是指集合中元素的个数，以及其中len(q)＝len(c)＝m，α1为参加本次dtw计算的两城市用地数据的距离衰减系数；最后定义sim1＝d1·
α1即经衰减后的两城市间用地变化相似度。4.根据权利要求1所述的方法，其特征在于：所述步骤20具体包括：首先求两多维时间序列每个对应维度上较大的标准差，即两城市相同种类用地的平均波动的较大值作为这一维度的独立阈值；当两城市相同时间切片上同种类用地的数值的欧式距离都小于其所属维度的阈值时，就判定为公共子串长度加1；然后，再在dtw
d
计算时使用欧式距离、余弦距离或平方欧式距离单独计算两城市对应时间切片上的数值距离；其后，由公共子串长度计算出距离衰减系数，并与dtw
d
计算得出的两多维时间序列的距离相乘，得到经衰减后的两多维时间序列相似度，即两城市间的用地变化相似度，以达成衰减系数在多维时间序列dtw
d
计算中的使用；完成所有序列的两两比较得到城市间用地变化相似度矩阵。5.根据权利要求4所述的方法，其特征在于，包括如下步骤：设n为多维时间序列的维数，即城市用地种类数，设m为多维时间序列的时间切片数量，即城市用地数据的时间长度；q和c为分别为两两对比的城市；定义ts
i
＝max{stdq
i
,stdc
i
}，stdq
i
为城市q第i类用地的标准差，stdc
i
为城市c第i类用地的标准差，在计算公共子串长度时,考虑集合于是定义公共子串长度l2＝#v同样令路π＝[π0,
…
,π
p
],其中每个π
k
＝(a
k
,b
k
),那么令π
min
满足其中x取euc,cos,sqeuc中的任一种；计算其中len(q)＝len(c)＝m
最后定义sim2＝d2·
α2即经衰减后的两城市间用地变化相似度。6.根据权利要求1所述的方法，其特征在于，所述步骤20具体包括：将每个多维时间序列的各自维度求标准差再取平均值，求得每个多维时间序列的平均波动，即城市各类用地的平均变动幅度；使用dtw
d
计算对应时间切片上数值的欧氏距离，若求得的距离小于两序列中较大的平均标准差，则公共子串长度加1，再由公共子串长度计算出距离衰减系数；其后，再次使用dtw
d
、并在其中使用欧式距离、余弦距离或平方欧式距离计算两对应时间切片上的用地面积数值距离，得出两多维时间序列的距离；最后，将距离衰减系数与两多维时间序列的距离相乘，得到经衰减后的任意两城市间的用地相似度，以达成衰减系数在多维时间序列dtw
d
计算中的使用；完成所有序列的两两比较得到城市间用地变化相似度矩阵。7.根据权利要求6所述的方法，其特征在于，包括如下步骤：设n为多维时间序列的维数，即城市用地种类数，设m为多维时间序列的时间切片数量，即城市用地数据的时间长度；q和c为分别为两两对比的城市；定义其中，为城市q第i类用地面积的均值，i＝1,2,
…
,n，j＝1,2,
…
,m；stdq
i
定义成其中stdq
i
为城市q第i类用地的标准差；然后即城市q用地面积的平均波动；同样的定义其中，为城市c第i类用地的均值；stdc
i
定义成其中stdc
i
为城市c第i类用地的标准差；然后即城市c用地面积的平均波动；最后将阈值ts定义成max{μ(q),μ(c)}；接着定义q
·
,j
＝(q
1,j
,q
2,j
,
…
,q
n,j
),以及c
·
,j
＝(c
1,j
,c
2,j
,
…
,c
n,j
)令路π＝[π0,
…
,π
p
],其中每个π
k
＝(a
k
,b
k
),并记公共子串长度其中#是指集合中元素的个数,d
euc
是欧氏距离；同时，令π
min
满足其中x取euc,cos,sqeuc同样计算
其中len(q)＝len(c)＝m最后定义sim3＝d3·
α3即经衰减后的两城市间用地变化相似度。8.根据权利要求1至7任一项所述的方法，其特征在于，所述步骤30中，所述层次聚类算法通过枚举最小化成对聚类间最远样本距离complete、最小化成对聚类间平均样本距离值average与最小化成对聚类间最近样本距离值single三个参数，得出多种分类可能并进行对比；轮廓系数s的计算公式为：其中，a为某个样本与其所在簇内其他样本的平均距离，b为某个样本与下一个最邻近簇样本的平均距离。

技术总结
本发明公开了一种基于DTW和距离衰减的多维时间序列聚类方法，涉及城市发展研究技术领域。首先在使用DTW


技术研发人员：郭鑫 苏海龙 王德 宋强 纪立虎 郭峰
受保护的技术使用者：郭鑫
技术研发日：2023.04.20
技术公布日：2023/8/28