CN114298184A - 一种基于机器学习和优化算法的新能源公交工况构建方法 - Google Patents

一种基于机器学习和优化算法的新能源公交工况构建方法 Download PDF

Info

Publication number
CN114298184A
CN114298184A CN202111560908.4A CN202111560908A CN114298184A CN 114298184 A CN114298184 A CN 114298184A CN 202111560908 A CN202111560908 A CN 202111560908A CN 114298184 A CN114298184 A CN 114298184A
Authority
CN
China
Prior art keywords
matrix
working condition
clustering
characteristic parameters
bus
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202111560908.4A
Other languages
English (en)
Inventor
何洪文
黄汝臣
赵旭阳
孟祥飞
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Institute of Technology BIT
Original Assignee
Beijing Institute of Technology BIT
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Institute of Technology BIT filed Critical Beijing Institute of Technology BIT
Priority to CN202111560908.4A priority Critical patent/CN114298184A/zh
Publication of CN114298184A publication Critical patent/CN114298184A/zh
Pending legal-status Critical Current

Links

Images

Landscapes

  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明提供了一种基于机器学习和优化算法的新能源公交工况构建方法，其采用了多种方法联合，对采集到的数据进行预处理，降低了外界因素对采集的数据质量的负面影响。方法中所采用的主成分分析法对运动学片段的特征参数降维，在尽可能多地保留原始行驶数据信息的基础上，降低了计算的复杂度，消除了各特征参数之间的相关性，保证了分析结果的可靠性。通过改进的粒子群优化算法计算得到k‑means聚类算法的初始聚类中心，降低了k‑means算法对初始聚类中心的敏感性，提高了聚类结果的准确性。该方法最终构建的公交工况与原始行驶数据的统计特征差异很小，说明构建的公交工况能够很好地反映实际新能源公交线路的真实交通状况和驾驶特征。

Description

一种基于机器学习和优化算法的新能源公交工况构建方法
技术领域
本发明属于新能源汽车工况构建技术领域，尤其涉及一种基于机器学习和优化算法实现并适用于新能源公交车辆的工况构建方法。
背景技术
现阶段，针对新能源汽车的行驶特征分析与测试中均需要借助相应的标准工况数据，目前针对国内新能源公交汽车所主要应用的工况数据为“中国典型城市公交工况(Chinese Typical City Bus Driving Cycle,CCBC)”。然而，由于我国城市化水平的不断提高以及城市交通路网体系的逐渐完善，城市道路交通场景的多样性特征越来越明显，不同城市乃至不同线路的交通环境和驾驶特征均不尽相同，导致公交车的实际运行结果与以CCBC工况为基准认证的结果存在愈加明显的差异，因此有必要针对我国新能源公交的实际情况，构建能够反映本地区特定新能源公交路线交通状况和行驶特征的公交工况，为新能源公交车的燃油经济性和排放性提供更为精确的评价基准。
现有的工况构建方法主要分为两类：基于聚类算法的短行程构建方法和基于马尔科夫的构建方法，譬如现有技术：①.“基于K-均值聚类算法的行驶工况构建方法”(秦大同、詹森、漆正刚、陈淑江；《吉林大学学报(工学版)》，2016,46(02):383-389)；②.“基于LLEKM和马尔科夫链的城市轻型车行驶工况构建”(张惠玲、孔德学、余涛、敖谷昌、邵毅明；《长安大学学报(自然科学版)》，2021,41(05):118-126)。其中，基于聚类算法的短行程构建方法原理简单、易于实现，因此应用更加广泛，但由于聚类结果受初始聚类中心的影响较大，而通常情况下的聚类中心都是随机选取的，因而对所构建的公交工况精确程度易产生负面影响，使其不能很好地反映真实的交通场景和行驶特征；基于马尔科夫的构建方法实现难度较大，所构建工况的准确度在很大程度上依赖于所建立的马尔科夫状态转移矩阵的精确性，同时，马尔科夫状态转移矩阵的精确程度与计算成本之间存在矛盾，这也限制了基于马尔科夫的构建方法的广泛应用。此外，由于城市公交线路的各个站点都是固定的，每两个相邻公交站点之间的区间距离也都是一定的，这也对基于马尔科夫的构建方法造成了巨大的挑战。
发明内容
针对上述本领域中存在的技术问题，本发明提供了一种基于机器学习和优化算法的新能源公交工况构建方法，具体包括以下步骤：
步骤1、针对特定公交路线上行驶的新能源公交车辆，利用工况数据采集设备获取其原始行驶数据，并对所述原始行驶数据执行必要的预处理；
步骤2、将所述路线上每两个相邻公交站点之间的车速数据片段视为一个运动学片段，根据每个所述运动学片段的起点与终点时刻将预处理后的所述数据划分为多个运动学片段；
步骤3、从各运动学片段中提取与所述路线上的交通状况以及驾驶行为相关的特征参数，由提取的全部运动学片段的特征参数组成样本观测矩阵，并对所述样本观测矩阵执行标准化处理；
步骤4、利用主成分分析法对标准化处理后的所述样本观测矩阵执行降维，得到消除了所述特征参数之间相关性的若干主成分，用于描述不同运动学片段对应的交通状况以及驾驶行为；
步骤5、利用k-means聚类算法对所述主成分进行聚类，将各运动学片段分别划分为对应于不同交通场景的多个工况类型，完成该特定公交路线的完整工况构建；
步骤6、利用该特定公交路线上实际行驶的新能源公交车辆原始行驶数据，对所构建的完整工况进行对比验证，并实现定期更新。
进一步地，步骤1中所执行的预处理具体包括但不限于：
a.剔除发车等待、始发站站内行驶、终点站站内行驶、往返公交公司以及设备故障、采集中断等各种特定、特殊情况下的无效数据；
b.对因数据采集设备信号丢失造成的采样数据不连续，进而导致车速值短暂空缺的情况，利用插值法填充丢失的车速数据；
c.对因怠速时车辆抖动导致的车速小于1km/h的数据进行怠速处理；
d.根据采集时间，对车速数据进行时间间隔为1秒的插值；
e.对于短时间内车速变化较为剧烈，即加速度过大或过小等异常情况，利用滑动平均滤波法对车速数据进行降噪。
进一步地，步骤2划分运动学片段的方法具体包括：将每两个相邻的公交站点之间的车速片段视为一个运动学片段，第一个运动学片段的起点为在始发站的发车时刻，最后一个运动学片段的终点为到终点站的停车时刻，除此之外，其余所有运动学片段的起点和终点均分别为每两个相邻的公交站点从进站停车开始到离站停车结束的中间时刻。
进一步地，步骤3所提取的特征参数具体包括：速度特征参数、加速度特征参数、急动度特征参数、时间特征参数和距离特征参数5类；其中，速度特征参数包括：最高车速、平均车速、车速标准差；加速度特征参数包括：最大正加速度、平均正加速度、正加速度标准差、最小负加速度、平均负加速度、负加速度标准差、加速度标准差；急动度特征参数包括：最大正急动度、平均正急动度、正急动度标准差、最小负急动度、平均负急动度、负急动度标准差、急动度标准差；时间特征参数包括：怠速时间比例、加速时间比例、减速时间比例、运动学片段时长；距离特征参数为每个运动学片段相应的距离；
对于特定公交路线上m个运动学片段上分别包含的n个特征参数，构建以下形式的样本观测矩阵Am×n：
Figure BDA0003420599270000031
其中，aij(i＝1,2,…,m；j＝1,2,…,n)表示第i个片段样本的第j个特征参数；
为消除由于各个特征参数之间的量纲不统一带来的负面影响，对Am×n执行标准化处理后得到以下标准化矩阵Ym×n：
Figure BDA0003420599270000032
Figure BDA0003420599270000033
其中，yij(i＝1,2,…,m；j＝1,2,…,n)表示第i个样本经过标准化之后的第j个特征参数。
进一步地，步骤4中对样本观测矩阵执行降维具体包括执行以下过程：
4.1)计算所述标准化矩阵Ym×n的Pearson相关系数矩阵Rn×n：
Figure BDA0003420599270000034
Figure BDA0003420599270000035
其中，rij表示Pearson相关系数矩阵R的第i行第j列的元素(i＝1,2,…,n；j＝1,2,…,n)，
Figure BDA0003420599270000036
表示标准矩阵Y第i列的平均值，
Figure BDA0003420599270000037
表示标准矩阵Ym×n第j列的平均值，yki表示标准矩阵Ym×n的第k行第i列的元素，ykj表示标准矩阵Ym×n的第k行第j列的元素；
4.2)计算Pearson相关系数矩阵R的非负特征值，记为λi(i＝1,2,…,n)，将各个特征值按照从大到小的顺序排列，得到λ1>λ2>…>λn≥0；
4.3)计算各个主成分对应的贡献率μi：
Figure BDA0003420599270000041
4.4)将贡献率μi逐个累加，若累计贡献率达到预定阈值，则停止累加，此时确定出待提取的主成分的个数，即待提取的主成分的个数等于参与累加的贡献率的个数，记为w；
4.5)计算各特征值λi对应的特征向量li，形式为行向量，各特征向量组成的矩阵转置后即主成分载荷矩阵
Figure BDA0003420599270000042
Figure BDA0003420599270000043
其中，lij表示载荷矩阵的各元素(i＝1,2,…,n；j＝1,2,…,n)；
4.6)根据标准化矩阵Y和特征向量li计算各主成分的得分向量hi：
Figure BDA0003420599270000044
其中，hi为列向量(i＝1,2,…,n)yi为标准化矩阵Ym×n的列向量(i＝1,2,…,n)；
将各主成分的得分向量合并后转置，得到主成分得分矩阵
Figure BDA0003420599270000045
Figure BDA0003420599270000046
取得分矩阵
Figure BDA0003420599270000047
中的前w列用于后续步骤5的聚类。
进一步地，步骤5中对于基于k-means聚类算法的聚类，首先利用改进的粒子群优化算法计算初始聚类中心，包括：
①利用最值距离乘积法取出粒子群优化的初始聚类中心，所有的初始聚类中心组成一个粒子，聚类中心的个数等于聚类簇数，重复α次，产生α个粒子，完成粒子群的初始化；
②对每一个粒子，根据初始聚类中心，按照欧式距离最小化原则进行一次聚类，并计算每个粒子的适应度值，寻找每个粒子的最优位置和粒子群的全局最优位置；
③对于粒子群优化算法的惯性因子和学习因子基于以下公式进行计算：
Figure BDA0003420599270000051
其中，ωi表示第i个粒子的惯性因子(i＝1,2,…,α)，fiti表示第i个粒子的适应度，fitbest表示粒子群的全局最优值，t表示当前的迭代次数，tmax表示最大迭代次数；
④更新各个粒子的速度与位置：
Figure BDA0003420599270000052
其中，vi表示第i个粒子的速度(i＝1,2,…,α)，xi表示第i个粒子的位置，上标t表示当前的迭代次数，ρ1和ρ2均为[0,1]之间的随机数，pi表示第i个粒子的最优位置，pg表示粒子群的最优位置；
⑤将更新之后的粒子用做聚类中心，重复步骤②至步骤④；
⑥判断当前是否满足条件即迭代次数是否达到设定的阈值tmax，若满足条件，则结束粒子群优化，若不满足，则转向步骤③。
进一步地，步骤5中基于k-means聚类算法的聚类具体包括依次执行：
5.1)建立样本类别矩阵Qm×k，并初始化其为全零矩阵；
5.2)对所有样本数据进行聚类，并更新样本类别矩阵Qm×k，更新公式为：
Figure BDA0003420599270000053
其中，qij表示样本类别矩阵中的元素(i＝1,2,…,m；j＝1,2,…,k)，xi表示第i个样本，Cj表示第j类簇；
5.3)判断当前样本类别矩阵Qm×k是否满足不再发生变化或聚类次数达到预设的上限值，若满足条件，则结束k-means聚类，输出聚类结果，若不满足，则转向步骤5.2。
进一步地，步骤5中在完成聚类和不同工况类型的划分以后，计算每一个工况片段的特征参数，判断每一个工况片段所具体对应的交通场景，包括：计算每两个相邻站点之间的公交区间对应不同交通场景的概率，取概率最大的交通场景为该公交区间所属；若在多于一种交通场景下的概率相等，则分别计算相应交通场景下的该公交区间的所有运动学片段与相应聚类中心的平均距离，取平均距离最小的交通场景为该公交区间所属；
针对该特定公交路线上的全部公交区间完成交通场景的确定，即实现了对该特定公交路线的工况构建。
进一步地，步骤6中对对所构建的完整工况进行对比验证，可具体采用以下方式：
a.分别计算构建的新能源公交工况和原始行驶数据的特征参数，计算差异率；
b.分别绘制构建的新能源公交工况和原始行驶数据的速度-加速度联合概率分布图，以及二者之间的速度-加速度联合概率误差分布图。
上述本发明所提供的基于机器学习和优化算法的新能源公交工况构建方法，相对现有技术至少具有以下有益效果：
1、其采用了多种方法联合，对采集到的数据进行预处理，降低了外界因素对采集的数据质量的负面影响。
2、方法中所采用的主成分分析法对运动学片段的特征参数降维，在尽可能多地保留原始行驶数据信息的基础上，降低了计算的复杂度，消除了各特征参数之间的相关性，保证了分析结果的可靠性。
3、通过改进的粒子群优化算法计算得到k-means聚类算法的初始聚类中心，降低了k-means算法对初始聚类中心的敏感性，提高了聚类结果的准确性。
4、该方法最终构建的公交工况与原始行驶数据的统计特征差异很小，说明构建的公交工况能够很好地反映实际新能源公交线路的真实交通状况和驾驶特征。
附图说明
图1为本发明所提供方法的整体流程图；
图2为Pearson相关系数矩阵热力图；
图3为累计贡献率的pareto图；
图4为改进粒子群优化算法的实施流程图；
图5为用前3个主成分展示的最终聚类结果；
图6为聚类得到的3种不同工况片段库；
图7为基于本发明所构建的郑州市727路新能源公交工况；
图8为基于本发明所构建的工况与原始行驶数据的验证对比结果。
具体实施方式
下面将结合附图对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
本发明提供的一种基于机器学习和优化算法的新能源公交工况构建方法，如图1所示，具体包括以下步骤：
步骤1、针对特定公交路线上行驶的新能源公交车辆，利用工况数据采集设备获取其原始行驶数据，并对所述原始行驶数据执行必要的预处理；
步骤2、将所述路线上每两个相邻公交站点之间的车速数据片段视为一个运动学片段，根据每个所述运动学片段的起点与终点时刻将预处理后的所述数据划分为多个运动学片段；
步骤3、从各运动学片段中提取与所述路线上的交通状况以及驾驶行为相关的特征参数，由提取的全部运动学片段的特征参数组成样本观测矩阵，并对所述样本观测矩阵执行标准化处理；
步骤4、利用主成分分析法对标准化处理后的所述样本观测矩阵执行降维，得到消除了所述特征参数之间相关性的若干主成分，用于描述不同运动学片段对应的交通状况以及驾驶行为；
步骤5、利用改进的粒子群优化算法和k-means聚类算法对所述主成分进行聚类，将各运动学片段分别划分为对应于不同交通场景的多个工况类型，完成该特定公交路线的完整工况构建；
步骤6、利用该特定公交路线上实际行驶的新能源公交车辆原始行驶数据，对所构建的完整工况进行对比验证，并实现定期更新。
在本发明的一个优选实施方式中，步骤1中所执行的预处理具体包括但不限于：
a.剔除发车等待、始发站站内行驶、终点站站内行驶、往返公交公司以及设备故障、采集中断等各种特定、特殊情况下的无效数据；
b.对因数据采集设备信号丢失造成的采样数据不连续，进而导致车速值短暂空缺的情况，利用插值法填充丢失的车速数据；
c.对因怠速时车辆抖动导致的车速小于1km/h的数据进行怠速处理；
d.根据采集时间，对车速数据进行时间间隔为1秒的插值；
e.对于短时间内车速变化较为剧烈，即加速度过大或过小等异常情况，利用滑动平均滤波法对车速数据进行降噪。
在本发明的一个优选实施方式中，针对郑州新能源公交路线727路，提取原始数据并经过预处理后的一共有28条从始发站到终点站的完整车速-时间曲线，总时长共约20小时，总里程共约410千米。通过将每两个相邻的公交站点之间的车速片段视为一个运动学片段，第一个运动学片段的起点为在始发站的发车时刻，最后一个运动学片段的终点为到终点站的停车时刻，除此之外，其余所有运动学片段的起点和终点均分别为每两个相邻的公交站点从进站停车开始到离站停车结束的中间时刻，共得到672个运动学片段。
在本发明的一个优选实施方式中，步骤2划分运动学片段的方法具体包括：将每两个相邻的公交站点之间的车速片段视为一个运动学片段，第一个运动学片段的起点为在始发站的发车时刻，最后一个运动学片段的终点为到终点站的停车时刻，除此之外，其余所有运动学片段的起点和终点均分别为每两个相邻的公交站点从进站停车开始到离站停车结束的中间时刻。
在本发明的一个优选实施方式中，步骤3所提取的特征参数具体包括：速度特征参数、加速度特征参数、急动度特征参数、时间特征参数和距离特征参数5类；其中，速度特征参数包括：最高车速、平均车速、车速标准差；加速度特征参数包括：最大正加速度、平均正加速度、正加速度标准差、最小负加速度、平均负加速度、负加速度标准差、加速度标准差；急动度特征参数包括：最大正急动度、平均正急动度、正急动度标准差、最小负急动度、平均负急动度、负急动度标准差、急动度标准差；时间特征参数包括：怠速时间比例、加速时间比例、减速时间比例、运动学片段时长；距离特征参数为每个运动学片段相应的距离。具体如表1所示：
表1运动学片段特征参数
编号 特征参数 单位 编号 特征参数 单位
p01 最高车速 km/h p12 平均正急动度 m/s<sup>3</sup>
p02 平均车速 km/h p13 正急动度标准差 m/s<sup>3</sup>
p03 车速标准差 km/h p14 最小负急动度 m/s<sup>3</sup>
p04 最大正加速度 m/s<sup>2</sup> p15 平均负急动度 m/s<sup>3</sup>
p05 平均正加速度 m/s<sup>2</sup> p16 负急动度标准差 m/s<sup>3</sup>
p06 正加速度标准差 m/s<sup>2</sup> p17 急动度标准差 m/s<sup>3</sup>
p07 最小负加速度 m/s<sup>2</sup> p18 怠速时间比例 ％
p08 平均负加速度 m/s<sup>2</sup> p19 加速时间比例 ％
p09 负加速度标准差 m/s<sup>2</sup> p20 减速时间比例 ％
p10 加速度标准差 m/s<sup>2</sup> p21 运动学片段时长 s
p11 最大正急动度 m/s<sup>3</sup> p22 运动学片段距离 m
对于特定公交路线上m个运动学片段上分别包含的n个特征参数，构建以下形式的样本观测矩阵Am×n：
Figure BDA0003420599270000081
其中，aij(i＝1,2,…,m；j＝1,2,…,n)表示第i个片段样本的第j个特征参数；
为消除由于各个特征参数之间的量纲不统一带来的负面影响，对Am×n执行标准化处理后得到以下标准化矩阵Ym×n：
Figure BDA0003420599270000082
Figure BDA0003420599270000083
其中，yij(i＝1,2,…,m；j＝1,2,…,n)表示第i个样本经过标准化之后的第j个特征参数。
在本发明的一个优选实施方式中，步骤4中对样本观测矩阵执行降维具体包括执行以下过程：
4.1)计算所述标准化矩阵Ym×n的Pearson相关系数矩阵Rn×n：
Figure BDA0003420599270000091
Figure BDA0003420599270000092
其中，rij表示Pearson相关系数矩阵R的第i行第j列的元素(i＝1,2,…,n；j＝1,2,…,n)，
Figure BDA0003420599270000093
表示标准矩阵Y第i列的平均值，
Figure BDA0003420599270000094
表示标准矩阵Ym×n第j列的平均值，yki表示标准矩阵Ym×n的第k行第i列的元素，ykj表示标准矩阵Ym×n的第k行第j列的元素，计算得到的Pearson相关系数矩阵热力图如图2所示；
4.2)计算Pearson相关系数矩阵R的非负特征值，记为λi(i＝1,2,…,n)，将各个特征值按照从大到小的顺序排列，得到λ1>λ2>…>λn≥0；
4.3)计算各个主成分对应的贡献率μi：
Figure BDA0003420599270000095
具体计算得到的各主成分的特征值和贡献率如表2所示：
表2各主成分的特征值和贡献率
主成分 特征值 贡献率/％ 主成分 特征值 贡献率/％
1 7.2380 32.8999 12 0.1402 0.6372
2 4.7680 21.6729 13 0.1285 0.5840
3 2.7256 12.3890 14 0.1202 0.5463
4 2.4080 10.9453 15 0.0897 0.4077
5 1.4028 6.3764 16 0.0749 0.3404
6 0.7423 3.3741 17 0.0653 0.2968
7 0.6786 3.0845 18 0.0461 0.2093
8 0.5836 2.6529 19 0.0169 0.0768
9 0.3860 1.7547 20 0.0103 0.0468
10 0.2100 0.9544 21 0.0027 0.0121
11 0.1616 0.7344 22 0.0009 0.0041
4.4)将贡献率μi逐个累加，若累计贡献率达到预定阈值，在此实施例中选择为80％，则停止累加，此时确定出待提取的主成分的个数，即待提取的主成分的个数w等于参与累加的贡献率的个数，累计贡献率的pareto图如图3所示，前5个主成分的累计贡献率超过了80％，为84.28％，则可认为只提取前5个主成分便可以包含80％以上的原始行驶数据信息，因此前5个主成分可以替代22个特征参数，消除了特征参数之间的相关性；4.5)计算各特征值λi对应的特征向量li，形式为行向量，各特征向量组成的矩阵转置后即主成分载荷矩阵
Figure BDA0003420599270000096
Figure BDA0003420599270000101
其中，lij表示载荷矩阵的各元素(i＝1,2,…,n；j＝1,2,…,n)；
4.6)根据标准化矩阵Y和特征向量li计算各主成分的得分向量hi：
Figure BDA0003420599270000102
其中，hi为列向量(i＝1,2,…,n)yi为标准化矩阵Ym×n的列向量(i＝1,2,…,n)；
将各主成分的得分向量合并后转置，得到主成分得分矩阵
Figure BDA0003420599270000103
Figure BDA0003420599270000104
取得分矩阵
Figure BDA0003420599270000105
中的前w列用于后续步骤5的聚类。
在本发明的一个优选实施方式中，步骤5中对于基于k-means聚类算法的聚类，首先利用改进的粒子群优化算法计算初始聚类中心，本领域技术人员也应当知晓，即使不额外采用该改进的粒子群算法，本发明仍然可以解决相应的技术问题并实现相应的有益效果。如图4所示，包括：
①利用最值距离乘积法取出粒子群优化的初始聚类中心，所有的初始聚类中心组成一个粒子，聚类中心的个数等于聚类簇数，重复α次，产生α个粒子，完成粒子群的初始化；
②对每一个粒子，根据初始聚类中心，按照欧式距离最小化原则进行一次聚类，并计算每个粒子的适应度值，寻找每个粒子的最优位置和粒子群的全局最优位置；
欧式距离最小化原则的计算公式为：
Figure BDA0003420599270000106
Figure BDA0003420599270000107
式中，Err表示距离的平方误差，k表示聚类的簇的个数，Ci表示聚类的簇(i＝1,2,…,k)，x表示样本，δi表示各簇的中心点，||x-δi||2表示各样本与中心点的欧氏距离；
各粒子的适应度函数为：
Figure BDA0003420599270000111
式中，fitε表示各粒子的适应度(ε＝1,2,…,α)，
Figure BDA0003420599270000112
表示第i类的样本数(i＝1,2,…,k)，xij表示第i类的第j个样本；
不同于常规粒子群优化算法中的惯性因子线性变化、学习因子始终为常数，本发明对粒子群优化过程中的惯性因子ω、自身学习因子c1和全局学习因子c2根据各粒子的适应度值进行自适应调整：不同粒子的惯性因子ω因适应度值的不同而不同，粒子的适应度值越小，则对应的ω越大，表明粒子此时所处的位置较差，需要提高全局搜索能力，加快收敛速度，反之，若粒子的适应度值越大，则对应的ω越小，表明粒子此时所处的位置较好，需要提高局部搜索能力；自身学习因子c1在迭代初期取值较大，此时粒子受个体最优的影响更大，可以避免种群在迭代前期陷入局部最优，全局学习因子c2在迭代后期取值较大，此时粒子受全局最优的影响更大，保证了种群在迭代后期向全局最优位置收敛，进而加快了收敛速度。因此，对于粒子群优化算法的惯性因子和学习因子基于以下公式进行计算：
Figure BDA0003420599270000113
其中，ωi表示第i个粒子的惯性因子(i＝1,2,…,α)，fiti表示第i个粒子的适应度，fitbest表示粒子群的全局最优值，t表示当前的迭代次数，tmax表示最大迭代次数；
③更新各个粒子的速度与位置：
Figure BDA0003420599270000114
其中，vi表示第i个粒子的速度(i＝1,2,…,α)，xi表示第i个粒子的位置，上标t表示当前的迭代次数，ρ1和ρ2均为[0,1]之间的随机数，pi表示第i个粒子的最优位置，pg表示粒子群的最优位置；
④将更新之后的粒子用做聚类中心，重复步骤5.2至步骤5.4；
⑤判断当前是否满足条件即迭代次数是否达到设定的阈值tmax，若满足条件，则结束粒子群优化，若不满足，则转向步骤5.3。
在本发明的一个优选实施方式中，所述最值距离乘积法具体包括：
①从样本集合
Figure BDA0003420599270000115
中随机选取一个样本元素作为第一个初始聚类中心δ1，将δ1加入到粒子群优化的初始聚类中心集合Φ中，并将该样本元素从集合
Figure BDA0003420599270000116
中删除；
②计算更新后的集合
Figure BDA0003420599270000117
中所有元素到δ1的距离，选取距离最大的元素作为第二个初始聚类中心δ2；
③将δ2加入到集合Φ中，并将对应的样本元素从集合
Figure BDA0003420599270000118
中删除；
④分别计算更新后的集合
Figure BDA0003420599270000122
中每个元素到集合Φ中各个元素的距离，存储为数组S，数组S的个数与集合
Figure BDA0003420599270000123
中的元素个数相等；
⑤计算集合
Figure BDA0003420599270000124
中的每个元素对应的数组S中的最大值与最小值的乘积，取最大的乘积值对应的元素；
⑥判断当前是否满足条件：选取的初始聚类中心的个数等于待选取的聚类中心个数k，若满足条件，则输出粒子群优化的初始聚类中心集合Φ，若不满足，则转向步骤③。
进一步地，步骤5中基于k-means聚类算法的聚类具体包括依次执行：
5.1)建立样本类别矩阵Qm×k，并初始化其为全零矩阵；
5.2)对所有样本数据进行聚类，并更新样本类别矩阵Qm×k，更新公式为：
Figure BDA0003420599270000121
其中，qij表示样本类别矩阵中的元素(i＝1,2,…,m；j＝1,2,…,k)，xi表示第i个样本，Cj表示第j类簇；
5.3)判断当前样本类别矩阵Qm×k是否满足不再发生变化或聚类次数达到预设的上限值，若满足条件，则结束k-means聚类，输出聚类结果，若不满足，则转向步骤5.2。为方便可视化，用前3个主成分展示的最终聚类结果如图5所示。
在本发明的一个优选实施方式中，步骤5中在完成聚类和不同工况类型的划分以后，畅行工况、正常工况和拥堵工况这三种不同的交通场景，计算每一个工况片段的特征参数，选取部分特征参数如表3所示：工况片段库一的最高车速和平均车速最高、怠速时间比例最短，为畅行工况；工况片段库三的最高车速和平均车速最低、怠速时间比例最长，为拥堵工况；工况片段库二的最高车速、平均车速最低和怠速时间比例均介于工况片段库一和工况片段库三之间，为正常工况；畅行工况片段库、正常工况片段库、拥堵工况片段库分别如图6a、图6b、图6c所示；
表3工况片段库的部分特征参数
特征参数 单位 工况片段库一 工况片段库二 工况片段库三
最高车速 km/h 60.30 57.36 55.20
平均车速 km/h 29.19 19.03 15.18
怠速时间比例 ％ 10.66 24.31 28.07
最大正加速度 m/s<sup>2</sup> 1.74 1.74 1.72
平均正加速度 m/s<sup>2</sup> 0.42 0.65 0.48
加速时间比例 ％ 47.72 40.04 37.06
最小负加速度 m/s<sup>2</sup> -2.26 -2.37 -1.90
平均负加速度 m/s<sup>2</sup> -0.47 -0.70 -0.49
减速时间比例 ％ 42.74 37.29 36.18
判断每一个工况片段所具体对应的交通场景，包括：计算每两个相邻站点之间的公交区间对应不同交通场景的概率，取概率最大的交通场景为该公交区间所属；若在多于一种交通场景下的概率相等，则分别计算相应交通场景下的该公交区间的所有运动学片段与相应聚类中心的平均距离，取平均距离最小的交通场景为该公交区间所属；
针对该特定公交路线上的全部公交区间完成交通场景的确定，即实现了对该特定公交路线的工况构建，如图7所示。
在本发明的一个优选实施方式中，步骤6中对对所构建的完整工况进行对比验证，可具体采用以下方式：
a.分别计算构建的新能源公交工况和原始行驶数据的特征参数，计算差异率，如表4所示，对于本实施例，各项的差异率均在3％以内；
表4部分特征参数及其差异率
特征参数 单位 原始行驶数据 构建的公交工况 差异率/％
平均车速 km/h 20.4022 19.9691 2.12
车速标准差 km/h 17.4011 17.1755 1.30
平均正加速度 m/s<sup>2</sup> 0.5048 0.4913 2.67
正加速度标准差 m/s<sup>2</sup> 0.4014 0.3954 1.49
平均负加速度 m/s<sup>2</sup> -0.5401 -0.5466 1.20
负加速度标准差 m/s<sup>2</sup> 0.4531 0.4533 0.04
怠速时间比例 ％ 0.2261 0.2204 2.52
b.分别绘制原始行驶数据和构建的公交工况的速度-加速度联合概率分布图如图8a、图8b所示，以及二者之间的速度-加速度联合概率误差分布图如图8c所示，图8a和图8b的相似度很高，误差较大的地方发生在加速度较低处，且最大误差不超过0.01，以上关于差异率和速度-加速度联合概率分布图的分析说明了构建的工况能很好地代表原始行驶数据，并能很好地反映郑州市727路新能源公交路线的交通状况和行驶特征。
应理解，本发明实施例中各步骤的序号的大小并不意味着执行顺序的先后，各过程的执行顺序应以其功能和内在逻辑确定，而不应对本发明实施例的实施过程构成任何限定。
尽管已经示出和描述了本发明的实施例，对于本领域的普通技术人员而言，可以理解在不脱离本发明的原理和精神的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由所附权利要求及其等同物限定。

Claims (9)

1.一种基于机器学习和优化算法的新能源公交工况构建方法，其特征在于：具体包括以下步骤：
步骤1、针对特定公交路线上行驶的新能源公交车辆，利用工况数据采集设备获取其原始行驶数据，并对所述原始行驶数据执行必要的预处理；
步骤2、将所述路线上每两个相邻公交站点之间的车速数据片段视为一个运动学片段，根据每个所述运动学片段的起点与终点时刻将预处理后的所述数据划分为多个运动学片段；
步骤3、从各运动学片段中提取与所述路线上的交通状况以及驾驶行为相关的特征参数，由提取的全部运动学片段的特征参数组成样本观测矩阵，并对所述样本观测矩阵执行标准化处理；
步骤4、利用主成分分析法对标准化处理后的所述样本观测矩阵执行降维，得到消除了所述特征参数之间相关性的若干主成分，用于描述不同运动学片段对应的交通状况以及驾驶行为；
步骤5、利用k-means聚类算法对所述主成分进行聚类，将各运动学片段分别划分为对应于不同交通场景的多个工况类型，完成该特定公交路线的完整工况构建；
步骤6、利用该特定公交路线上实际行驶的新能源公交车辆原始行驶数据，对所构建的完整工况进行对比验证，并实现定期更新。
2.如权利要求1所述的方法，其特征在于：步骤1中所执行的预处理具体包括：
a.剔除发车等待、始发站站内行驶、终点站站内行驶、往返公交公司以及设备故障、采集的各种无效数据；
b.对因数据采集设备信号丢失造成的采样数据不连续，进而导致车速值短暂空缺的情况，利用插值法填充丢失的车速数据；
c.对因怠速时车辆抖动导致的车速小于1km/h的数据进行怠速处理；
d.根据采集时间，对车速数据进行时间间隔为1秒的插值；
e.对于短时间内车速异常情况，利用滑动平均滤波法对车速数据进行降噪。
3.如权利要求1所述的方法，其特征在于：步骤2划分运动学片段的方法具体包括：将每两个相邻的公交站点之间的车速片段视为一个运动学片段，第一个运动学片段的起点为在始发站的发车时刻，最后一个运动学片段的终点为到终点站的停车时刻，除此之外，其余所有运动学片段的起点和终点均分别为每两个相邻的公交站点从进站停车开始到离站停车结束的中间时刻。
4.如权利要求1所述的方法，其特征在于：步骤3所提取的特征参数具体包括：速度特征参数、加速度特征参数、急动度特征参数、时间特征参数和距离特征参数5类；其中，速度特征参数包括：最高车速、平均车速、车速标准差；加速度特征参数包括：最大正加速度、平均正加速度、正加速度标准差、最小负加速度、平均负加速度、负加速度标准差、加速度标准差；急动度特征参数包括：最大正急动度、平均正急动度、正急动度标准差、最小负急动度、平均负急动度、负急动度标准差、急动度标准差；时间特征参数包括：怠速时间比例、加速时间比例、减速时间比例、运动学片段时长；距离特征参数为每个运动学片段相应的距离；
对于特定公交路线上m个运动学片段上分别包含的n个特征参数，构建以下形式的样本观测矩阵Am×n：
Figure FDA0003420599260000021
其中，aij(i＝1,2,…,m；j＝1,2,…,n)表示第i个片段样本的第j个特征参数；
为消除由于各个特征参数之间的量纲不统一带来的负面影响，对Am×n执行标准化处理后得到以下标准化矩阵Ym×n：
Figure FDA0003420599260000022
Figure FDA0003420599260000023
其中，yij(i＝1,2,…,m；j＝1,2,…,n)表示第i个样本经过标准化之后的第j个特征参数。
5.如权利要求4所述的方法，其特征在于：步骤4中对样本观测矩阵执行降维具体包括执行以下过程：
4.1)计算所述标准化矩阵Ym×n的Pearson相关系数矩阵Rn×n：
Figure FDA0003420599260000024
Figure FDA0003420599260000031
其中，rij表示Pearson相关系数矩阵R的第i行第j列的元素(i＝1,2,…,n；j＝1,2,…,n)，
Figure FDA0003420599260000032
表示标准矩阵Y第i列的平均值，
Figure FDA0003420599260000033
表示标准矩阵Ym×n第j列的平均值，yki表示标准矩阵Ym×n的第k行第i列的元素，ykj表示标准矩阵Ym×n的第k行第j列的元素；
4.2)计算Pearson相关系数矩阵R的非负特征值，记为λi(i＝1,2,…,n)，将各个特征值按照从大到小的顺序排列，得到λ1>λ2>…>λn≥0；
4.3)计算各个主成分对应的贡献率μi：
Figure FDA0003420599260000034
4.4)将贡献率μi逐个累加，若累计贡献率达到预定阈值，则停止累加，此时确定出待提取的主成分的个数，即待提取的主成分的个数等于参与累加的贡献率的个数，记为w；
4.5)计算各特征值λi对应的特征向量li，形式为行向量，各特征向量组成的矩阵转置后即主成分载荷矩阵
Figure FDA0003420599260000035
Figure FDA0003420599260000036
其中，lij表示载荷矩阵的各元素(i＝1,2,…,n；j＝1,2,…,n)；
4.6)根据标准化矩阵Y和特征向量li计算各主成分的得分向量hi：
Figure FDA0003420599260000037
其中，hi为列向量(i＝1,2,…,n)yi为标准化矩阵Ym×n的列向量(i＝1,2,…,n)；
将各主成分的得分向量合并后转置，得到主成分得分矩阵
Figure FDA0003420599260000038
Figure FDA0003420599260000039
取得分矩阵
Figure FDA00034205992600000310
中的前w列用于后续步骤5的聚类。
6.如权利要求1所述的方法，其特征在于：步骤5中对于基于k-means聚类算法的聚类，首先利用改进的粒子群优化算法计算初始聚类中心，包括：
①利用最值距离乘积法取出粒子群优化的初始聚类中心，所有的初始聚类中心组成一个粒子，聚类中心的个数等于聚类簇数，重复α次，产生α个粒子，完成粒子群的初始化；
②对每一个粒子，根据初始聚类中心，按照欧式距离最小化原则进行一次聚类，并计算每个粒子的适应度值，寻找每个粒子的最优位置和粒子群的全局最优位置；
③对于粒子群优化算法的惯性因子和学习因子基于以下公式进行计算：
Figure FDA0003420599260000041
其中，ωi表示第i个粒子的惯性因子(i＝1,2,…,α)，fiti表示第i个粒子的适应度，fitbest表示粒子群的全局最优值，t表示当前的迭代次数，tmax表示最大迭代次数；
④更新各个粒子的速度与位置：
Figure FDA0003420599260000042
其中，vi表示第i个粒子的速度(i＝1,2,…,α)，xi表示第i个粒子的位置，上标t表示当前的迭代次数，ρ1和ρ2均为[0,1]之间的随机数，pi表示第i个粒子的最优位置，pg表示粒子群的最优位置；
⑤将更新之后的粒子用做聚类中心，重复步骤②至步骤④；
⑥判断当前是否满足条件即迭代次数是否达到设定的阈值tmax，若满足条件，则结束粒子群优化，若不满足，则转向步骤③。
7.如权利要求1所述的方法，其特征在于：步骤5中基于k-means聚类算法的聚类具体包括依次执行：
5.1)建立样本类别矩阵Qm×k，并初始化其为全零矩阵；
5.2)对所有样本数据进行聚类，并更新样本类别矩阵Qm×k，更新公式为：
Figure FDA0003420599260000043
其中，qij表示样本类别矩阵中的元素(i＝1,2,…,m；j＝1,2,…,k)，xi表示第i个样本，Cj表示第j类簇；
5.3)判断当前样本类别矩阵Qm×k是否满足不再发生变化或聚类次数达到预设的上限值，若满足条件，则结束k-means聚类，输出聚类结果，若不满足，则转向步骤5.2。
8.如权利要求1所述的方法，其特征在于：步骤5中在完成聚类和不同工况类型的划分以后，计算每一个工况片段的特征参数，判断每一个工况片段所具体对应的交通场景，包括：计算每两个相邻站点之间的公交区间对应不同交通场景的概率，取概率最大的交通场景为该公交区间所属；若在多于一种交通场景下的概率相等，则分别计算相应交通场景下的该公交区间的所有运动学片段与相应聚类中心的平均距离，取平均距离最小的交通场景为该公交区间所属；
针对该特定公交路线上的全部公交区间完成交通场景的确定，即实现了对该特定公交路线的工况构建。
9.如权利要求1所述的方法，其特征在于：步骤6中对对所构建的完整工况进行对比验证，具体采用以下方式：
a.分别计算构建的新能源公交工况和原始行驶数据的特征参数，计算差异率；
b.分别绘制构建的新能源公交工况和原始行驶数据的速度-加速度联合概率分布图，以及二者之间的速度-加速度联合概率误差分布图。
CN202111560908.4A 2021-12-20 2021-12-20 一种基于机器学习和优化算法的新能源公交工况构建方法 Pending CN114298184A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202111560908.4A CN114298184A (zh) 2021-12-20 2021-12-20 一种基于机器学习和优化算法的新能源公交工况构建方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202111560908.4A CN114298184A (zh) 2021-12-20 2021-12-20 一种基于机器学习和优化算法的新能源公交工况构建方法

Publications (1)

Publication Number Publication Date
CN114298184A true CN114298184A (zh) 2022-04-08

Family

ID=80968115

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202111560908.4A Pending CN114298184A (zh) 2021-12-20 2021-12-20 一种基于机器学习和优化算法的新能源公交工况构建方法

Country Status (1)

Country Link
CN (1) CN114298184A (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114662620A (zh) * 2022-05-24 2022-06-24 岚图汽车科技有限公司 一种市场用户的汽车耐久载荷数据处理方法和装置
CN116719831A (zh) * 2023-08-03 2023-09-08 四川中测仪器科技有限公司 一种面向健康监测的标准数据库建立与更新方法
CN116776229A (zh) * 2023-08-17 2023-09-19 深圳市城市交通规划设计研究中心股份有限公司 面向排放碳因子的汽车行驶典型工况划分方法
CN117669998A (zh) * 2024-02-01 2024-03-08 聊城大学 一种考虑乘客载荷变化的公交工况构建方法

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114662620A (zh) * 2022-05-24 2022-06-24 岚图汽车科技有限公司 一种市场用户的汽车耐久载荷数据处理方法和装置
CN116719831A (zh) * 2023-08-03 2023-09-08 四川中测仪器科技有限公司 一种面向健康监测的标准数据库建立与更新方法
CN116719831B (zh) * 2023-08-03 2023-10-27 四川中测仪器科技有限公司 一种面向健康监测的标准数据库建立与更新方法
CN116776229A (zh) * 2023-08-17 2023-09-19 深圳市城市交通规划设计研究中心股份有限公司 面向排放碳因子的汽车行驶典型工况划分方法
CN116776229B (zh) * 2023-08-17 2023-12-26 深圳市城市交通规划设计研究中心股份有限公司 面向碳排放因子的汽车行驶典型工况划分方法
CN117669998A (zh) * 2024-02-01 2024-03-08 聊城大学 一种考虑乘客载荷变化的公交工况构建方法

Similar Documents

Publication Publication Date Title
CN114298184A (zh) 一种基于机器学习和优化算法的新能源公交工况构建方法
CN111081016B (zh) 一种基于复杂网络理论的城市交通异常识别方法
Xie et al. A hybrid method combining Markov prediction and fuzzy classification for driving condition recognition
CN107146013A (zh) 一种基于灰色预测和支持向量机的分类型电动汽车需求时空分布动态预测方法
CN109635914B (zh) 基于混合智能遗传粒子群的优化极限学习机轨迹预测方法
CN109887279B (zh) 一种交通拥堵预测方法及系统
CN110634299B (zh) 基于多源轨迹数据的城市交通状态精细划分与识别方法
CN110979342B (zh) 一种用于车辆全局能量管理控制的工况信息获取方法
CN113436433B (zh) 一种高效的城市交通离群值检测方法
CN112785077A (zh) 基于时空数据的出行需求预测方法及系统
CN111815948A (zh) 基于工况特征的车辆行驶工况预测方法
CN113222385A (zh) 一种电动汽车行驶工况构建与评价方法
CN114492919A (zh) 一种基于交通流量的电动汽车充电负荷预测系统及方法
CN112884014A (zh) 一种基于路段拓扑结构分类的交通速度短时预测方法
CN114446064A (zh) 一种分析高速公路服务区流量的方法、装置、存储介质及终端
CN115169714A (zh) 城市地铁进出站客流量预测方法
CN111191824A (zh) 一种动力电池容量衰减预测方法及系统
CN113642768A (zh) 一种基于工况重构的车辆行驶能耗预测方法
CN112950926A (zh) 一种基于大数据和深度学习的城市主干道路速度预测方法
CN109858559B (zh) 基于交通流宏观基本图的自适应交通分析路网简化方法
CN110728459A (zh) 出行方式识别系统、方法、装置及模型训练方法、装置
CN112948965A (zh) 一种基于机器学习和统计验证的汽车行驶工况的构建方法
CN109741597B (zh) 一种基于改进深度森林的公交车路段运行时间预测方法
CN115100847B (zh) 面向低渗透率网联车轨迹数据的排队服务时间估计方法
Yugendar et al. Driving cycle estimation and validation for Ludhiana City, India

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination