基于馬爾可夫鏈的車輛規(guī)劃行駛工況構(gòu)建方法研究*

王偉 曲輔凡 楊鈁 李文博

（1.中汽研汽車檢驗(yàn)中心（天津）有限公司，天津 300300；2.中國第一汽車股份有限公司研發(fā)總院，長春 130013）

主題詞：工況構(gòu)建 聚類分析 馬爾可夫鏈 在線地圖

1 前言

車輛行駛工況是車輛性能提升的共性基礎(chǔ)，對于車輛能耗優(yōu)化、控制策略標(biāo)定及新技術(shù)的開發(fā)具有重要意義。當(dāng)前，應(yīng)用最廣泛的車輛行駛工況是標(biāo)準(zhǔn)循環(huán)工況，但受地域差異、道路差異影響，標(biāo)準(zhǔn)循環(huán)工況與車輛實(shí)際工況差異較大[1]。因此有必要研究一種能夠規(guī)劃、預(yù)測車輛實(shí)際行駛狀態(tài)的工況構(gòu)建方法。

國內(nèi)外學(xué)者在車輛工況構(gòu)建方法方面進(jìn)行了大量的研究，工況構(gòu)建包含行駛數(shù)據(jù)壓縮重構(gòu)和工況合成2個部分[2]。數(shù)據(jù)壓縮重構(gòu)過程將運(yùn)行數(shù)據(jù)視為離散片段，提取速度、加速度、行駛距離等典型參數(shù)，通過主成分分析（Principal Component Analysis，PCA）法對數(shù)據(jù)進(jìn)行降維，保留原始數(shù)據(jù)的核心信息[3]，如尹安東等以合肥市純電動汽車的行駛工況為對象，采用主成分分析法構(gòu)建行駛工況，對不同類型工況的汽車能量消耗量進(jìn)行估算，但是該方法忽略了車輛實(shí)際行駛過程中的不確定性，能耗預(yù)測結(jié)果有待驗(yàn)證。工況合成方法包含聚類分析法、馬爾可夫鏈法、V-A 矩陣法等算法[4-5]。秦大同等采用K 均值（K-Means）聚類分析法構(gòu)造車輛能耗與行駛特征的關(guān)聯(lián)關(guān)系，合成某山地城市的典型行駛工況，但是忽略了實(shí)時交通狀態(tài)的影響，難以代表實(shí)際行駛特征[6-7]。Hereijgers 等提出一種同時考慮汽車行駛速度和道路坡度的多維馬爾可夫鏈的工況合成方法，構(gòu)建的工況更符合真實(shí)路況數(shù)據(jù)，但是該方法以車速作為依據(jù)進(jìn)行整體分析，劃分為加速、減速、勻速、怠速（或停車）4 種狀態(tài)[8-9]，導(dǎo)致該方法合成的工況受初始隨機(jī)車速的影響較大，難以反映車輛實(shí)際行駛狀態(tài)。

綜上，當(dāng)前工況構(gòu)建方法均為針對歷史數(shù)據(jù)的工況構(gòu)建方法，沒有融合道路交通實(shí)時運(yùn)行狀態(tài)，無法反映車輛未來行駛工況。本文提出一種基于馬爾可夫鏈的車輛規(guī)劃行駛工況構(gòu)建方法，采集車輛行駛數(shù)據(jù)，采用馬爾可夫鏈構(gòu)建工況片段庫，將在線地圖規(guī)劃路徑與工況片段融合，最終生成規(guī)劃行駛工況。

2 基于馬爾可夫鏈的在線地圖規(guī)劃工況構(gòu)建

工況片段庫基于實(shí)際行駛數(shù)據(jù)構(gòu)建，用于在線地圖行駛工況規(guī)劃。本文提出的工況構(gòu)建流程如圖1所示。

圖1 在線地圖規(guī)劃行駛工況構(gòu)建流程

2.1 車輛行駛工況片段庫構(gòu)建

2.1.1 行駛數(shù)據(jù)清洗

道路和環(huán)境等因素的干擾會造成采集數(shù)據(jù)的異常，甚至缺失，故需對采集的數(shù)據(jù)進(jìn)行清洗。異常數(shù)據(jù)主要體現(xiàn)在停車時間或怠速時間、加速度和車速等方面。數(shù)據(jù)清洗后需采用遞推平滑濾波算法對原始數(shù)據(jù)進(jìn)行濾波處理，降低異常因素對車速數(shù)據(jù)的影響。

由于車速數(shù)據(jù)連續(xù)采集的時間較長，行駛特征難以準(zhǔn)確表達(dá)，因此將數(shù)據(jù)分段處理。城市道路行駛車輛因受到交通條件的限制，存在停車、加速、減速、勻速狀態(tài)間的頻繁切換，因此將從一次停車或者怠速的起始時刻至下一次停車或者怠速的起始時刻定義為一個交通片段[10-11]。通常，每個交通片段至少包含1個加速狀態(tài)和1個減速狀態(tài)，時長不少于10 s，加速度一般處于-4～4 m/s2范圍內(nèi)。將采集的數(shù)據(jù)分為數(shù)量為N的交通片段，選擇如表1 所示的特征參數(shù)描述每個片段的加速、減速、勻速、停車（或怠速）等行駛特征，其中平均行駛速度為行駛過程中的平均車速。

表1 片段特征參數(shù)

2.1.2 數(shù)據(jù)降維

經(jīng)過上述處理后得到的交通片段數(shù)據(jù)量大、維度高，且交通片段各特征參數(shù)間具有一定的相關(guān)性，直接聚類存在計算量大、聚類效果差等問題，因此需要將特征參數(shù)降維[12-13]。本文采用主成分分析法對車輛行駛數(shù)據(jù)進(jìn)行降維處理，將原始變量線性組合成相互獨(dú)立的新變量。

將所有片段的特征參數(shù)記為特征矩陣X：

式中，n為交通片段的數(shù)量；p為每個片段特征參數(shù)的數(shù)量。

對X進(jìn)行標(biāo)準(zhǔn)化處理，得到標(biāo)準(zhǔn)化矩陣S，其元素Sij為：

處理后的矩陣S并不改變各參量的相關(guān)性。根據(jù)矩陣S計算出相關(guān)系數(shù)矩陣及其特征值和特征參數(shù)，基于特征值計算出累計貢獻(xiàn)率。選擇特征值大于1，且累積貢獻(xiàn)率在80%以上的主成分，用于行駛工況構(gòu)建。

2.1.3 K-Means聚類

采用K-Means聚類算法對主成分進(jìn)行聚類，具體步驟為[14-15]：

a.根據(jù)實(shí)際問題確定分類數(shù)量k，并在每一類中確定初始聚類中心；

b.計算各樣本到聚類中心的距離，距離較近的歸入一類，歐式距離的計算公式為：

式中，xik為第i個樣本的第k個主成分；xjk為聚類中心j的第k個主成分；q為主成分的數(shù)量。

c.計算每一類的中心位置，并將該位置確定為新的聚類中心；

d.按照新的聚類中心重新進(jìn)行分類，重復(fù)步驟b、步驟c 的操作，隨著重復(fù)次數(shù)的增多，聚類中心不再發(fā)生大的偏移，聚類結(jié)果趨于穩(wěn)定。

按上述步驟，獲得了類別數(shù)量為n的片段樣本、每類樣本包含的片段數(shù)量及特征參數(shù)的均值，用于馬爾可夫鏈工況構(gòu)建。

2.1.4 基于馬爾可夫鏈的片段庫構(gòu)建

基于測試數(shù)據(jù)構(gòu)建的片段庫為在線地圖規(guī)劃行駛工況構(gòu)建提供基礎(chǔ)，因此需要選用合理的方法構(gòu)建片段庫。馬爾可夫鏈能夠基于歷史交通數(shù)據(jù)較好地預(yù)測未來行駛工況，且車輛行駛特征符合馬爾可夫性，因此采用馬爾可夫鏈構(gòu)建片段庫。采用臨近法將車輛行駛狀態(tài)離散為有限的狀態(tài)[16]：

在當(dāng)前時刻行駛狀態(tài)v(k)為vi時，下一時刻行駛狀態(tài)v(k+1)為vj的概率為[17]：

式中，Pi,j為狀態(tài)轉(zhuǎn)移概率矩陣P的第i行第j列元素；概率轉(zhuǎn)移矩陣P為：

狀態(tài)轉(zhuǎn)移矩陣的元素采用極大似然法求得：

式中，F(xiàn)ij為車輛從狀態(tài)vi轉(zhuǎn)移到vj的次數(shù)；Fk為當(dāng)前狀態(tài)的總統(tǒng)計量。

基于馬爾可夫鏈構(gòu)建工況的過程為：設(shè)定當(dāng)前狀態(tài)為v1，生成隨機(jī)數(shù)λ∈(0,1)，通過λ與狀態(tài)轉(zhuǎn)移矩陣P的相應(yīng)元素大小進(jìn)行對比，確定下一狀態(tài)v2，如此重復(fù)計算，直至累計距離達(dá)到目標(biāo)車速。

構(gòu)建的片段庫用于在線地圖工況構(gòu)建，而在線地圖片段長度、平均車速等具有不確定性，因此構(gòu)建的片段庫應(yīng)當(dāng)盡可能覆蓋全部行駛車速和距離。

2.2 在線地圖規(guī)劃工況融合

構(gòu)建在線地圖規(guī)劃行駛工況需要基于地圖應(yīng)用程序編程接口（Application Programming Interface，API）數(shù)據(jù)生成基礎(chǔ)規(guī)劃工況，然后將構(gòu)建的片段庫與基礎(chǔ)規(guī)劃工況融合，形成能預(yù)測實(shí)際行駛特征的工況。

2.2.1 在線地圖基礎(chǔ)規(guī)劃工況

在地圖中規(guī)劃路徑，地圖API可獲取行駛距離和預(yù)計行駛時間，并可獲取在規(guī)劃路徑片段的距離和通過時間、交通燈的位置，通過以下計算得出路徑片段的目標(biāo)車速。

首先根據(jù)交通燈位置以及車輛通過該片段的預(yù)估時間，確定在交通燈處是否停車以及停車時間。之后利用片段長度與通過時間計算片段平均車速，并根據(jù)交通片段長度與擁堵系數(shù)進(jìn)行修正。最后識別連續(xù)交通片段中車速相等的勻速段，觀察勻速段內(nèi)是否包含停車點(diǎn)。如果未包含，則目標(biāo)車速為該勻速段的車速。如果包含，則統(tǒng)計該片段內(nèi)從0加速的次數(shù)n1、減速到0的次數(shù)n2和停車時間，計算目標(biāo)車速Vt：

式中，a為平均加速度；d為平均減速度；T為通過時間；Ns為停車次數(shù)；ti為每次停車時間；D為勻速段的距離。

計算目標(biāo)車速內(nèi)的速度突變點(diǎn)需識別車速突變點(diǎn)的初始速度和終止速度，計算速度突變點(diǎn)的加速、減速過程需求車速及距離：

式中，V為目標(biāo)車速；X為行駛距離；X0為初始距離；V0為初始位置的初始車速；t為加速、減速過程的時間。

通過上述計算得出基于地圖API的基礎(chǔ)規(guī)劃工況，僅包含勻速和停車穩(wěn)態(tài)工況，在細(xì)分的片段中不能進(jìn)一步規(guī)劃出車輛實(shí)際行駛工況。

2.2.2 馬爾可夫鏈片段融合

基于在線地圖需求信息，采用貪心算法對每個交通小片段選取最優(yōu)的工況片段：

式中，Δl為地圖交通小片段長度與k1個工況片段總長度的差值；k1為工況片段數(shù)量；d(i)為第i個工況片段的距離；L為地圖交通小片段的長度。

以交通小片段的通過時間與工況片段的差值最小為尋優(yōu)目標(biāo)，且差值Δt應(yīng)不小于0：

式中，t(i)為第i個工況片段的時間；T1為交通小片段的通過時間。

由于地圖交通小片段的長度和時間存在不確定性，工況片段尋優(yōu)可能出現(xiàn)3種情況：

a.Δl≥0 且Δt>0，則在Δt時間內(nèi)補(bǔ)充車速為Δl/Δt的勻速段；

b.Δl>0 且Δt=0，則在時間T1內(nèi)，將工況片段速度增大L/(L-Δl)倍；

c.Δl=0且Δt=0，則無需處理。

通過上述實(shí)際采集的數(shù)據(jù)所構(gòu)造的工況片段，以規(guī)劃路徑片段的距離和平均車速作為約束條件，從片段庫中選取、組合片段，最終實(shí)現(xiàn)了生成的工況與規(guī)劃路徑的距離和時間相等。

3 驗(yàn)證分析

在試驗(yàn)車內(nèi)安裝車速傳感器、加速度傳感器及數(shù)據(jù)采集設(shè)備（未采集交叉口、交通燈及車道特征等信號），采集同一駕駛員在某城市駕駛車輛連續(xù)3 個月的車輛行駛數(shù)據(jù)，行駛道路涉及城市工況、郊區(qū)工況、快速路工況等。降噪后的數(shù)據(jù)如圖2所示，可以看出濾波后的數(shù)據(jù)平滑度更高。將濾波后的數(shù)據(jù)劃分為11 587個片段，計算出各片段的特征參數(shù)如表2 所示。再按照本文提出的方法構(gòu)建在線地圖基礎(chǔ)規(guī)劃工況、在線地圖規(guī)劃行駛工況，采集車輛實(shí)際行駛工況。以某款電動車為研究對象，將3種工況導(dǎo)入新能源汽車數(shù)字化虛擬仿真平臺（VPAT2021），進(jìn)行驗(yàn)證分析。

表2 片段特征參數(shù)

各片段的特征參數(shù)在主成分分析后得到的相應(yīng)特征值、貢獻(xiàn)率及累計貢獻(xiàn)率如表3所示，前3個主成分的特征值均大于1，累計貢獻(xiàn)率為84.27%。因此選取主成分M1、M2、M3 來表征15 個特征參數(shù)所包含的行駛信息。

表3 主成分分析結(jié)果

對3 項(xiàng)主成分進(jìn)行K-Means 聚類分析，得出4 類片段樣本集，各類樣本的部分特征參數(shù)如表4所示，由表4中速度和加速度可知，類別1～類別4 分別代表城郊工況、高速工況、城市工況、嚴(yán)重?fù)矶鹿r。

表4 聚類特征參數(shù)（部分）

基于聚類結(jié)果和行駛數(shù)據(jù)，采用式（6）計算車輛在不同類之間切換的概率矩陣P1，如表5所示。

表5 類狀態(tài)轉(zhuǎn)移概率P1 %

在聚類樣本內(nèi)部，車輛存在怠速（或停車）、加速、勻速、減速4種行駛狀態(tài)，不同類所代表的車速、加減速度等行駛狀態(tài)不同，行駛狀態(tài)轉(zhuǎn)移概率也不同。以類別1為例，采用式（6）計算出各行駛狀態(tài)轉(zhuǎn)移概率矩陣P2如表6所示。

表6 類別1行駛狀態(tài)轉(zhuǎn)移概率P2 %

基于概率轉(zhuǎn)移矩陣，利用馬爾可夫鏈合成片段的步驟為：

a.選擇類別3 為初始類狀態(tài)，選擇怠速或停車狀態(tài)為初始行駛狀態(tài)；

b.根據(jù)類狀態(tài)轉(zhuǎn)移矩陣P1和隨機(jī)數(shù)λ選擇接下來的類狀態(tài)；

c.根據(jù)該類的行駛狀態(tài)轉(zhuǎn)移矩陣P2和隨機(jī)數(shù)λ選擇行駛狀態(tài)；

d.重復(fù)步驟b、步驟c，組合構(gòu)建的片段，直到滿足片段距離要求。

計算平均車速是否達(dá)到目標(biāo)平均車速，將達(dá)到要求的片段列入候選工況庫，然后對比候選工況庫的平均誤差，選擇平均誤差最小的為典型片段?；隈R爾可夫鏈構(gòu)建的片段示例如圖3所示。

圖3 基于馬爾可夫鏈的車速片段

在地圖上規(guī)劃車輛實(shí)際行駛的出發(fā)點(diǎn)和終點(diǎn)，生成在線地圖工況，從工況片段庫選取片段，融合在線地圖工況，形成規(guī)劃行駛工況，如圖4所示，與車輛實(shí)際行駛狀態(tài)相比，構(gòu)造的工況特征值如表7所示。

圖4 規(guī)劃行駛工況

由表7可見，基于馬爾可夫鏈片段庫構(gòu)建的工況與實(shí)際行駛工況的平均誤差僅為4.29%，最大相對誤差為9.09%。在線地圖基礎(chǔ)規(guī)劃工況以勻速為主，因此最終與實(shí)際行駛工況的平均誤差達(dá)到316.40%，明顯高于前者。故基于馬爾可夫鏈片段庫的車輛規(guī)劃工況構(gòu)建方法能夠較好地反映車輛實(shí)際行駛狀態(tài)。

在VPAT2021中搭建某款電動車的模型，車輛基本參數(shù)如表8所示。

表8 車輛基本參數(shù)

進(jìn)行工況能耗仿真對比分析，仿真過程如圖5 所示，仿真結(jié)果如表9 所示，基于馬爾可夫鏈片段庫規(guī)劃的車輛工況與實(shí)際行駛工況的誤差更小，百公里電耗為12.73 kW·h，誤差僅為4.09%。

圖5 VPAT2021能耗仿真過程

表9 能耗仿真對比

4 結(jié)束語

本文基于采集的實(shí)際交通數(shù)據(jù)，經(jīng)過數(shù)據(jù)清洗、劃分，計算出多個交通片段的特征參數(shù)，通過PCA法和KMeans聚類法對數(shù)據(jù)進(jìn)行分析，利用馬爾可夫鏈構(gòu)建了片段庫，與在線地圖基礎(chǔ)規(guī)劃工況融合，構(gòu)建了在線地圖規(guī)劃行駛工況。

將在線地圖規(guī)劃行駛工況與車輛實(shí)際行駛工況進(jìn)行對比，采用特征參數(shù)和仿真能耗2 個維度評價，發(fā)現(xiàn)采用馬爾可夫鏈片段庫與在線地圖融合的方式能更準(zhǔn)確地反映車輛實(shí)際行駛特征，工況特征參數(shù)平均誤差為4.29%，能耗仿真誤差為4.09%。