在城市交通流特性研究中, 机动车出行轨迹包含详细的交通流微观参数, 通过对城市交通路网中所有运行车辆的出行轨迹进行提取、聚类和整合分析, 可以系统科学地再现所有车辆的运行场景, 有效获取城市网络交通流的宏观运行状态, 进而为分析城市交通需求的结构和时空分布特性提供数据支撑.

与其他交通信息采集技术相比, 车牌自动识别系统工作具有连续性强、数据精确度高, 检测样本量大等优点^[1], 国内外越来越多的专家学者运用车牌识别数据进行城市交通特性的深入分析.如孙剑等^[1]提出了车牌自动识别环境下的车辆出行矩阵估计新方法；李晓莉等^[2]提出了基于异常值数据表现以及行程时间分布特征的异常值剔除方法；畅玉皎等^[3]基于车牌识别数据, 通过k-means聚类数据挖掘方法, 提取了路网中的通勤特征车辆, 并研究了通勤车辆出行行为的时空特性.上述研究主要集中在车牌识别数据的概率统计及宏观特性提取上, 对于出行轨迹、出行行为等交通流微观参数的挖掘提取较少涉及.

出行轨迹提取的主要难点是出行轨迹的补全.当前, 轨迹补全的方法主要有2种：最短路径法(The shortest path, SP)及基于限速的拓扑寻优法(topological optimization based on speed constraint, SC+TO).其中, 最短路径法认为出行者在路径决策时总是选择最短路径；基于限速的拓扑寻优法认为出行者的速度不会过高或过低, 轨迹速度与路网最高通行速度及最低通行速度的差的绝对值之和越小；同时, 轨迹在静态拓扑(距离、转弯等)上的指标值越小, 则该轨迹被选择的可能性越大.如Yang等^[4]基于最短路径法的思想, 结合粒子滤波器及路径流量估计重建了车辆的运行轨迹；Quddus等^[5]通过改进的地图匹配算法, 利用最短路径算法完成了稀疏车辆轨迹的补全.王龙飞^[6]则基于车牌识别数据, 利用SC+TO法结合TOPSIS法联合决策完成了车辆出行轨迹的提取及重构.最短路径法算法复杂度较小, 有利于进行模型的快速求解, 但在实际的路径决策时, 受到道路拥堵、信息获取延迟等因素的影响.出行者可能不选择最短路径, 从而造成该算法的准确度较低, 漏检轨迹节点数越多, 表现越差.SC+TO法融合了车牌识别数据的动态指标及静态指标, 多样化的决策属性可以更好地描述出行者的决策行为, 该算法准确度较SP法有了较大的提高, 但算法没有充分挖掘车牌识别数据的历史属性, 忽略了出行者路径决策的主观行为, 如出行者总是偏爱某些路径, 即使它不是最短或最快的.

针对上述路径提取补全算法的局限性, 提出了基于K则最短路径算法(K shortest paths, KSP)及灰色关联法(GRA算法)的出行轨迹补全算法, 通过对车牌自动识别数据的深入挖掘分析, 完成了对城市交通网络中机动车出行轨迹的提取, 并基于贵阳市南明区的实际数据验证了算法的有效性.

1 数据预处理

主要通过车牌识别数据、交叉口点位数据及路段时空数据3种数据来实现城市所有机动车出行轨迹的提取及补全, 受多种因素影响, 数据往往存在冗杂数据及错误数据, 需要进行数据的预处理.

1.1 车牌识别数据

车牌识别数据即所有车辆的检测记录集合, 其包含多种属性, 该算法仅利用其中的部分属性数据, 需要将其他的属性数据作为“冗杂”数据剔除掉, 以提高算法的运行效率.算法利用到的属性数据为车牌号、检测时间、设备点位地址、设备点位编号、进口道方向及车道编号.

在车牌数据采集过程中, 受到设备识别率、信号传输干扰、套牌车等的影响, 如表 1所示, 车牌识别数据存在以下异常.

其中, 由于设备漏检及数据逻辑关系混乱需要借助车辆的轨迹数据进行识别剔除, 故本节的预处理主要针对数据乱码、检测错误2种数据异常情况.对于数据乱码及检测错误, 由于缺乏其他有效的信息对其进行还原, 故采取删除的处理策略, 删除方法如下：

1) 车牌字符长度约束

正常车牌的字符长度为7~9, 字符长度不在该区间的检测车牌号均为异常车牌.

2) 车牌日检测次数约束

将日车牌识别数据根据车牌属性分组, 统计各分组的频数, 即为各车牌Plate的日检测次数N.以贵阳市2016年2月20日的车牌识别数据为例.为车牌日检测次数超过800次的车牌号如表 2所示.由表 2可得, 正常车牌日检测次数的上限为826, 超过上限的均为异常数据.为了计算方便, 可设置车牌日检测次数上限阈值ψ₊为1 000.同时由于检测次数为1的车牌识别数据对于机动车出行轨迹的提取无意义, 可设置车牌日检测次数下限阈值ψ_-为1, 车牌日检测次数在(ψ_-, ψ₊)之间的数据即为正常车牌识别数据.

将同时满足车牌字符长度及车牌日检测次数约束的车牌识别数据保留, 其他数据均删除即可完成车牌识别数据的清洗.

1.2 交叉口点位数据

提取的机动车出行轨迹定义如下：

$ p = \left( {{\rm{nod}}{{\rm{e}}_1} \to {\rm{nod}}{{\rm{e}}_2} \to \cdots {\rm{nod}}{{\rm{e}}_i} \cdots \to {\rm{nod}}{{\rm{e}}_n}} \right). $

(1)

式中：p为机动车某条出行轨迹, node_i, node_i+1为邻接的两交叉口.

由式(1)可以看出, 出行轨迹的基本单位为相邻交叉口之间的路段.

1.2.1 交叉口实际点位聚类

同一个交叉口, 可能存在多个不同时期安装的视频检测设备, 但其表征的均是该交叉口点位node_i, 故需要对视频检测设备基于交叉口进行聚类, 并对每类进行ID编号, 例如图 1的交叉口1-3, 其ID号为1, 2, 3.

1.2.2 交叉口虚拟点位设置

点位实际路网往往存在无视频检测设备的交叉口, 导致出行轨迹在该交叉口完全缺失, 需要在此设置虚拟设备点位, 以便于后期的轨迹补全.如图 1所示, 交叉口4无实际卡口设备, 可在此设置虚拟设备点位, 并给予相应的ID编号, 与上述的交叉口实际点位共同组成交叉口检索矩阵C.

1.3 路段时空数据

路段时空数据包括路段行程时间和路段距离数据.假设车辆经过两相邻交叉口a、b的时刻为t_a、t_b, 则路段ab的行程时间为

$ t = {t_b} - {t_a}. $

(2)

对车牌号及时间戳排序后的所有车牌识别数据进行遍历, 即可得到各路段所有经过车辆的行程时间数据.由于驾驶员行为差异(绕行及中途停车等)、车牌识别有误等原因, 路段行程时间存在异常值, 需要对行程时间中的异常值进行处理, 处理策略采取李晓莉等^[2]提出的行程时间异常值处理方法, 该方法首先通过设置行程时间的上下限值、重复检测判断阈值来过滤掉明显的异常数据；然后设置统计时窗, 对所有统计时窗内的行程时间数据依次以均值和两倍标准差、中位值和平均绝对偏差为条件进行循环过滤, 直至行程时间样本数据不再发生变化, 则异常值处理完成, 处理后得到路网各路段所有经过车辆的行程时间数据.处理前、后的路段行程时间数据分别如图 2(a)、(b)所示.

2 机动车出行轨迹提取

基于车牌识别数据机动车出行轨迹提取主要包括出行链的获取及分离、出行轨迹补全.

2.1 出行链获取及分离

对车牌识别数据的所有记录行按车牌号、检测时间戳排序后, 即可得到每辆车经过的设备点位集合, 为设备点位匹配上相应的交叉口编号ID, 则每个车牌号的ID编号及时间戳的集合即为每辆车的出行链.出行链分离的目的是将同一辆车的出行链拆分为一定次数的出行, 通过车牌识别数据的时空关系来完成出行链分离, 即采取相邻节点间速度结合交叉口邻接矩阵来进行出行链分离位置的确定, 具体步骤如下：

Step 1：将聚类后的交叉口点位映射到路网拓扑上, 建立交叉口邻接矩阵D, 若两交叉口a, b直接相连, 则邻接矩阵元素赋值1, 即D(a, b)=1, 否则D(a, b)=0.

Step 2：提取出行链相邻节点的车流方向信息(包括进口道方向及车道编号)及检测时间差Δt, 计算车道方向匹配的最短行驶路径s及其距离d, 得到出行链相邻节点的速度v=d/Δt.

Step 3：获取最短行驶路径s的速度上限阈值v_u、下限阈值v_l.

Step 4：出行链分离

出行链相邻节点间的4类情况, 其特性如表 3所示.

对所有出行链执行上述4个步骤, 即可确定所有出行链的分离位置及漏检位置, 在分离位置进行分割操作即可完成出行链的分离.出行链分离完成后, 为获得完整的机动车出行轨迹, 需要设计相应的策略在出行链的漏检位置处进行出行轨迹的补全.

v_l, v_u的确定方法如下：

1) 设定统计时窗t_w.交通流在一天之内的波动性导致了行程时间随时间的波动, 因此设定统计时窗, 认为同一时窗内的行程时间数据具有相同的模式.结合城市交通在时间上的变化特点, 将00:00-06:00的t_w设置为30 min, 06:00-24:00的t_w设置为15 min^{[3, 7-8]}.

2) 对统计时窗内的路段行程时间异常值进行处理后, 取其最大值和最小值作为路段当前时间窗下的行程时间上下限值, 最短行驶路径s的路径行程时间上限t_u, 下限为t_l.计算方式如下：

设最短行驶路径包含的路段为(r₁, r₂, …, r_n), 其路段长度分别为(l₁, l₂, …, l_n), 其中具有行程时间数据的路段长度分别为(l′₁, l′₂, …, l′_m), 路段行程时间在对应时间窗的上下限值为[(t₁^l, t₁^u), (t₂^l, t₂^u), …, (t_m^l, t_m^u)], 则t_l、t_u分别为

$ {t_{\rm{l}}} = \frac{{\sum\limits_{i = 1}^n {{l_i}} }}{{\sum\limits_{i = 1}^m {{{l'}_i}} }}\sum\limits_{i = 1}^m {t_i^1} ,{t_{\rm{u}}} = \frac{{\sum\limits_{i = 1}^n {{l_i}} }}{{\sum\limits_{i = 1}^m {{{l'}_i}} }}\sum\limits_{i = 1}^m {t_i^{\rm{u}}} . $

(3)

式中：t_i^l为路段i的行程时间上限, t_i^u为路段i的行程时间下限.

3) 计算速度上下限阈值分别为

$ {v_1} = \frac{d}{{{t_{\rm{u}}}\xi }},{v_{\rm{u}}} = \frac{{\xi d}}{{{t_{\rm{l}}}}}. $

(4)

式中：ξ为调整系数, 一般设置为1.15, 主要是为了预防机动车实际出行路径非最短路径及行程时间过度滤波等因素造成的速度上下限区间偏小.

2.2 机动车出行轨迹补全

对机动车出行链进行分离之后即可得到出行轨迹, 但由于车流密度大、车速过快、异常天气、设备故障及拍摄角度不对等原因, 造成交叉口设备漏检车辆；同时, 当前城市的视频检测设备无法覆盖所有交叉口, 而完整的出行轨迹要求其相邻节点即为相邻交叉口, 故上述出行链直接分离得到的出行轨迹往往是缺失的, 需要设计相应策略进行补全.

2.2.1 构建出行轨迹可行解集合

为降低算法的运算复杂度, 在进行多目标决策前, 需要构建出行轨迹可行解集合, 构建过程如下：

假设漏检部分的上条记录点为R₁, 下条记录点为R₂.

Step 1：出行者受多种因素的影响, 实际出行轨迹可能不遵循最短路径, 故需要寻求出行OD点之间的多个备选优化路径, 形成最短路径组.运用K则最短路径算法(KSP)^[9]基于距离搜索R₁点与R₂点之间的前K条最短路径, 并根据路径长度升序排列, 构建出行轨迹可行解集合(p₁, p₂, …, p_K).在出行者实际的路径决策过程中, 其待选择路径数一般不超过10, 故KSP算法中的K取10.

Step 2：根据R₁点及R₂点记录的进口道方向、车道编号提取出与R₁连接的路段r_s及与R₂连接的路段r_e, 同时提取待选路径的首尾路段分别为r′_s和r′_e.

Step 3：将出行轨迹可行解集合中$ \left. \begin{array}{l} r{'_{\rm{s}}}{\rm{ = }}{r_{\rm{s}}}\\ r{'_{\rm{e}}}{\rm{ = }}{r_{\rm{e}}} \end{array} \right\}$的路径提取出来, 构建新的出行轨迹可行解集合(p₁^*, p₂^*, …, p_n^*).

2.2.2 出行轨迹决策指标

出行轨迹的补全需要合理模拟出行者的路径决策过程, 出行者在实际出行决策过程中, 其考虑的因素主要有路径距离、路径行程时间、信号交叉口个数(易受信号灯影响增加路径行程时间)^[10]、转弯次数等.通过对上述要素进行多目标寻优, 可以得到两点之间的最优路径, 但出行者由于获取信息的滞后性及固有习惯等因素, 实际选择的路径可能非最优路径^[10-12].针对这种情况, 另外设计2个动态指标来进行补全路径的校正, 分别为行程时间相符程度与轨迹偏好程度.

根据漏检部分上、下条记录的时间戳, 可以获取车辆由R₁点到R₂点的实际路径行程时间, 在确定的设备密度及合理的行程时间扩样方法基础上, 可以计算得到出行轨迹可行解集合中涉及的所有路段的行程时间特征值, 进而得到各轨迹可行解的计算路径行程时间, 轨迹的计算路径时间与实际路径行程时间越接近, 则该轨迹为实际轨迹的概率越大.同时通过统计车牌历史数据中采取各出行轨迹方案的车辆数, 即可得到出行轨迹可行解集合中各路径被选择的概率, 其表征了出行者对于各方案的偏好程度.

综上, 设计了路径距离、行程时间相符程度、轨迹偏好程度、信号交叉口个数、车辆转弯次数5个指标来进行出行轨迹的补全决策.由于出行者对于各指标的增长敏感性是逐渐降低的, 故基于Max-Min方法^[4]采取负指数函数的形式来进行各指标的标准化, 具体如下：

1) 路径距离：

$ {x_1}\left( {{p_i}} \right) = \left\{ {\begin{array}{*{20}{c}} {\exp \left( { - \alpha \left( {\frac{{{L_i} - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{L_j}} \right)}}{{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{L_j}} \right) - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{L_j}} \right)}}} \right)} \right),}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{L_j}} \right) \ne {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{L_j}} \right);}\\ {1,}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{L_j}} \right) = {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{L_j}} \right).} \end{array}} \right. $

(5)

式中：L_i为方案p_i的路径距离；L_j为方案p_j的路径距离; M为方案集合.

2) 行程时间相符程度：

$ {x_2}\left( {{p_i}} \right) = \left\{ {\begin{array}{*{20}{c}} {\exp \left( { - \alpha \left( {\frac{{\left| {{t_i} - {t_ * }} \right| - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {\left| {{t_j} - {t_ * }} \right|} \right)}}{{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {\left| {{t_j} - {t_ * }} \right|} \right) - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {\left| {{t_j} - {t_ * }} \right|} \right)}}} \right)} \right),}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {\left| {{t_j} - {t_ * }} \right|} \right) \ne {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {\left| {{t_j} - {t_ * }} \right|} \right);}\\ {1,}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {\left| {{t_j} - {t_ * }} \right|} \right) = {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {\left| {{t_j} - {t_ * }} \right|} \right).} \end{array}} \right. $

(6)

式中：t_i为方案p_i的计算路径行程时间，t_j为方案p_j的计算路径行程时间, 其为各路段行程时间特征值(取算数平均值)的和；t_*为相邻记录实际路径行程时间.

3) 轨迹偏爱程度：

$ {x_3}\left( {{p_i}} \right) = \left\{ {\begin{array}{*{20}{c}} {\exp \left( { - \alpha \left( {\frac{{{n_{pi}} - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{n_{pj}}} \right)}}{{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{n_{pj}}} \right) - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{n_{pj}}} \right)}}} \right)} \right),}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{n_{pj}}} \right) \ne {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{n_{pj}}} \right);}\\ {1,}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{n_{pj}}} \right) = {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{n_{pj}}} \right).} \end{array}} \right. $

(7)

式中：n_pi为历史数据中方案p_i的车辆数, n_pj为历史数据中方案p_j的车辆数.

4) 信号交叉口个数：

$ {x_4}\left( {{p_i}} \right) = \left\{ {\begin{array}{*{20}{c}} {\exp \left( { - \alpha \left( {\frac{{{N_i} - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{N_j}} \right)}}{{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{N_j}} \right) - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{N_j}} \right)}}} \right)} \right),}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{N_j}} \right) \ne {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{N_j}} \right);}\\ {1,}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{N_j}} \right) = {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{N_j}} \right).} \end{array}} \right. $

(8)

式中: N_i为方案p_i的信号交叉口个数, N_j为方案p_j的信号交叉口个数.

5) 车辆转弯次数：

$ {x_5}\left( {{p_i}} \right) = \left\{ {\begin{array}{*{20}{c}} {\exp \left( { - \alpha \left( {\frac{{{D_i} - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{D_j}} \right)}}{{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{D_j}} \right) - {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{D_j}} \right)}}} \right)} \right),}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{D_j}} \right) \ne {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{D_j}} \right);}\\ {1,}&{{\rm{Ma}}{{\rm{x}}_{j \in M}}\left( {{D_j}} \right) = {\rm{Mi}}{{\rm{n}}_{j \in M}}\left( {{D_j}} \right).} \end{array}} \right. $

(9)

式中：D_i为方案p_i的车辆转弯次数，D_j为方案p_j的车辆转弯次数.

2.2.3 基于GRA的出行轨迹补全决策算法

出行轨迹补全本质上是多目标决策问题, 多目标决策的经典方法主要有层次分析法(AHP)、数据包络分析法(DEA)、逼近于理想解的排序方法(TOPSIS)、灰色关联法(GRA)等^[13-14], 考虑到出行轨迹补全算法的数据源特性和算法特点, 选用灰色关联法^[15]来完成机动车出行轨迹的补全.计算得到出行轨迹可行解集合中所有轨迹的各决策指标值后, 基于GRA算法决策出最优的补全路径, 其具体步骤如下：

Step 1:构建标准化决策矩阵M′

$ \begin{array}{l} \mathit{\boldsymbol{M'}} = \left[ {\begin{array}{*{20}{c}} {{{x'}_{11}}}&{{{x'}_{12}}}& \cdots &{{{x'}_{15}}}\\ {{{x'}_{21}}}&{{{x'}_{22}}}& \cdots &{{{x'}_{25}}}\\ \cdots&\cdots &{{{x'}_{ij}}}& \cdots \\ {{{x'}_{n1}}}&{{{x'}_{n2}}}& \cdots &{{{x'}_{n5}}} \end{array}} \right]\\ {{x'}_{ij}} = \frac{{{x_{ij}}}}{{\sqrt {\sum\limits_{i = 1}^n {{{\left( {{x_{ij}}} \right)}^2}} } }}\left( {1 \le i \le n,1 \le j \le 5} \right). \end{array} $

(10)

式中：(x_i1, x_i2, x_i3, x_i4, x_i5)分别为出行轨迹i的路径距离、行程时间相符程度、轨迹偏爱程度、信号交叉口个数、车辆转弯次数的指标值.

Step 2:构建加权标准化决策矩阵M^*

假设5种属性的决策权重为(ω₁, ω₂, ω₃, ω₄, ω₅), 则加权决策矩阵M^*为

$ \begin{array}{l} {\mathit{\boldsymbol{M}}^ * } = \mathit{\boldsymbol{M'}}{\rm{Diag}}\left( {{\omega _1},{\omega _2}, \cdots ,{\omega _5}} \right) = \\ \left[ {\begin{array}{*{20}{c}} {{\omega _1}{{x'}_{11}}}&{{\omega _2}{{x'}_{12}}}& \cdots &{{\omega _5}{{x'}_{15}}}\\ {{\omega _1}{{x'}_{21}}}&{{\omega _2}{{x'}_{22}}}& \cdots &{{\omega _5}{{x'}_{25}}}\\ \cdots&\cdots &{{\omega _j}{{x'}_{ij}}}& \cdots \\ {{\omega _1}{{x'}_{n1}}}&{{\omega _2}{{x'}_{n2}}}& \cdots &{{\omega _5}{{x'}_{n5}}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {{y_{11}}}&{{y_{12}}}& \cdots &{{y_{15}}}\\ {{y_{21}}}&{{y_{22}}}& \cdots &{{y_{25}}}\\ \cdots&\cdots &{{y_{ij}}}& \cdots \\ {{y_{n1}}}&{{y_{n2}}}& \cdots &{{y_{n5}}} \end{array}} \right]. \end{array} $

(11)

式中：y_ij为ω_j和x′_ij的乘积, 是出行轨迹i的加权指标值.其中灰色关联法权重确定的方法参考文献[15-16].

Step 3:构建差异度矩阵

差异度矩阵为

$ \begin{array}{l} \mathit{\boldsymbol{ \boldsymbol{\varDelta} }} = \left[ {\begin{array}{*{20}{c}} {{\Delta _{01}}\left( 1 \right)}&{{\Delta _{01}}\left( 2 \right)}& \cdots &{{\Delta _{01}}\left( 5 \right)}\\ {{\Delta _{02}}\left( 1 \right)}&{{\Delta _{02}}\left( 2 \right)}& \cdots &{{\Delta _{02}}\left( 5 \right)}\\ \cdots&\cdots &{{\Delta _{0i}}\left( j \right)}& \cdots \\ {{\Delta _{0n}}\left( 1 \right)}&{{\Delta _{0n}}\left( 2 \right)}& \cdots &{{\Delta _{0n}}\left( 5 \right)} \end{array}} \right].\\ {\Delta _{0i}}\left( j \right) = \left| {{{y'}_0}\left( j \right) - {y_{ij}}} \right|. \end{array} $

(12)

式中：Δ_0i(j)为指标i的最大值y′₀(j)与y_ij差的绝对值, 其中y′₀(j)的计算方式如下：

$ \begin{array}{l} {{y'}_0}\left( j \right) = \max {\left[ {{y_{1j}},{y_{2j}}, \cdots {y_{ij}}, \cdots {y_{nj}}} \right]^{\rm{T}}}\\ 1 \le i \le n,1 \le j \le 5. \end{array} $

(13)

Step 4:计算轨迹i的灰色关联度为

$ \begin{array}{l} {\psi _i} = \\ \frac{1}{5}\sum\limits_{j = 1}^5 {\left( {\frac{{\mathop {\min }\limits_i \mathop {\min }\limits_j {\Delta _{0i}}\left( j \right) + \varsigma \mathop {\max }\limits_i \mathop {\max }\limits_j {\Delta _{0i}}\left( j \right)}}{{{\Delta _{0i}}\left( j \right) + \varsigma \mathop {\max }\limits_i \mathop {\max }\limits_j {\Delta _{0i}}\left( j \right)}}} \right)} . \end{array} $

(14)

式中：$ \zeta$为分辨系数, 一般取值为0.5.

Step 5：计算出行轨迹可行解集合中所有轨迹的灰色关联度, 其中, 灰色关联度最大的轨迹方案即为最终轨迹补全方案.

3 实验结果及分析 3.1 数据集

以贵阳南明区某局部路网(1.6 km×1.8 km)进行计算分析, 该区域视频设备点位分布如图 3所示.算法数据集主要包含以下4部分：

1) 历史车牌识别数据P_h.选取图 3区域2016年2月22日到3月18日(工作日, 共20天)的车牌识别数据作为历史数据训练集, 共包含37 928 635条记录.

2) 当前车牌数据P_c.选取上述区域2016年3月22日的车牌识别数据作为实验数据集, 共包含1 833 941条记录, 记录属性与历史车牌数据相同.

3) 道路网数据Node_data, Link_data.运用Python及百度API爬取上述区域的道路网信息, 包含52个节点、82条路段.节点属性为交叉口地址和节点编号, 路段属性为路段编号及路段长度.

4) 设备点位数据Point_data.包含67个实际设备点位、13个虚拟设备点位.基于交叉口聚类后, 得到52个设备点位数据, 每条数据包含点位地址、点位编号、经度、纬度4项属性.

3.2 实验结果

对当前车牌数据集P_c执行预处理算法, 得到新的车牌数据集P′_c, 共包含1 046 586条记录, 删除比例为(1 833 941-1 046 586)/1 833 941=42.93%, 删除的记录主要包括异常数据、出租车数据及全天车牌检测次数为1的数据.对P′_c执行上述出行轨迹提取算法, 运算得到所有车辆当天的轨迹数据集P^*_c, 共包含1 620 872条记录, 去除掉出行长度仅为1条记录的车牌数据, 剩余1 506 985条记录.经过上述步骤, 共提取出有效出行路径208 271条, 其中存在漏检情况的路径有159 017条, 说明漏检情况是大量存在的.

3.2.1 算法有效性分析

出行轨迹提取算法的核心为漏检补全算法, 为验证补全算法的有效性, 从轨迹数据集P_c^*中随机选取m条不存在漏检且节点数不小于10的路径构成标准路径集P_L=(p₁, p₂, …, p_i, …, p_m).针对P_L中的每条路径p_i, 随机删除首末节点间j个连续的节点(j=1, 2, 3, …, n)来模拟漏检节点, 得到待补全路径集P_Lj=(p_j¹, p_j², …, p_jⁱ, …, p_j^m), 对P_Lj执行轨迹补全算法得到补全路径集P′_Lj=(p_j^1′, p_j^2′, …, p_j^i′, …, p_j^m′).建立准确率向量ζ_j=[ζ_j¹, ζ_j², …, ζ_jⁱ, …, ζ_j^m], ζ_i取值方法如下:

$ \zeta _j^i = \left\{ \begin{array}{l} 1,\;\;\;\;p_j^{i'} = {p_i};\\ 0,\;\;\;\;p_j^{i'} \ne {p_i}. \end{array} \right. $

(15)

则补全算法在删除节点数为j时的准确率为

$ {\psi _j} = \frac{{\sum\limits_{i = 1}^m {\zeta _j^i} }}{m} \times 100\% ,m\;取\;300. $

(16)

选取的对比算法如下：

(1) 最短路径法(SP).

(2) 基于限速的拓扑寻优法(SC+TO).

提出的路径补全算法(KSP+GRA)和2种对比算法的准确率表现如图 4所示.由图 4中可以看出, 3种轨迹补全算法的准确率均随删除节点数的增加而下降, 其中最短路径法下降趋势最为明显, 当漏检节点数为6时, 其算法准确率就下降到15%, 这是由其决策属性过于单一造成的.图中表示算法的准确率.随着漏检节点数的增加, 两节点间往往存在多条距离相近的较短路径, 此时出行者可能会根据其他路网信息(转弯、信号灯、路段交通状态)做出非最短路径的出行决策, 故最短路径法的准确率对漏检节点数敏感性较高, 而基于速度约束的拓扑寻优法由于添加了路网拓扑及速度约束属性, 其准确率及鲁棒性较SP有较大提高.相比2种对比算法, 路径补全算法(KSP+GRA)的准确率及鲁棒性均具有明显优势, 当漏检节点数为6的时候, 算法仍具有75.33%的准确率, 且漏检节点数越多, 其算法优势越明显.

3.2.2 漏检情况分析

存在漏检情况的某条出行轨迹的示意图如图 5所示, 由图可知, 出行轨迹上的点位分为3种：正常点位、漏检实际点位、漏检虚拟点位.其中漏检实际点位与漏检虚拟点位统称为漏检点位.漏检实际点位是由视频检测设备未捕捉及识别到车牌引起的, 而漏检虚拟点位则是由交叉口无视频检测设备造成的, 相当于机动车在虚拟点位的漏检率为100%.

测试区域当天所有出行轨迹的漏检情况如表 4所示, 表中，N_d为实际+虚拟点位漏检个数，N_d′为实际点位漏检个数，N_p为路径数，P_d为当前路径数占路径总数的百分比，其中N_d′为1或2的路径对应的P_d和大于80%, 而漏检点个数为1-2的漏检路径往往执行最短路补全策略就可以达到较好的补全效果.随着漏检点个数的增加, 漏检位置之间存在的可行路径数急剧增加, 此时SP算法及SC+TO算法等静态补全策略准确率将急剧下降, KSP+GRA算法的优势开始得到体现.KSP+GRA算法既包含了路径距离、信号交叉口个数、转弯数等静态属性, 又包含了轨迹偏爱程度及行程时间相符程度等动态属性, 使其面对不同的漏检情况都能达到良好的补全效果, 具有较强的鲁棒性.

仅考虑实际点位的漏检, 由表 4可以得到该区域当天的漏检记录数为

$ {N_{\rm{m}}} = \sum\limits_{i = 1}^n {\left( {{N_{{\rm{d'}}}} \times {N_{\rm{p}}}} \right)} = 266067. $

(17)

计算得到视频设备的漏检率(266 067/(266 067+1 046 586)=20.27%), 与视频设备的理论漏检率(1-0.90×9=19%)基本相同.

同时研究了漏检在时间轴上的分布情况, 如图 6所示, 实际+虚拟点位及单独实际点位的出行轨迹全天平均漏检点个数分别为1.769、1.686, 且漏检点个数的比例分布在全天各时段均较为稳定, 说明视频检测设备的识别率稳定且出行路径在时间和空间上分布较为均匀.结合算法有效性及路径漏检情况分析, 提出的基于车牌识别数据的机动车出行轨迹提取算法在测试区域的综合准确率在92%以上.

其中综合准确率的计算方法为

$ \rho = \sum\limits_{j = 1}^n {{\psi _j}{P_{dj}}} . $

(18)

式中：ρ为综合准确率, P_dj为测试区域内漏检节点数为j的路径所占总路径的百分比.

4 结语

提出了基于车牌识别数据的机动车出行轨迹提取算法, 并针对车辆出行轨迹中存在的漏检情况, 设计了基于KSP及GRA的出行轨迹补全决策算法, 该算法既充分利用了路网的静态拓扑属性, 同时科学地将路段行程时间、历史轨迹信息等动态信息属性融入到算法模型中, 克服了单纯依靠静态拓扑信息无法精确匹配实际轨迹的缺陷.最后, 以贵阳市南明区某区域为算例进行了算法的计算演示, 验证了该算法的有效性.