基于K-均值聚類算法的汽車行駛工況構(gòu)建

王震，宋志堯

（200093 上海市 上海理工大學(xué)）

0 引言

汽車行駛工況代表車輛在實(shí)際行駛時(shí)速度與時(shí)間的對(duì)應(yīng)函數(shù)關(guān)系。行駛工況作為一項(xiàng)重要的參考指標(biāo)，影響著車輛在排放方面的測(cè)試標(biāo)準(zhǔn)以及汽車各項(xiàng)性能指標(biāo)[1]。近幾年，我國(guó)汽車保有量快速增長(zhǎng)，道路交通狀況改變了很多，昔日以NEDC 工況為基準(zhǔn)所優(yōu)化標(biāo)定的汽車，實(shí)際油耗產(chǎn)生了越來(lái)越大的偏差[2]。由于地域不同，各個(gè)城市間的地理環(huán)境差異性較大，造成各個(gè)城市的交通發(fā)展情況不一[3]，需要深入研究我國(guó)的實(shí)際道路的測(cè)試工況作為車輛測(cè)試與評(píng)價(jià)的重要依據(jù)[4]。

本文以某城市實(shí)際運(yùn)行載荷數(shù)據(jù)，對(duì)原始載荷數(shù)據(jù)中存在的異常數(shù)據(jù)進(jìn)行預(yù)處理并采用T4253H 平滑算法[5]對(duì)原始速度時(shí)間序列進(jìn)行濾波處理。通過(guò)片段劃分并構(gòu)造片段多維特征參數(shù)，采用主成分分析方法實(shí)現(xiàn)參數(shù)降維，采用聚類效果評(píng)價(jià)方法確定聚類數(shù)目，之后利用K-means 聚類方法[6]將運(yùn)行片段分為4 種工況，基于馬爾科夫鏈原理[7]確定工況間的順序關(guān)系，最終結(jié)合怠速時(shí)間特征與各工況總時(shí)間比例確定汽車行駛工況。

1 數(shù)據(jù)預(yù)處理

本文采集的某城市原始載荷數(shù)據(jù)共140 h，3 800 km，包含不同道路類型、不同駕駛習(xí)慣的載荷時(shí)間歷程。先對(duì)原始載荷歷程進(jìn)行平滑過(guò)渡處理，使得處理后的時(shí)間序列連續(xù)，并且濾波后可以消除由于數(shù)據(jù)采集終端存在的零點(diǎn)漂移，即白噪聲造成一定的車速偏差，可以通過(guò)濾波處理后來(lái)降低數(shù)據(jù)誤差。

本文擬采用T4253H 濾波首先將數(shù)據(jù)按照跨度4，2，5，3 的順序進(jìn)行中心移動(dòng)平均處理，對(duì)于原始數(shù)據(jù){fi}中的N 個(gè)非平穩(wěn)數(shù)據(jù){yi}，認(rèn)定每隔4，2，5，3 個(gè)跨度的區(qū)間內(nèi)的數(shù)據(jù)都是平穩(wěn)的，其均值接近于常量。所以可以每隔4，2，5，3 的跨度取相鄰數(shù)據(jù)的均值，來(lái)表示每隔這么多個(gè)跨度的一個(gè)均值，通常用該均值表示數(shù)據(jù)的測(cè)量結(jié)果或信號(hào)數(shù)據(jù)。例如每隔5 個(gè)取一個(gè)均值的表達(dá)式為

同理，每隔4 個(gè)取值，即y4=1/4(y1+y2+y3+y4),以此類推，可得一般表達(dá)式為

式（2）中，k 可取4，2，5，3；i 取1，2，…，k。顯然，這樣所得到的{fk=yk}，原始數(shù)據(jù)經(jīng)過(guò)處理之后會(huì)更加平滑且無(wú)異常數(shù)據(jù)。但實(shí)際上，不同的數(shù)據(jù)距離平滑數(shù)據(jù)的遠(yuǎn)近會(huì)產(chǎn)生不等權(quán)的情況，因此對(duì)于不同的復(fù)雜數(shù)據(jù)，需要通過(guò)加權(quán)平均來(lái)表示平滑數(shù)據(jù)，所以這里就把Hanning 作為權(quán)重再做移動(dòng)平均處理。其一般算式為

式中：mi——權(quán)系數(shù)，且；p，q——小于4，2，5，3 的任一正整數(shù)。由此可以對(duì)數(shù)據(jù)進(jìn)行平滑處理。

此外，還需要對(duì)原始數(shù)據(jù)中存在的異常數(shù)據(jù)進(jìn)行處理。對(duì)于異常加減速的數(shù)據(jù),如圖1 所示，此種異常加減速數(shù)據(jù)影響片段劃分，應(yīng)將異常數(shù)據(jù)剔除；對(duì)于怠速時(shí)間大于180 s 的片段，如圖2 所示，通過(guò)編程將最長(zhǎng)怠速時(shí)間設(shè)置為180 s；由于設(shè)備異常出現(xiàn)無(wú)波動(dòng)的數(shù)據(jù)，如圖3 所示，速度時(shí)間歷程無(wú)波動(dòng)可能由于設(shè)備異常，該片段會(huì)使汽車運(yùn)動(dòng)學(xué)片段失真且對(duì)工況構(gòu)建影響較大，對(duì)其進(jìn)行刪除處理；對(duì)于緩行片段，如圖4所示。最高車速小于5 km/h 的緩行片段視為怠速工況處理。

圖1 異常加減速數(shù)據(jù)Fig.1 Abnormal acceleration and deceleration data

圖2 長(zhǎng)時(shí)間怠速數(shù)據(jù)Fig.2 Long-term idling data

圖3 異常無(wú)波動(dòng)數(shù)據(jù)Fig.3 Abnormal no-fluctuation data

圖4 緩行片段數(shù)據(jù)Fig.4 Slow driving data

2 片段劃分與特征參數(shù)

2.1 片段劃分

對(duì)預(yù)處理后的速度時(shí)間歷程進(jìn)行片段劃分。將速度從零開始，經(jīng)過(guò)加減速到零為止定義為運(yùn)行段；將速度持續(xù)為零或最高速度低于5 km/h 的片段定義為怠速段；運(yùn)動(dòng)學(xué)片段是汽車從怠速狀態(tài)開始到下一個(gè)怠速狀態(tài)開始之間的車速區(qū)間。具體片段劃分示意圖如圖5 所示。最終通過(guò)編程將載荷數(shù)據(jù)劃分出2 517 個(gè)運(yùn)動(dòng)學(xué)片段。

圖5 片段劃分示意圖Fig.5 Schematic diagram of segmentation

2.2 特征參數(shù)構(gòu)造

汽車在實(shí)際道路上的行駛過(guò)程需要通過(guò)一些基本的特征參數(shù)來(lái)反映各個(gè)運(yùn)動(dòng)學(xué)片段的特征。本文結(jié)合汽車運(yùn)行載荷數(shù)據(jù)，對(duì)所有運(yùn)行片段構(gòu)造如表1 中15 個(gè)特征參數(shù)，并對(duì)所有運(yùn)行片段計(jì)算出總體特征參數(shù)矩陣。其中，平均速度指包含怠速段的運(yùn)動(dòng)學(xué)片段的平均速度，平均行駛速度指不包含怠速段的運(yùn)行段的平均速度。

表1 特征參數(shù)說(shuō)明Tab.1 Characteristic parameter escription

3 行駛工況構(gòu)建

3.1 主成分分析

設(shè)有n 個(gè)片段，每個(gè)片段有p 個(gè)特征參數(shù)，記為X=（x1，x2，…，xp）。設(shè)隨機(jī)變量X 的均值為μ，協(xié)方差矩陣為∑。主成分是將p 個(gè)特征參數(shù)線性組合的問(wèn)題。對(duì)X 進(jìn)行線性變換可以生成新的綜合指標(biāo)即主成分，記為y1，y2，…，yp。樣本數(shù)據(jù)可表示為

主成分分析的數(shù)學(xué)模型為

從相關(guān)性矩陣出發(fā)求解主成分，在對(duì)特征參數(shù)矩陣進(jìn)行降維處理之前，由于各特征參數(shù)的量綱不同，造成特征參數(shù)值得分散程度會(huì)有很大差異，從而影響分析結(jié)果，因此需要對(duì)特征參數(shù)矩陣做標(biāo)準(zhǔn)化處理。

主成分的目的是減少指標(biāo)的個(gè)數(shù)，累積貢獻(xiàn)率達(dá)到80%以上。這樣的目的是保留大部分信息，同時(shí)減少指標(biāo)變量。根據(jù)總方差解釋（表2）提取前4 個(gè)主成分，累積貢獻(xiàn)率達(dá)到82.88%。

表2 總方差解釋Tab.2 Total variance explanation

根據(jù)因子載荷量與特征值的關(guān)系可以得出主成分系數(shù)矩陣。根據(jù)主成分系數(shù)矩陣乘以標(biāo)準(zhǔn)化后的特征參數(shù)矩陣即可得到主成分得分矩陣，如表3 所示。主成分得分矩陣可用于后續(xù)聚類分析。

表3 主成分得分矩陣Tab.3 Principal component score matrix

3.2 K-means 聚類分析

3.2.1 聚類數(shù)目確定

本文采用肘部法則（Elbow method）確定聚類數(shù)目（k 值）。肘部法則主要是隨著聚類數(shù)目k 值的增大，尋找損失函數(shù)的拐點(diǎn)，損失函數(shù)為各個(gè)類別畸變程度之和，也就是誤差平方和（Sum of squared errors，SSE）。

式中：Ci——第 i 個(gè)簇；p——Ci中的樣本點(diǎn)；mi——Ci的質(zhì)心；SSE——聚類效果的好壞。

當(dāng)聚類數(shù)目為2～6 時(shí)，計(jì)算出所有情況的誤差平方和，如圖6 所示。當(dāng)聚類數(shù)目為4 時(shí)出現(xiàn)拐點(diǎn)，說(shuō)明聚類分為4 類較為合理。

圖6 誤差平方和計(jì)算結(jié)果Fig.6 Calculation results of sum of squared errors

3.2.2 聚類結(jié)果

將降維后的主成分得分矩陣進(jìn)行聚類分析，預(yù)先確定聚類數(shù)目之后，采用K-means 聚類算法將運(yùn)行學(xué)片段進(jìn)行分類，基于前3 個(gè)主成分繪制出聚類云圖如圖7 所示。

圖7 聚類結(jié)果示意圖Fig.7 Schematic diagram of clustering results

聚類后第1 類包含526 個(gè)片段，第2 類包含637 個(gè)片段，第3 類包含460 個(gè)片段，第4 類包含894 個(gè)片段。并根據(jù)聚類后的結(jié)果，計(jì)算出每一類典型特征參數(shù)的平均值代表每類工況的參數(shù)特征，如表4 所示。

表4 聚類后各類均值特征參數(shù)Tab.4 Various types of mean feature parameters after clustering

聚類后的4 種工況特征值具有一定的差異性：類別3 屬于高速行駛片段；類別2 屬于中速行駛片段；類別1 與類別4 片段最高車速與行駛時(shí)間較為接近，但類別1 的怠速時(shí)間較長(zhǎng)，類別4 怠速時(shí)間較短，行駛里程較高?？梢哉J(rèn)為類別1 為擁堵道路片段，類別4 為暢通道路片段。

3.3 典型工況分析與提取

基于聚類后各類的結(jié)果，對(duì)每類工況所包含的時(shí)間與里程比例進(jìn)行統(tǒng)計(jì)分析，作為拼接片段的參考依據(jù)?；诘∷贂r(shí)間比例確定每類工況怠速時(shí)間，基于各類工況總時(shí)間比例確定出最終要拼接的各類的時(shí)間。聚類后四種工況的總時(shí)間比例及怠速時(shí)間比例如圖8、圖9 所示。

圖8 各工況總時(shí)間比例Fig.8 Proportion of total time under various working conditions

圖9 各工況怠速時(shí)間比例Fig.9 Proportion of idle time under various working conditions

將每類工況到各自聚類中心較近的幾個(gè)片段作為候選片段。由于通過(guò)K-means 聚類分析后，聚類中心處的片段與均值特征參數(shù)對(duì)應(yīng)的片段一致，因此，根據(jù)表4 各類工況的特征參數(shù)的平均值確定出每類工況的最優(yōu)片段。

3.4 行駛工況構(gòu)建

本文采用馬爾科夫鏈原理確定各工況狀態(tài)轉(zhuǎn)移關(guān)系。馬爾可夫鏈?zhǔn)谴硪唤M離散隨機(jī)變量的集合。具體指對(duì)概率空間（Ω，F(xiàn)，P）內(nèi)以一維可數(shù)集為指數(shù)集的隨機(jī)變量集合X={Xn：n>0}，若隨機(jī)變量的取值都在可數(shù)集內(nèi)：X=si，si∈S，且隨機(jī)變量的條件概率滿足如下關(guān)系：

則X 被稱為馬爾可夫鏈。對(duì)于一個(gè)靜態(tài)的馬爾可夫鏈，隨著馬爾可夫鏈的增長(zhǎng)，整個(gè)模型事件分布保持不變。基于馬爾可夫鏈的這個(gè)性質(zhì)，通過(guò)模型事件轉(zhuǎn)移矩陣構(gòu)建的行駛工況能夠代表整個(gè)原始數(shù)據(jù)的行駛工況，利用最大似然估計(jì)法來(lái)估計(jì)每類工況的轉(zhuǎn)移矩陣。

基于馬爾科夫鏈原理對(duì)4 類工況的拼接順序進(jìn)行確定。由于本文在數(shù)據(jù)預(yù)處理過(guò)程中將異常數(shù)據(jù)進(jìn)行剔除，并基于預(yù)處理后的數(shù)據(jù)進(jìn)行片段劃分，劃分的各片段之間依次連續(xù)?；诰垲惙治龅慕Y(jié)果，可以確定每個(gè)片段所屬的工況類別，從而構(gòu)成馬爾科夫鏈模型，而4 種工況類別構(gòu)成狀態(tài)空間，通過(guò)馬爾科夫鏈原理計(jì)算出各工況間的狀態(tài)轉(zhuǎn)移矩陣，如表5 所示。

表5 工況狀態(tài)轉(zhuǎn)移概率矩陣Tab.5 State transition probability matrix

從所有工況中隨機(jī)選取一個(gè)工況作為初始工況，基于狀態(tài)轉(zhuǎn)移概率分布構(gòu)造隨機(jī)數(shù)，若該隨機(jī)數(shù)大于第r 工況下的轉(zhuǎn)移概率，小于第q 工況下的轉(zhuǎn)移概率，則下一個(gè)工況為第r 個(gè)工況，以此進(jìn)行工況拼接。由于隨機(jī)選取工況4 的概率最高，為35.5%，因此得出各工況拼接關(guān)系依次為工況4、工況1、工況2、工況3。最終結(jié)合各工況間時(shí)間比例關(guān)系構(gòu)建出汽車行駛工況，如圖10所示，總時(shí)間為1 280 s。

圖10 構(gòu)建的汽車行駛工況Fig.10 Constructed vehicle driving conditions

4 汽車行駛工況驗(yàn)證

為驗(yàn)證所構(gòu)建的行駛工況是否合理，需計(jì)算汽車行駛工況的特征參數(shù)與經(jīng)處理后數(shù)據(jù)源的特征參數(shù)之間的相對(duì)誤差。相對(duì)誤差計(jì)算公式如下：

式中：ε——相對(duì)誤差；xi——構(gòu)建的汽車行駛工況特征參數(shù)值；xj——經(jīng)處理后數(shù)據(jù)源的特征參數(shù)值。選取較為典型的特征參數(shù)并計(jì)算構(gòu)建好的行駛工況的特征值與預(yù)處理之后的特征值進(jìn)行對(duì)比。如圖11 所示，相對(duì)誤差基本上在10%以內(nèi)，在合理范圍之內(nèi)。

圖11 行駛工況特征參數(shù)誤差對(duì)比Fig.11 Comparison of characteristic parameter errors in driving conditions

將原始預(yù)處理后的數(shù)據(jù)與構(gòu)建的汽車行駛工況進(jìn)行速度與加速度聯(lián)合分布概率差異對(duì)比，如圖12 所示。差異值基本在±2%之間，比較合理。

圖12 速度-加速度聯(lián)合分布概率差異值Fig.12 Probability difference value of joint velocity-acceleration distribution

5 結(jié)論

對(duì)汽車實(shí)際道路行駛數(shù)據(jù)進(jìn)行了預(yù)處理，基于程序?qū)崿F(xiàn)運(yùn)動(dòng)學(xué)片段的自動(dòng)劃分，并構(gòu)造表征片段的特征參數(shù)，基于主成分原理實(shí)現(xiàn)參數(shù)降維，采用K-means 聚類方法將運(yùn)動(dòng)學(xué)片段分為4 類，且聚類效果良好?；隈R爾科夫鏈原理構(gòu)建各工況狀態(tài)轉(zhuǎn)移概率矩陣，結(jié)合工況時(shí)間比例關(guān)系構(gòu)建汽車行駛工況。

根據(jù)特征參數(shù)相對(duì)誤差檢驗(yàn)及速度-加速度聯(lián)合分布概率差異值進(jìn)行對(duì)比，確定最終組合后的行駛工況，相對(duì)誤差在10%以內(nèi)，速度-加速度聯(lián)合概率分布差異值在±2%以內(nèi)，效果良好。