2. 西南交通大学 生命科学与工程学院, 四川 成都 610031;
3. 浙江大学 计算机科学与技术学院, 浙江 杭州 310027;
4. 科学计算与智能信息处理广西高校重点实验室(广西师范学院), 广西 南宁 530023;
5. 四川大学 计算机学院, 四川 成都 610065
2. School of Life Science and Engineering, Southwest Jiaotong University, Chengdu 610031, China;
3. College of Computer Science and Technology, Zhejiang University, Hangzhou 310027, China;
4. Science Computing and Intelligent Information Processing of Guangxi Higher Education Key Laboratory (Guangxi Teachers Education University), Nanning 530023, China;
5. College of Computer Science, Sichuan University, Chengdu 610065, China
随着我国智慧城市建设飞速发展,智能交通系统(intelligent transportation system,简称ITS)的应用变得日趋普及和大众化.文献[1]指出了大规模数据驱动下ITS系统中轨迹的重要性和可能形成的崭新应用.目前,各种交通信息采集技术(如RFID、无线传感器、视频监控等)已被广泛地运用于高速公路、城市交通的路段和卡口,每天都会收集海量的实时交通数据.同时,问题和困境突显:(1) 智能交通的潜在价值还没有得到有效挖掘;(2) 对交通信息的感知和收集能力有限;(3) 对存在于各个ITS中的海量数据无法高效存储及运用、有效分析,缺乏对交通态势的预测和研判能力,对公众的实时交通信息服务很难满足.如何对海量的移动对象轨迹数据进行准确、高效的分析处理及预测,做出即时和正确的交通疏通,进而为有效改善实际交通拥堵状况提供更加智能的基于位置的服务,成为目前亟待解决的关键问题.一个典型应用是:Wang等人[2]提出了一种基于大规模出租车轨迹数据的交通拥堵可视化分析方法,通过道路速度计算对交通拥堵自动检测.
基于位置的智能服务[3]作为一种新产生的技术,为解决轨迹数据挖掘问题提供了新思路,其目标在于挖掘对象的空间位置所蕴含的价值,分析对象之间的位置关系,从而提供个性化的位置服务.例如,通过记录用户以往的就餐地点,分析其饮食、消费习惯,为用户提供满意的餐厅.为了准确地提供上述基于位置的服务,实时甚至提前获取用户的位置变得尤为重要,因此,移动对象轨迹预测逐渐成为当前研究热点.关于移动对象轨迹预测的研究,为相关科研工作提供了理论基础.例如,对于GPS轨迹数据的预处理工作可以辅助提高具体应用的时间效率,提供了更加规范的数据,避免了繁琐的前期工作.此外,轨迹预测工作有利于分析人的行为特征,为社会学家分析个体及群体行为模式提供辅助支持.同时,移动对象轨迹预测具有较高的应用价值.例如,由于出租车数量有限,人们往往难以在高峰期打到车,利用轨迹预测算法可以有效地调度出租车.
位置大数据(location big data)[4]是构成位置社会感知的重要资源,具有相当大的体量,其所使用的数据集在规模和复杂程度上均已达到了“大”数据的层次,代表性实例见表 1.
![]() |
Table 1 Instances of location big data 表 1 位置大数据实例 |
当前,针对大数据环境下的移动对象轨迹预测研究的方法和模型较少,然而,仅一条完整的轨迹往往就包含成百上千的位置点,亟需高效的位置大数据预处理方法.此外,对轨迹进行查询往往需要为位置点构建索引,已有的索引结构,如DISC-tree[7],均需对路网和轨迹点构建双层树状索引,代价极高.
1 相关工作国内外研究者针对移动对象频繁模式挖掘及轨迹预测已经展开了相关研究:在使用时空数据挖掘技术挖掘频繁轨迹模式方面,Monreale[8]提出了WhereNext方法,从历史轨迹数据中抽取对象的运动模式,借助其表示用户经常频繁访问的地点,同时利用T-pattern tree查询发现最佳匹配轨迹;Ying等人[9]基于地理和轨迹的语义特征预测移动对象下一时刻位置点信息,该方法通过挖掘同类用户的常见行为特性来预测其未来位置.
马尔可夫链模型在挖掘频繁模式方面应用广泛.Song等人[10]提出基于状态的移动对象运动模型,运用马尔可夫转移概率解释移动对象在不同状态间的转换.Ishikawa等人[11]将整个交通路网划分为大小不同的网格并判断移动对象所处位置,以R-tree作为辅助索引,使用马尔可夫链描述移动对象在网格间转移的概率,但所提算法主要适用于针对特定区域的分析和精确查询.Asahara[12]借助混合马尔可夫链模型预测行人运动,综合考虑行人的个人特征和历史状态进行预测.实验结果表明,与马尔可夫链和隐马尔可夫模型相比,混合马尔可夫链模型具有较高的预测准确率.Gambs[13]扩展了Mobility Markov Chain(MMC)模型用于移动对象位置预测,本质是高阶马尔可夫模型,预测精度在70%~95%之间,但计算开销较大.Qiao等人[14]将HMM应用于移动对象轨迹预测,但未考虑大数据环境下算法的运行时间性能问题.此外,Qiao等人[15]利用高斯混合模型对移动对象复杂运动模式建模,统计不同运动模式的概率分布,进而实现准确、高效的位置预测.
基于HMM的轨迹查询预测算法通常将空间划分为不相交的区域,利用区域代表真实的轨迹点,精简了轨迹数据表示过程,提高了预测计算的速度,适用于位置大数据的分析预测.但是,基于HMM的预测模型不可避免地存在答案丢失[16]和精度依赖问题(预测精度依赖于历史轨迹数据).此外,上述算法基于向量空间进行预测,仅适用于不同路口间的路段,当移动对象运动至道路交汇处时,该方法无法预测具体的移动方向.再者,HMM模型中有很多参数需要人为设置,不同轨迹数据上预测精度差异较大.针对上述方法的不足,本文主要提出了两项创新性工作:(1) 新的位置大数据处理方法,利用基于密度的聚类方法减少HMM中状态数量,并对轨迹数据进行分段处理,以提取局部位置数据特征;(2) 新的参数自适应轨迹预测算法,通过计算机轨迹点最小间隔,自适应调整区域划分网格大小和轨迹分段的大小.
2 基本概念及问题描述针对移动对象位置大数据进行轨迹预测主要步骤包括:(1) 地图预处理,实现网格化分区;(2) 位置轨迹预处理,通过对轨迹数据的简化、抽象和聚集操作,建立高效轨迹预测模型;(3) 局部位置数据特征提取,提取出移动对象的高阶连续运动模式;(4) 位置大数据建模,对轨迹数据的时间和空间维度降维,构建全局模型;(5) 特征关联及预测,运用抽取的轨迹模式借助轨迹预测技术挖掘对象的运动趋势.
本文将轨迹预测问题具体描述为:对于某一特定用户,已知给定轨迹{T(t)|t=1,2,3,…,n,表示不同时间点}和该用户的预测模型参数λ={π,A,B|π表示初始化概率,A表示隐状态转移概率矩阵,B表示隐状态与观察状态间转移概率矩阵},求解n+1时刻移动对象位置T(n+1).通常,移动对象的轨迹数据由二维坐标点组成,将多个点按时间先后顺序连接起来形成一条完整的动态运动轨迹.轨迹位置大数据除了具有大数据的4V特性以外,还具有时序性、不确定性和动态变化性.以GPS数据为例,轨迹数据是离散分布在地图上的点,而且体量很大,如果对每条轨迹利用链表进行表达,那么将会占用大量的内存空间,而且不同轨迹之间的相似度将难以描述.因此,为了有效地表示轨迹,本文将平面区域划分成不相交的区域,即,网格.这样就可以使用网格序列代表点序列,完成对轨迹的表示.网格表示法简化了轨迹,优化了内存的使用,同时提高了轨迹相似度计算的效率,但是不可避免地带来了答案丢失问题(具体参见第4.1节介绍).
为解决上述问题,本文引入隐马尔可夫模型,首先,通过模型定义完成对轨迹模型的建立;然后,通过解答隐马尔可夫模型的3个问题,即,评估问题、解码问题、学习问题,达到轨迹预测的目的.下面给出基于隐马尔可夫模型进行轨迹预测的基本概念,由马尔可夫链引出隐马尔可夫模型.
定义1(轨迹序列). 给定欧氏空间,轨迹序列S={s1,s2,…,sn}为按时间排序的离散轨迹点,si表示第i个轨迹点,si=(xi,yi,ti),1≤i≤n.
定义2(预测轨迹序列). 已知轨迹S={s1,s2,…,sk},利用轨迹预测模型,得到预测轨迹序列Tp={Tp1,Tp2,…, Tpn},其中,n>k,Tpi表示第i个轨迹点.
定义3(马尔可夫轨迹链). 已知轨迹随机过程Ω={X(t),t∈T}的状态空间Σ是有限集或可列集,对于T内任意n+1个轨迹时间序列t1<t2<t3<…<tn<tn+1和Σ内任意n+1个状态j1,j2,j3,…,jn,jn+1,如果条件概率:
$\begin{align} & P\left\{ X\left( {{t}_{n}}_{+1} \right)={{j}_{n}}_{+1}|X\left( {{t}_{1}} \right)={{j}_{1}},X\left( {{t}_{2}} \right)={{j}_{2}},X\left( {{t}_{3}} \right)={{j}_{3}},\ldots ,X\left( {{t}_{n}} \right)={{j}_{n}} \right\} \\ & =P\left\{ X\left( {{t}_{n}}_{+1} \right)={{j}_{n}}_{+1}|X\left( {{t}_{n}} \right)={{j}_{n}} \right\} \\ \end{align}$ | (1) |
恒成立,则称此过程为马尔可夫轨迹链.公式(1)称为马尔可夫性质,或称无后效性.
定义4(轨迹状态概率). pi(t)=P{X(t)=aj,t∈T}称为马尔可夫轨迹链X(t)=aj的状态概率.
定义5(轨迹状态转移概率). pij(t,t')=P{X(t')=aj|X(t)=ai},称为马尔可夫轨迹链在X(t)=ai的条件下X(t')=aj的状态转移概率;由pij组成的矩阵称为马尔可夫轨迹链的状态转移矩阵,描述由状态ai转移到aj的概率.
定义6(隐马尔可夫轨迹模型)[14]. 用一个六元组H={S,H,R,π,A,B}来描述,假设X(n)表示轨迹隐状态序列, O(n)表示轨迹观测状态序列,则:
• S:表示轨迹序列,S={s1,s2,…,sn}为有序的离散轨迹点{si=(xi,yi,t)|1≤i≤n}所构成的序列;
• H:轨迹隐状态集合,通过对训练轨迹数据按固定长度分段形成.隐状态由{hi,i=1,2,3,...,N}表示,N表示隐藏状态数量;
• R:轨迹观测状态集合,由空间网格划分单元组成,由{oi,i=1,2,3,...,M}表示,M表示观测状态数量;
• П={πi}:轨迹初始状态概率分布,πi=P(q1=i)表示初始化时选择某个状态的概率;
• A={aij}:隐状态转移概率矩阵,其中,hij=P{X(t+1)=αj|X(t)=ai},即,下层Markov轨迹链状态转移矩阵;
• B={bik}:轨迹观测状态与隐状态转移概率矩阵,也称为混淆矩阵,其中,bik=P{O(t)=ok|X(t)=hi},描述某个时刻由隐藏状态ai得到观测状态ok的概率.
隐状态转移矩阵和混淆矩阵与时间无关,即,当系统演化时,A和B并不随时间改变.当N和M固定时,本文使用λ={π,A,B}表示模型的参数.
图 1是隐马尔可夫轨迹模型的一个实例,椭圆代表训练轨迹通过划分后形成的隐状态,二维移动空间平面网格划分单元作为观测状态.其中,轨迹隐状态集合H={s1,s2,s3,s4,s5},状态数N=5,轨迹观测状态集合R={g1,g2,g3, g4,g5,g6,g7,g8,g9},M=9.隐状态转移概率矩阵和混淆矩阵如下:
$\begin{align} & A=\{{{h}_{ij}}\}=\left\{ \begin{matrix} \begin{matrix} 0 \\ 0 \\ \end{matrix} & amp;\begin{matrix} 0.3 \\ 0 \\ \end{matrix} & amp;\begin{matrix} 0.7 \\ 0.4 \\ \end{matrix} & amp;\begin{matrix} \begin{matrix} 0 & amp;0 \\ \end{matrix} \\ \begin{matrix} 0.2 & amp;0 \\ \end{matrix}.4 \\ \end{matrix} \\ 0 & amp;0 & amp;0 & amp;\begin{matrix} 0.45 & amp;0.55 \\ \end{matrix} \\ 0 & amp;0 & amp;0 & amp;\begin{matrix} 0 & amp;1 \\ \end{matrix} \\ 0 & amp;0 & amp;0 & amp;\begin{matrix} 0 & amp;0 \\ \end{matrix} \\ \end{matrix} \right\}, \\ & B=\{{{b}_{ik}}\}= \\ & \left\{ \begin{matrix} 0.4 & amp;0 & amp;0 & amp;0.6 & amp;0 & amp;0 & amp;0 & amp;0 & amp;0 \\ 0 & amp;0 & amp;0 & amp;0 & amp;0.5 & amp;0 & amp;0.5 & amp;0 & amp;0 \\ 0 & amp;0 & amp;0 & amp;0 & amp;0.35 & amp;0.6 & amp;0 & amp;0.05 & amp;0 \\ 0 & amp;0 & amp;1 & amp;0 & amp;0 & amp;0 & amp;0 & amp;0 & amp;0 \\ 0 & amp;0 & amp;0 & amp;0 & amp;0 & amp;0.3 & amp;0 & amp;0 & amp;0.7 \\ \end{matrix} \right\} \\ \end{align}$
![]() | Fig. 1 Instance of trajectory hidden Markov model图 1 隐马尔可夫轨迹模型示例 |
为了计算矩阵A,需要将一条完整轨迹按固定长度分段,检查一段轨迹上的所有点pk是否都落在S中:如果“是”,则记下该隐状态的标识符;否则,记为“空”(用*表示).因此,子轨迹能够表示为一个状态序列如s1s3**s5.然后,根据该序列计算轨迹状态转移概率.矩阵B的计算通过检查pk是否属于gi和si或只属于gi.在使用隐马尔可夫模型对位置大数据进行建模前,需要对轨迹数据进行预处理工作,包括去噪、分段等处理,使轨迹数据得以化简,进而更加符合大数据环境下隐马尔可夫模型建模要求,下一节将详细阐述.
3 位置大数据预处理面对体量和规模巨大、时空关系复杂的位置大数据,常常需要对轨迹数据进行预处理,本文提出新型大数据环境下位置密度分区和局部位置数据特征提取方法.位置密度分区通过抽取对象的特征,将特征相似的多个对象组成类,达到对数据分类与筛选的目的,同一个类中的成员对象具有相似的属性.DBSCAN(density-based spatial clustering of applications with noise)算法[17]是基于密度对数据进行处理的经典方法.轨迹数据在时空上的分布具有不确定性,某些区域轨迹点较密集,而一些区域较分散,因此特别适合于利用基于密度的分区方法对轨迹数据进行预处理.位置大数据分析通常需要基于路网或地图数据展开,因此需要将连续的空间地图进行离散化处理,划分为多个区域,依位置密度分布进行分区是一种较好的策略.
本文对DBSCAN聚类算法进行改进,用于对位置大数据进行分区划分.此外,位置大数据包含了大量信息,用途广泛,在使用这些数据之前,有诸多问题需要解决:(1) 连续的轨迹数据由定位终端按一定时间间隔采样而成,由于移动对象速度变化、定位设备精确度不高等因素的影响,生成的轨迹数据不完全符合真实情况;(2) 由于数据采集设备自身的不稳定性,采集到的数据集中往往会包含一些噪声点.
3.1 位置密度分区本节使用基于密度的方法对位置大数据分区,将轨迹点与由信号不稳定等因素产生的噪声点区分,其本质是轨迹聚类,输入参数包括邻域半径ε和最小轨迹点数Θ.下面给出位置密度分区算法中的主要概念.
定义7(ε邻域). 给定轨迹点p,将半径为ε内的区域称为该轨迹点的ε邻域.
定义8(核心对象). 给定轨迹点p的ε邻域内轨迹点数量大于Θ,则对象p被称为核心对象.
位置密度分区法的基本思想是:通过遍历位置大数据中每个轨迹点的邻域来生成簇.假设一个轨迹点p的邻域内包含多于Θ个轨迹点,则创建一个新簇,将p作为该簇的核心对象.然后,递归遍历核心对象直接密度可达的对象,过程中包含簇合并操作.当没有新的点可以添加到任何簇时,过程结束.算法如下所示[14]:
算法1. 基于密度的轨迹分区算法——TraCluster.
输入:位置大数据集合S,聚簇半径ε,最少轨迹点数目Θ.
输出:达到聚类密度要求的簇集合C={C1,C2,…,Cn}.
1. n=0; | //初始化簇的个数为0 |
2. for each unvisited point p in S | |
3. report p as visited; | //将 p标记为已访问 |
4. R= Neighbours( p, ε); | |
5. if size( R)< θ then | |
6. report p as noise; | //如果满足 size( R)< θ,则将 p标记为噪声 |
7. else | |
8. create a new cluster C k; | //建立新簇 C k |
9. ExpandCluster( p, N, C k, ε, θ); | |
10. end if | |
11. end for |
在位置密度分区过程中,遍历S中所有轨迹点.如算法1所示.
• 第1行初始化簇个数;
• 第2行开始遍历轨迹点;
• 第3行将轨迹点p标记为已访问;
• 第4行~第11行计算轨迹点p与其他所有轨迹点的距离,将距离小于ε的轨迹点存入集合R中,如果R小于Θ,则标记p为噪声;否则,以p为核心建立新簇,并调用ExpandCluster算法递归访问R中轨迹点.
ExpandCluster算法如下所示:
算法2. ExpandCluster(p,R,Ck,ε,Θ).
输入:p的邻域R,簇Ck,半径ε,最少轨迹点数目Θ.
输出:所有达到密度要求的簇.
1. add p to cluster C k; | //首先将核心点加入簇 C k |
2. for each point p' in R | |
3. report p' as visited; | |
4. R'= Neighbours( p', ε); | //对 R邻域内的所有点在进行半径检查 |
5. if size( R')≥ θ then | |
6. add p' to cluster C k; | //将 p'加入簇 C k |
7. ExpandCluster( p', R', C k, ε, θ); |
在ExpandCluster算法中,深度优先访问R中轨迹点:
• 首先,将点p加入簇Ck中(第1行);
• 将p标记为已访问(第2行、第3行);
• 判断R中点p'是否为核心对象,如果是,则将p'的加入到Ck中(第4行~第6行);执行递归运算深度优先判断p'的ε邻域(第7行).
借助空间索引技术DISC-tree[7],本文所提位置密度分区算法时间复杂度为O(nlogn),n表示数据集中对象的数目.
3.2 局部位置数据特征提取为了提高模型对位置大数据预处理的效率,需要对轨迹数据进行分段处理,轨迹分片称为“段”.本文提出轨迹分段算法TraSegment[14],用于提取局部位置数据特征.轨迹分段算法仍然使用ε邻域中半径ε作为参数.如算法3所示,轨迹分段的基本思想是:
(1) 进行初始化工作(第1行).
(2) 遍历TraList中的轨迹点(第2行),判断其是否已经访问过:如果访问过,则跳过;否则,设置其段id为c,并标记为已访问(第3行~第8行).从j=i+1点开始遍历,如果点i与点j的距离小于ε,则标记点j的段id为c,并标记为已访问(第9行~第17行).
算法3. 局部位置数据特征提取算法——TraSegment.
输入:包含n个位置点的轨迹链表TraList,分段半径ε.
输出:所有达到分割要求的簇.
1. c=1; | //设置段 id从1开始 |
2. for i=0 to TraList.length | //遍历轨迹点 |
3. if TraList[ i]. visit==true | |
4. continue; | |
5. else | |
6. TraList[ i]. id= c; | |
7. TraList[ i]. visit=true; | //将 p标记为已访问 |
8. end if | |
9. for j= i+1 to TraList. length | |
10. if TraList[ i]. visit==true | |
11. continue; | |
12. end if | |
13. if Distance( TraList[ i], TraList[ j])< ε | |
14. TraList[ j]. id= c; | |
15. TraList[ j]. visit=true; | //将 p标记为已访问 |
16. end if | |
17. end for | |
18. c++; | |
19. end for |
轨迹分段的效果如图 2所示,可以看到:通过应用分段算法,连续的轨迹被划分成为不同的片段.轨迹分段主要作用是:(1) 进一步进行位置大数据预处理操作,提高预测算法运行效率;(2) 位置大数据中局部位置数据的特征提取,提取轨迹预测模型中包含的隐状态.
![]() | Fig. 2 Results of trajectory segmentation图 2 轨迹分段结果 |
本文提出的基于HMM的参数自适应轨迹预测算法以隐马尔可夫模型为基础,根据轨迹数据特点,抽象出与HMM模型相对应的轨迹隐状态和观测状态,通过解决HMM模型的3个基本问题完成轨迹预测.
轨迹模型建立过程中,首先使用轨迹分段算法对位置大数据进行分段,将用于训练的轨迹划分成不同的段,用{Ci,i=1,2,3,...,M}表示.同时,为简化轨迹,仍采用平面区域划分方法将轨迹表示成为网格序列,用{gi,i=1,2,3,..., N}表示.如图 3所示,轨迹T1被表示为{g1,g2,g6,g10}和{C1,C2,C3},前者为网格序列表示,后者为段序列表示,分别对应HMM中的观测状态和隐状态;M和N分别为隐状态和观测状态数量.
![]() | Fig. 3 Answer-Loss problem图 3 答案丢失问题 |
由图 3可知,将移动空间按照网格划分后,原本相近的两条轨迹被划分到了不同的网格单元中.在完整的空间区域中,可以通过计算两条轨迹中轨迹点的欧氏距离来判断轨迹相似度.但是在经网格划分后的空间中,欧式距离这一参数将不再具有意义,取而代之的是网格的邻近程度.轨迹T1和轨迹T2在网格划分下相似度较低,因为两条轨迹分属于不同的网格,邻近程度较低,使得原本相似度较高的两条轨迹在网格划分下表现出较低相似度,这就是答案丢失问题.为了解决这一问题,引入HMM模型,轨迹T1和T2同属于{C1,C2,C3}这一段序列.在HMM模型中,使用段序列表示轨迹,通过判断不同轨迹是否属于相同的段来计算轨迹相似度,这样就避免了直接使用网格序列表示轨迹带来的答案丢失问题.
应用所提出的轨迹模型,当输入轨迹由网格序列表示时,可以根据HMM的解码问题求解出对应的隐状态序列,即得到由段序列表示的轨迹.解码问题是HMM中的一个重要问题.在很多情况下,人们对隐状态更感兴趣,因为其包含了一些不能被直接观察到的有价值的信息.对于轨迹预测模型而言,当已知一个网格表示的轨迹序列Rm(t)={gi|i=1,2,3,...,N}时,如何根据已知的λ={π,A,B}求出对应最可能的段序列Cn(t)={Ci|i=1,2,3,...,M},这一问题就可以转化为HMM中的解码问题,也是本文解决的一个关键问题.已知:
$P(R|\lambda )=\underset{{{C}_{1}},{{C}_{2}},...,{{C}_{M}}}{\mathop{\text{Max}}}\,P(R|{{C}_{i}},\lambda )P({{C}_{i}}|\lambda )$ | (2) |
在进行轨迹预测时,首先将输入的轨迹投影到网格划分的区域中.如图 4所示,将轨迹l1转换为网格序列{g1,g2,g6,g10}表示,然后使用训练轨迹数据建立预测模型,应用评估问题,针对不同轨迹进行预测,即:将轨迹带入模型计算每条轨迹出现的概率,取概率最大轨迹作为预测轨迹T(t),t=n+1时刻对应的隐状态,并用轨迹分段的中心点代表预测点.
![]() | Fig. 4 Example of trajectory projection图 4 轨迹投影举例 |
针对评估问题,其具体描述为:已知模型参数λ={π,A,B},如何计算某一给定网格表示的轨迹序列Rm(t)={gi| i=1,2,3,...,N}在模型中出现的概率,即P(O|λ).从HMM模型的定义中可以知,对于某一给定的观测序列,其对应的隐状态序列存在NL个,其中,L为给定观察序列的长度.假设网格序列R对应的段序列集为Ci=(C1,C2,…,CT), T=NL,则:
$P(R|\lambda )=\sum\limits_{{{C}_{1}},{{C}_{2}},...,{{C}_{T}}}{P(R|C,\lambda )}P(C|\lambda )$ | (3) |
上式的计算时间复杂度为O(NL).为了高效地解决轨迹预测的评估问题,本文使用前向算法Forward[18].
4.2 自适应参数选择算法基于HMM的轨迹预测模型的建立,与轨迹所在区域的大小密切相关.轨迹预测模型的精确度受3个主要因素的影响,分别是HMM模型参数λ={π,A,B}、区域划分网格大小gridSize、轨迹分段的大小segSize.在实际应用中,用于预测的已知轨迹可能并不符合预测模型的理想输入.由于移动对象的运动速度可以是随机变化的,所以任意两个轨迹点的间隔大小并不一致,这就导致了由输入产生的隐状态链不连续的问题.为了解决这一问题,本文首先计算得到轨迹点最小间隔minDis,然后对包含轨迹的区域进行整体缩放,使得minDis小于模型的聚类最小间隔.这样得到的轨迹最小间隔满足了模型要求,但是缩放使其他间隔加大,所以要通过添加轨迹中间点填补空缺轨迹.上述过程被称为参数自适应选择,如算法4所示,基本思想为:第1行根据输入位置点计算轨迹点最小间隔;第2行、第3行将参数segSize,gridSize设置为最小间隔;第4行~第10行判断轨迹间隔是否过大,如超过要求,则应用线性插值法插入轨迹点,进行轨迹补全;第11行、第12行将补全后的位置点投影到网格上表示;第13行返回结果.将算法4得到的网格序列作为输入数据带入第4.4节介绍的轨迹预测算法,便可以获得一条连续的预测轨迹序列.自适应参数选择算法主要解决真实应用中用于预测的轨迹间隔随机变化导致预测准确率下降的问题,根据输入轨迹自动选取参数组合,避免了隐状态不连续、状态停留问题.
算法4. 自适应参数选择算法——ParaSelection(S).
输入:任意一条轨迹序列S={s1,s2,…,sn}.
输出:轨迹在网格上的投影点集.
1. min=MinDistance(S); //取输入轨迹点最小间隔
2. segSize=min;
3. gridSize=min;
4. distance=0;
5. for (i=0; i<S.length-1; i++)
6. distance=Distance(S[i],S[i+1]);
7. if (distance>segSize×2)
8. Insert();
9. end if
10. end for
11. for each p in S
grid.add(TransfromIntoGrid(p)); //格式化网格序列
12. end for
13. return grid;
4.3 模型存在问题解决[14](1) 隐状态链不连续问题解决策略
HMTP模型中包含一条马尔可夫链,用来描述段的状态转移.如图 5所示,圆点表示训练数据,三角形表示预测时的前n个轨迹点.段簇的转移代表着移动对象的轨迹序列,当网格大小不大于段的大小时,预测所需的轨迹点所处的网格对应轨迹形成的段,因此,从图 5中可以发现,轨迹对应的隐状态是连续的.
![]() | Fig. 5 Trajectory state transfer图 5 轨迹状态转移 |
然而当网格过大时,模型对于轨迹的描述就变得不尽如人意了.如图 6所示,网格尺寸相对于段过大,当判别轨迹点p对应的隐状态时出现偏差.因为对于p所处的网格R来说,其对应的概率最大段为C,因此,落入R中的轨迹点都将被视为属于C,而忽略p本应属于的段,导致隐状态链不连续.
![]() | Fig. 6 Prediction bias caused by zooming in the grid图 6 放大网格导致的预测偏差 |
由于预测算法使用前向算法进行计算,在状态转移矩阵时,不连续的隐状态之间的转移概率为0,因此前向算法的输出结果必为0,导致预测失败.解决这一问题的一种方法是:在HMM建立过程中分梯度多次使用历史数据,即,以周期{T|T>gridSize}为间隔提取轨迹点,与连续轨迹点一同进行状态转移矩阵生成计算,得到的矩阵中不连续的两个状态的转移概率.这样,前向算法的计算得以进行.
(2) 状态停留问题解决策略
状态停留问题指的是前n个轨迹点中,存在多个连续点属于同一段的情况.由状态转移矩阵的定义可知, pii=0,1≤i≤N,导致预测失败.因此,可以通过实验获得pii的经验值,使得0<pii<max{pij},1≤i≤N,1≤j≤N.
4.4 基于HMM的自适应轨迹预测算法计算出HMM模型的具体参数λ={π,A,B}以后,下一步进行轨迹预测.
预测问题具体描述为:已知给定轨迹{T(t),t=1,2,3,…,n}和λ,求解T(n+1).对此,可以通过枚举下一轨迹点可能属于的网格,组成轨迹{T(t),t=1,2,3,…,n+1},带入评估问题,求解出概率最大的预测轨迹,其T(n+1)即为预测点.但是这种方法求得的预测点为网格,预测精度较低.因此,改进评估问题中的前向算法,使得输出为概率最大的段中心点,这样得到的预测点精度相比网格有很大提高,轨迹预测算法如下所示.
算法5. 基于HMM的自适应轨迹预测算法.
输入:任意一条轨迹序列S={s1,s2,…,sn}.
输出:预测点坐标(xn+1,yn+1, tn+1).
1. TraCluster( S); | //基于密度的轨迹分区预处理 |
2. TraSegment( S); | //轨迹分段 |
3. ParaSelection( S); | //初始化操作 |
4. Transform( S); | //将位置大数据转化为网格链表示 |
5. if ( S. length≤0) | |
6. return false; | |
7. Forward( S, tail); | //使用前向算法计算 |
8. Viterbi( S, bestSeq); | //使用维特比算法计算输入轨迹对应的隐状态 |
9. for i=0 to m | //遍历预测时刻隐状态 i, m表示隐状态转移概率矩阵维度 |
10. if ( i in bestSeq) | //防止位置点倒退,去除已遍历过的段 |
11. continue; | |
12. for j=0 to m | //遍历 tail中隐状态 j |
13. next[ i]+= tail[ j]× M[ j, i]; | //计算由状态 j转换为 i的概率, M表示隐状态转移概率矩阵 |
14. max= i while next is Maximum; | // i值即为预测的隐状态的段号 |
15. return ( x n +1, y n +1, t n +1); | //第 i个段中心点的坐标 |
轨迹预测算法基本思想是:
•针对每条轨迹序列,第1行进行位置密度分区操作,第2行进行局部位置数据特征提取,第3行调用自适应参数选择算法;
•第4行进行轨迹投影将输入原始轨迹转化为由网格表示的网格链;
•第5行、第6行判断输入轨迹长度是否符合要求;
•第7行、第8行分别使用前向算法和维特比算法得到前向变量的最后一列tail,用以保存当前状态到各个状态的状态转移概率和最佳隐状态序列bestSeq;
•第10行、第11行判断隐状态是否已经在bestSeq中出现,如果出现,则舍弃该隐状态;
•第12行~第14行计算从S序列对应的最后一个隐状态转移到所有遍历的隐状态的概率;
•当概率取最大时,第15行返回对应的隐状态的中心点,作为预测的目标点.
5 实验与性能分析 5.1 实验环境及数据集描述本实验所使用的位置大数据集来自微软亚洲研究院的GeoLife项目[19],由182个用户5年间11 129天的GPS轨迹数据组成.数据集包含17 621条轨迹数据,包含23 667 828个轨迹点,总长度达1 292 951公里.数据分布在北京市主城区的各条街道上,很好地反映了移动对象的真实运动行为.这些GPS数据由不同的定位装置收集,包括GPS接收机、智能手机等.文中算法利用C#程序设计语言实现,使用Microsoft Visual Studio 2008开发环境,数据库为SQL Server 2008.实验硬件平台为:Intel(R) Core(TM) 2 Duo P8700 2.53GHz CPU,2GB内存,操作系统平台为Windows 7.实验中实现了基于HMM的自适应轨迹预测算法SATP、朴素的基于HMM未考虑参数自适应选择的预测算法Naïve、基于时间连续贝叶斯网络的轨迹预测算法PutMode[20]、基于高斯混合模型的轨迹预测算法GMTP[15].为方便比较不同算法的优劣,性能评价指标定义如下:
定义9(预测命中率). 已知轨迹序列T={T1,T2,…,Tk},预测轨迹序列Tp={Tp1,Tp2,…,Tpn},k<n,dist(p,q)表示空间中p和q两点间的欧氏距离,η表示距离阈值,则dist(Ti,Tpi)<η时表示预测命中,定义如下:
$H({{T}_{i}},T{{p}_{i}})=\left\{ \begin{array}{*{35}{l}} 1,\text{ }(dist({{T}_{i}},T{{P}_{i}})<\eta ) \\ 0,\text{ }(dist({{T}_{i}},T{{P}_{i}})>\eta ) \\ \end{array} \right.$ | (4) |
定义10(预测准确率). 已知轨迹序列T,预测轨迹序列Tp,则预测准确率定义为
$Accuracy=\frac{\sum\limits_{i=1}^{n}{H({{T}_{i}},T{{p}_{i}})}}{|Tp|}$ | (5) |
其中,|Tp|代表预测轨迹序列的长度.
定义11(预测偏离度). 已知轨迹序列T,预测轨迹序列Tp,则预测偏离度定义为
$Deviation=\frac{\sum\limits_{i=1}^{n}{dist({{T}_{i}},T{{p}_{i}})}}{|Tp|}$ | (6) |
由于实际应用中移动对象的移动速度随机变化,不是恒定不变的,导致轨迹间隔大小不一,本文提出了基于HMM的自适应轨迹预测算法SATP.下面对比在移动对象速度随机变化,即,轨迹间隔大小不一的情况下,不同算法的预测效果.实验采用8组不同轨迹数据集,每组数据包含2.5×106个轨迹点作为训练数据.
通过对比各种算法在不同位置大数据集上的预测精度可以发现(如图 7所示):
![]() | Fig. 7 Accuracy comparison in randomly changing speed图 7 变速情况下预测准确率比较 |
(1) 与Naïve和PutMode算法相比,SATP模型在不同数据集上进行轨迹预测具有较高的预测准确率,平均高出Naïve算法46.7%,高出PutMode算法27.2%,高出GMTP算法10.3%.因为SATP模型总结了基于HMM的轨迹预测模型中隐状态的特点,针对朴素算法参数恒定、无法灵活应对速度随机变化的移动对象这一不足提出改进,对输入轨迹动态判断,分析移动对象的速度特征,在每次预测前调整网格大小和分段大小这两个参数,使得模型更加适合对移动对象当前和未来位置的预测.
(2) 当面对不同轨迹数据时,Naïve,PutMode和GMTP算法的预测准确率不稳定,时高时低,而SATP模型相对稳定,准确率维持在较高水平.主要是因为其可以针对不同数据自适应选择参数,更好地适应不同轨迹数据.
图 8给出了不同预测算法在8组实验数据上的预测偏离度:当预测轨迹由速度随机变化的移动对象产生时, SATP模型的参数也相应改变,其中分段大小这一参数也随之改变.根据预测偏离度的定义可知,分段大小会对偏离度产生一定影响,因此,SATP模型的预测偏离度低于Naïve和PutMode算法.此外可以发现,GMTP算法的预测偏离度略高于SATP.其原因在于:GMTP算法利用高斯过程回归预测移动对象最可能的运动轨迹,可以较为准确地预测对象未来位置点.此外可以发现,PutMode算法的预测偏差最高.原因在于其采用不确定性轨迹圆柱的圆心代表预测点,位置预测的精度较低,因此会产生较大的预测偏差.
![]() | Fig. 8 Deviation comparison in randomly changing speed图 8 变速情况下预测偏离度比较 |
对于速度变化较低的移动对象,通过实验对比4种算法的预测效果.实验结果如图 9和图 10所示.通过观察可以发现:在同样的数据集上,尽管SATP和Naïve算法的预测准确率接近,但是SATP模型的预测偏离度低于Naïve算法,轨迹预测的精度明显好于PutMode和GMTP算法.实验结果表明:对速度相对恒定的移动对象进行轨迹预测时,SATP和Naïve算法的准确率相近.因为匀速运动的移动对象的轨迹间隔变化较小,不需要更改分段大小或网格大小就可以获得较高的预测准确率.由于SATP模型动态改变分段大小,使得在距离上的划分粒度更为精细,因此预测点与实际点的距离较小,预测偏离误差低于恒参的Naïve算法.
![]() | Fig. 9 Accuracy comparison in constant speed图 9 恒定速度情况下预测准确率比较 |
![]() | Fig. 10 Deviation comparison in constant speed图 10 恒定速度情况下预测偏离度比较 |
综上所述,SATP模型综合考虑移动对象真实状态特征,根据输入数据自主调整模型参数,改善预测效果.与未考虑参数自适应调节的Naïve,PutMode和GMTP算法相比,该模型在预测准确率和偏差上具有明显的优势.
5.4 模型建立和预测时间效率分析在实时预测系统中,轨迹预测模型的建立与预测效率十分重要,为了验证轨迹预测模型建立和预测的时间效率,选取10 000条训练轨迹进行实验,其中每条轨迹包含约10万个轨迹点.本节比较4种算法的运行时间性能,结果如图 11所示.
![]() | Fig. 11 Execution time of creating TP model and prediction图 11 轨迹预测模型建立及预测时间 |
实验结果表明:
(1) GMTP算法略优于SATP算法.与SATP和PutMode算法相比,Naïve算法模型建立及预测时间较低,随训练轨迹增多而增大,保持线性关系.当训练轨迹数量达到10 000条(大约109个轨迹点)时,模型建立及预测时间为3s,实时性较好.
(2) SATP模型的建立及预测时间与Naïve算法相比有所增加,但总体消耗保持在可以接受的范围内,而且与训练轨迹数量成线性关系.SATP模型的建立及预测时间稍高于Naïve算法的原因在于:轨迹预处理阶段,由于模型参数动态变化,需要自适应选择参数,导致需要进行区域划分和轨迹分段,耗费了时间;反观Naïve算法,由于参数恒定,所以模型建立一次性完成,只进行一次轨迹分段与区域划分,因此时间开销较小.实验中为了展示算法的特点,选取极端情况下的数据进行测试,而真实情况下移动对象的速度变化相对缓慢.也就是说,SATP模型的参数变化频率较低,因此时间开销相对较小.与PutMode算法相比,其运行时间平均降低38.9%.其原因在于:PutMode算法进行预测时需要花费大量时间构建时间连续贝叶斯网络;此外,其轨迹聚类操作的时间开销较大,不适合位置大数据的预测.
(3) 在算法的真实应用中,用于训练的轨迹数据远小于实验数据,在位置大数据环境下(1万条轨迹,大约100万个轨迹点),SATP模型建立和预测时间均保持在20s以下,具有良好的实时性.
6 结论与展望本文旨在研究位置大数据环境下的新型高效、准确的轨迹预测方法,为了对位置大数据进行高效预处理,提出了基于密度的聚类方法用于位置密度分区和轨迹分段局部位置数据特征提取方法.为了适应移动对象速度动态变化的真实运动场景,提出了一种自适应参数选择算法,基于该算法设计实现了新型轨迹预测模型.位置大数据集上的实验,验证了本文所提自适应轨迹预测方法的时间性能优势及预测的精准性.最后,在基于HMM的轨迹预测模型基础上,设计实现了一个移动对象轨迹预测系统,综合运用前文所述理论依据,提供了对预测结果直观的展示,系统详细描述请参见http://userweb.swjtu.edu.cn/Userweb/qiaoshaojie/demo2.htm.
未来的研究工作包括:(1) 充分考虑客观因素对移动对象位置预测的影响,如,红绿灯、天气等因素,提高预测算法对环境因素的自适应性;(2) 为了进一步提高算法运行效率,将本文提出的轨迹预测算法并行化,将算法移植到Hadoop平台,为公安部门提供实时交通流监控和预测,辅助智能交通控制.
[1] | Zhang JP, Wang FY, Wang KF, Lin WH, Xu X, Chen C. Data-Driven intelligent transportation systems: A survey. IEEE Trans. on Intelligent Transportation Systems, 2011,12(4):1624-1639 . |
[2] | Wang ZC, Lu M, Yuan XR, Zhang JP, Wetering H. Visual traffic jam analysis based on trajectory data. IEEE Trans. on Visualization and Computer Graphics, 2013,19(12):2159-2168 . |
[3] | Meng XF, Ding ZM. Mobile Data Management: Concepts and Techniques. Beijing: Tsinghua University Press, 2009. 185-200 (in Chinese). |
[4] | Guo C, Liu JN, Fang Y, Luo M, Cui JS. Value extraction and collaborative mining methods for location big data. Ruan Jian Xue Bao/Journal of Software, 2014,25(4):713-730 (in Chinese with English abstract). |
[5] | Calabrese F, Pereira FC, Francisco C, Di Lorenzo G, Liu L, Ratti C. The geography of taste: Analyzing cell-phone mobility and social events. In: Proc. of the 8th Int’l Conf. on Pervasive Computing. Springer-Verlag, 2010. 22-37 . |
[6] | Calabrese F, Smoreda Z, Blondel VD, Ratti C. Interplay between telecommunications and face-to-face interactions: A study using mobile phone data. PLoS ONE, 2011,6(7):e20814 . |
[7] | Qiao SJ, Han N, Wang C, Zhu F, Tang CJ. A two-tiered dynamic index structure of moving objects based on constrained networks. Chinese Journal of Computers, 2014,37(9):1947-1958 (in Chinese with English abstract) . |
[8] | Monreale A, Pinelli F, Trasarti R, Giannotti F. WhereNext: A location predictor on trajectory pattern mining. In: Proc. of the 15th ACM SIGKDD Int’l Conf. on Knowledge Discovery and Data Mining. New York: ACM Press, 2009.637-646 . |
[9] | Ying JJ, Lee W, Weng T, Tseng VS. Semantic trajectory mining for location prediction. In: Proc. of the 19th ACM SIGSPATIAL Int’l Conf. on Advances in Geographic Information Systems. New York: ACM Press, 2011. 34-43 . |
[10] | Song MB, Ryu JH, Lee SK, Hwang CS. Considering mobility patterns in moving objects database. In: Proc. of the 2003 Int’l Conf. on Parallel Processing. Washington: IEEE, 2003. 597-604 . |
[11] | Ishikawa Y, Tsukamoto Y, Kitagawa H. Extracting mobility statistics from indexed spatio-temporal datasets. In: Proc. of the 2nd Int’l Workshop on Spatio-Temporal Database Management. 2004. 9-16. |
[12] | Asahara A, Maruyama K, Sato A, Seto K. Pedestrian-Movement prediction based on mixed Markov-chain model. In: Proc. of the 19th ACM SIGSPATIAL Int’l Conf. on Advances in Geographic Information Systems. New York: ACM Press, 2011. 25-33 . |
[13] | Gambs S, Killijian M, Cortez DP, Miguel N. Next place prediction using mobility Markov chains. In: Proc. of the 1st Workshop on Measurement, Privacy, and Mobility. New York: ACM Press, 2012. 1-6 . |
[14] | Qiao SJ, Shen DY, Wang XT, Han N, Zhu W. A self-adaptive parameter selection trajectory prediction approach via hidden Markov models. IEEE Trans. on Intelligent Transportation Systems, 2015,16(1):284-296 . |
[15] | Qiao SJ, Jin K, Han N, Tang CJ, Gesangduoji, Gutierrez LA. Trajectory prediction algorithm based on Gaussian mixture model. Ruan Jian Xue Bao/Journal of Software, 2015,26(5):1048-1063 (in Chinese with English abstract). |
[16] | Jensen CS, Lin D, Ooi BC, Zhang R. Effective density queries on continuously moving objects. In: Proc. of the 22nd Int’l Conf. on Data Engineering. Washington: IEEE, 2006. 71-71 . |
[17] | Ester M, Kriegel HP, Sander J, Xu X. A density-based algorithm for discovering clusters in large spatial databases with noise. In: Proc. of the 2nd Int’l Conf.on Knowledge Discovery and Data Mining. AAAI Press, 1996. 226-231 . |
[18] | Rabiner LR. A tutorial on hidden Markov models and selected applications in speech recognition. Proc. of the IEEE, 1989,77(2): 257-286 . |
[19] | Zheng Y, Xie X, Ma WY. Geolife: A collaborative social networking service among user, location and trajectory. IEEE Data Engineering Bulletin, 2010,33(2):32-40. |
[20] | Qiao SJ, Tang CJ, Jin HD, Long T, Dai SC, Ku YC, Chau M. PutMode: Prediction of uncertain trajectories in moving objects databases. Applied Intelligence, 2010,33(3):370-386 . |
[3] | 孟小峰,丁治明.移动数据管理:概念与技术.北京:清华大学出版社,2009.185-200. |
[4] | 郭迟,刘经南,方媛,罗梦,崔竞松.位置大数据的价值提取与协同挖掘方法. 软件学报,2014,25(4):713-730. |
[7] | 乔少杰,韩楠,王超,祝峰,唐常杰.基于路网的移动对象动态双层索引结构.计算机学报,2014,37(9):1947-1958 . |
[15] | 乔少杰,金琨,韩楠,唐常杰,格桑多吉,Gutierrez LA.一种基于高斯混合模型的轨迹预测算法. 软件学报,2015,26(5):1048-1063. |