聚類和馬爾科夫方法結(jié)合的城市汽車行駛工況構(gòu)建

姜 平 石 琴 陳無畏

合肥工業(yè)大學,合肥,230009

0 引言

汽車行駛工況分為標準工況和非標準工況。當今世界上很多國家都以標準等形式提出了不同車型在不同應(yīng)用時的標準工況,如日本的10-15工況循環(huán)、歐洲經(jīng)濟委員會的ECE-R15工況循環(huán)、美國聯(lián)邦城市及高速公路循環(huán)(CVS-C/H)。我國汽車排放測試規(guī)程采用歐洲的排放標準,往往帶來很不滿意的結(jié)果,因此,需要對待開發(fā)的目標車型進行工況調(diào)查和測試,并得到合適的行駛工況。所以,開發(fā)和設(shè)計汽車、考核汽車排放污染水平都需要有相應(yīng)的行駛工況。

目前,國內(nèi)學者大多采用單純的短行程和聚類分析法對汽車行駛工況進行研究[1-3]。國外很多學者采用“微路徑”以及隨機過程選擇方法來構(gòu)建行駛工況[4-5]。本文在構(gòu)建汽車行駛工況時,不同于國內(nèi)傳統(tǒng)方法。傳統(tǒng)方法未考慮汽車燃油消耗與平均速度、發(fā)動機工作狀態(tài)分布等與行駛工況有關(guān)的因素,只是籠統(tǒng)地對整個汽車行駛工況進行構(gòu)建。本文將聚類分析法、因子分析法與馬爾科夫過程相結(jié)合,按照平均速度和速度標準偏差對汽車行駛工況進行類的劃分,將每類行駛工況分別看作一個馬爾科夫過程,計算每類行駛工況的片段轉(zhuǎn)移概率和發(fā)動機比功率,提出構(gòu)建候選工況的新方法。

1 因子分析

因子分析的基本目的就是用少數(shù)幾個因子描述許多指標或因素之間的關(guān)系,把變量表示成各因子的線性組合,即將關(guān)系比較密切的幾個變量歸在同一類中,每一類變量就成為一個因子,從而以較少的幾個因子反映原資料的大部分信息。

發(fā)動機工作狀態(tài)分布一般根據(jù)發(fā)動機比功率、汽車行駛速度、瞬時加速度決定,而且它又與汽車排放、燃油消耗密切相關(guān),因此在構(gòu)建汽車行駛工況時,應(yīng)考慮發(fā)動機工作狀態(tài)分布。比功率綜合了速度、加速度、坡度以及風阻等參數(shù),其物理意義為發(fā)動機輸出功率與機動車質(zhì)量的比值,其擬合公式如下[6]：

式中,v為車輛行駛速度,m/s;a為車輛行駛瞬態(tài)加速度,m/s2;θ為道路坡度,(°);Pa為發(fā)動機比功率,k W/t。

參考國外發(fā)動機工作狀態(tài)分布[7],考慮國內(nèi)混合交通實際情況,將我國城市典型道路行駛工況燃油效率劃分為表1所示的分布。表1的工作狀態(tài)分布是一個14維變量(B1～B14),這些變量分別由速度分布、比功率來確定。每個變量類似于一個單元,每一個單元類似于一個“小抽屜”,在這個“小抽屜”中,存放著該工況下測點的平均燃油消耗。行駛工況是逐秒采集的,每秒所觀察的數(shù)據(jù)被分配在表1中。如B1表示該工況下汽車制動時平均燃油消耗量,B2表示該工況下汽車怠速時平均燃油消耗量,B3表示該工況下汽車比功率小于0且速度大于0、小于40km/h時的平均燃油消耗量,依此類推,其他變量表示在該工況下不同比功率和速度區(qū)間的平均燃油消耗量。

表1 行駛工況燃油消耗分布的劃分

本文應(yīng)用統(tǒng)計分析軟件SPSS15.0中的因子分析方法來減少變量,采用主成分分析法提取公因子,提取特征值大于1的成分為主成分。因子模型的旋轉(zhuǎn)方法選取方差最大變換法。對合肥市5條典型道路工作狀態(tài)分布進行主成分分析,特征碎石圖見圖 1。前 2個主成分的特征值為1.47,累計方差貢獻率為84.3%,即前2個主成分包含了原有變量84.3%的信息,所以可以取這2個因子來代替原有的變量。

圖1 特征值碎石圖

通過對因子模型進行旋轉(zhuǎn),得到旋轉(zhuǎn)后的因子載荷矩陣如表2所示。根據(jù)旋轉(zhuǎn)后的因子載荷大小來判別因子主要解釋的變量,本文選取載荷大于0.8的變量為因子的主要解釋變量。因子1主要解釋了變量 B9、B10、B11、B12、B13 、B14,由表1可知,這些變量主要是高速狀態(tài)下的變量。低速狀態(tài)變量因子載荷小,如 B0、B1、B2、B3、B4、B5、B6、B7,說明因子 1與平均速度有關(guān),可以解釋為平均速度。因子2主要解釋了變量B1、B3、B4、B5,對怠速狀態(tài)變量B2影響最小,載荷在低速時先增大,中間速度時減小,高速時先減小后增大,這表明因子2與車速偏離平均速度的程度有關(guān),可以解釋為速度的標準偏差。

表2 旋轉(zhuǎn)后的因子載荷矩陣

2 兩步聚類分析

兩步聚類分析是一種新型的分層聚類方法,目前一般應(yīng)用在數(shù)據(jù)挖掘與多元統(tǒng)計的交叉領(lǐng)域——模式分類中,其算法適合任何尺度的變量[8]。采用上述與燃油消耗相關(guān)的兩個主要因子(平均速度和速度標準偏差)進行兩步聚類分析。聚類準則采用貝葉斯信息準則 Akaike(AIC)。AIC利用N-L平均信息量最小的原則來確定系統(tǒng)的階次和估計相應(yīng)參數(shù)值,其定階準則為[9]：

式中,N為模型獨立未知參數(shù)的個數(shù);θML為相應(yīng)參數(shù)的估計值;L(θML)為在θML條件下的似然函數(shù)值。

在工程應(yīng)用中,對于不同階次 N,首先用極大似然估計法估計系統(tǒng)的參數(shù)值,然后計算似然函數(shù)值及其AIC(N)值,以使AIC(N)值達到最小時的階次作為該模型的階次。本文利用SPSS15.0統(tǒng)計軟件的兩步聚類模型,采用AIC準則對合肥市典型道路的行駛工況進行聚類分析,結(jié)果可以聚類為3類：第1類行駛工況(類1)平均速度偏低,為0～20km/h,第2類行駛工況(類2)平均速度中等,為 20～40km/h,第3類行駛工況(類3)平均速度偏高,大于 40km/h。

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

3.1 樣本的采集

本研究中試驗車為轎車,試驗車盡可能跟隨測試路段上的車流正常行駛。速度信號采樣頻率為1Hz。城市道路行駛工況的研究一般是從多條道路中篩選出有代表性的典型道路進行研究。典型道路的特征是能綜合反映整個城市的主要交通狀況,即反映整個城市道路的道路等級(快速道、主干道、次干道及支路,以及車道數(shù)、機動車非機動車混合等),交通強度(交通流量、飽和度等),交叉口密度(路段內(nèi)交叉形式和數(shù)量等)和時間(早晚交通流量變化明顯程度等)。本文選取了合肥市勝利路、屯溪路、明光路等5條道路進行樣本采集,使采集的數(shù)據(jù)能反映合肥市整體道路的交通情況,統(tǒng)計的汽車交通流量具有代表性。

3.2 數(shù)據(jù)的分離

數(shù)據(jù)分離就是把行駛工況的速度時間曲線劃分成許多反映加速度變化的模型事件,包括加速模型事件、減速模型事件、怠速模型事件和勻速模型事件,這些模型事件是構(gòu)建行駛工況的基礎(chǔ)。模型事件使用(a1,a2,δ,n)準則進行定義。為了確定分類的最優(yōu)值,需要考慮以下目標函數(shù)和約束條件。

目標函數(shù)CAP是一個加權(quán)函數(shù)：

式中,k1,k2,k3,k4,k5為目標函數(shù)的權(quán)值。

為實現(xiàn)上述目標,要滿足以下5個約束條件：①加速度模型事件含有連續(xù)的加速度,因此,加速標準差θ1應(yīng)該很小;②對于任何勻速事件,最大速度和最低速度差θ2應(yīng)該小;③對于任何加速或減速事件,最大速度和最低速度差θ3應(yīng)大;④任何勻速事件的平均加速度θ4接近0;⑤因為一個單獨點的特征依賴于相鄰點特征,所以一個好的分類標準應(yīng)該包括一小部分單秒數(shù)據(jù)θ5。

通過迭代法使目標函數(shù)CAP最小,得到最優(yōu)的分類準則——(5,1,3,2)分類原則[10],即加速事件定義為單一瞬時加速度不小于2.22m/s2或者連續(xù)觀測加速度大于0.71 m/s2,或持續(xù)3s或以上時間,而且要實現(xiàn)3.2km/h以上的速度遞增;減速片段定義為單一瞬時加速度不大于-2.22m/s2,或任何連續(xù)觀測加速度不大于-0.71m/s2,持續(xù) 3s或以上時間,而且實現(xiàn)3.2km/h以上的速度遞減。其余的被列為怠速或勻速事件。這種分類原則考慮了現(xiàn)實世界的行駛工況,勻速模式事件可以包括一些小的加速和減速模型事件。按照這個分類原則,行駛工況片段劃分結(jié)果如圖2所示。怠速模型事件指速度為0的模型事件。

圖2 速度時間曲線片段劃分結(jié)果

3.3 轉(zhuǎn)移矩陣的計算

馬爾科夫過程是一個隨機過程,記為Zτ,τ=1,2,…,T,Zτ代表每個狀態(tài)的模型事件。當前時刻τ狀態(tài)的概率只與前一個時刻τ-1狀態(tài)有關(guān)。換句話說,汽車下個時間τ+1的狀態(tài)(如加速、減速、勻速和怠速片段)的概率只與當前時間τ的加速、減速、勻速和怠速片段有關(guān)。

在馬爾科夫鏈中,假設(shè)X是方程X=pX時(p為轉(zhuǎn)移矩陣,X為p的特征向量)的唯一解,則X為馬爾科夫鏈的靜態(tài)分布。對于一個靜態(tài)的馬爾科夫鏈,隨著馬爾科夫鏈的增長,整個模型事件分布保持不變。應(yīng)用這個性質(zhì),根據(jù)模型事件轉(zhuǎn)移矩陣構(gòu)建的代表性行駛工況能代表整個試驗數(shù)據(jù)的行駛工況。

為了估計每類行駛工況的轉(zhuǎn)移矩陣,運用Lee等[11]的最大似然估計法進行估計。對于一個固定的馬爾科夫過程,應(yīng)用貝葉斯公式得到連續(xù)的有序狀態(tài)的概率為

假定每個過程都是重復觀察的,N rs為τ-1時刻狀態(tài)r轉(zhuǎn)移到τ時刻狀態(tài)s的事件數(shù),則N次重復觀察的公式為

根據(jù)極大似然函數(shù),得到狀態(tài)轉(zhuǎn)移概率方程

依據(jù)上述原理,得到合肥市典型道路不同類別的行駛工況的片段轉(zhuǎn)移概率,如表3所示。例如,在類1工況中,若當前狀態(tài)是減速片段,則下個狀態(tài)為加速片段的概率是53%,為勻速片段的概率是38%,為怠速片段的概率是9%。類3工況常為高速工況,通常沒有怠速,故只有3個狀態(tài)。

表3 片段轉(zhuǎn)移概率 %

3.4 候選工況的構(gòu)建

構(gòu)建行駛工況的基本思路是截取若干段實測的數(shù)據(jù)片段,拼成一個完整的速度—時間歷程。隨機選擇一個小于100s的片段作為開始片段,即初始狀態(tài)。初始狀態(tài)的起始速度可以為0。下個狀態(tài)的選擇,即片段的選擇一般根據(jù)速度、發(fā)動機比功率和轉(zhuǎn)移概率來選擇。依照以下步驟構(gòu)建候選工況：

(1)初始狀態(tài)一旦確定,它的終點速度和比功率就確定了,假設(shè)為當前狀態(tài)。下個狀態(tài)的初始速度與上個片段最終速度絕對值誤差小于0.5km/h,而且下個狀態(tài)初始發(fā)動機比功率與當前狀態(tài)的終點比功率相比不能超過15k W/t。

(2)根據(jù)轉(zhuǎn)移矩陣和上面提到的速度與比功率原則,選擇下個狀態(tài)。

(3)依據(jù)上述原則,不斷選擇下一個狀態(tài),直到候選工況達到理想的長度為止。

3.5 代表性工況的選取

為了從眾多候選工況中選取合適的代表性行駛工況,本文考慮2個參數(shù),即候選工況與試驗數(shù)據(jù)的速度—加速度頻率差 d 1,4種模型事件(加速、減速、怠速、勻速)時間比例差值的絕對值d2：

式中,fr(a,v)為候選工況速度v和加速度a的比值;p(a,v)為試驗數(shù)據(jù)中速度v和加速度為a的比值;fm為候選工況中第m種模型事件的時間比值;pm為試驗數(shù)據(jù)中第m種模型事件的時間比值。

速度 —時間頻率按照速度間隔5km/h、加速度間隔0.278m/s2進行計算。

綜合考慮 d1與 d2,選擇綜合加權(quán)測試法(composite performance measure,CPM)進行評價,即取 CPM最小的候選工況為代表性行駛工況：

式中,w1為d1的標準偏差的倒數(shù);w2為d2的標準偏差的倒數(shù)。

4 結(jié)果分析

以合肥市5條典型道路為例,采用因子分析和兩步聚類分析方法把汽車行駛工況分成3類,再分別對每類行駛工況按照馬爾科夫方法構(gòu)建代表性行駛工況,并與試驗數(shù)據(jù)進行比較。特征參數(shù)包括平均速度v m、速度標準方差σv、正加速度的平均值av+、負加速度的平均值av-、怠速時間比例p i、加速時間比例 p a、減速時間比例 p d、勻速時間比例p c、正的發(fā)動機比功率平均值P+、負的發(fā)動機比功率平均值P-、發(fā)動機比功率平均值Pa、發(fā)動機比功率方差σvsp。

每類行駛工況所構(gòu)建的代表性行駛工況與試驗數(shù)據(jù)的特征參數(shù)比較見表4?？梢钥闯?低速行駛工況怠速比例高,高速行駛工況怠速比例為0。類1、類2和類3所構(gòu)建的代表性行駛工況與試驗數(shù)據(jù)特征參數(shù)平均誤差分別為3.47%、1.26%和1.28%。表明構(gòu)建的行駛工況與試驗數(shù)據(jù)相符,比較合理。

表4 試驗數(shù)據(jù)與代表性行駛工況的特征參數(shù)比較

為了檢驗各類行駛工況模型事件分布是否與試驗數(shù)據(jù)分布一致,將3類代表性行駛工況的片段轉(zhuǎn)移概率與試驗數(shù)據(jù)的片段轉(zhuǎn)移概率進行χ2檢驗,3類檢驗概率P值都大于0.05,說明構(gòu)建行駛工況與試驗數(shù)據(jù)來自同一樣本。構(gòu)建后的代表性行駛工況的χ2檢驗統(tǒng)計值、d 1、d 2和CPM值如表5所示。

表5 代表性行駛工況周期、評價參數(shù)和χ2分布檢驗

類1、類2、類3三類工況構(gòu)建的代表性行駛工況分別如圖3、圖4和圖5所示。

圖3 類1代表性行駛工況

圖4 類2代表性行駛工況

圖5 類3代表性行駛工況

5 結(jié)論

(1)本文先采用兩步聚類分析方法把行駛工況分成低速、中速、高速3類行駛工況,然后采用馬爾科夫方法構(gòu)建行駛工況。構(gòu)建行駛工況時考慮了與燃油消耗相關(guān)的因素(國內(nèi)發(fā)動機工作狀態(tài)分布),因此構(gòu)成的代表性行駛工況能較好地估計汽車的實際燃油消耗。

(2)利用新的分類準則將行駛片段聯(lián)系為特殊的模型事件來定義行駛變量,將速度—時間曲線劃分為加速、減速、勻速和怠速片段,根據(jù)速度、比功率和轉(zhuǎn)移矩陣求取各片段間的轉(zhuǎn)移概率,建立了新的特征參數(shù)評價準則和代表性行駛工況的選取原則。

(3)利用組合方法構(gòu)建的代表性行駛工況與試驗數(shù)據(jù)的特征參數(shù)相比,平均相對誤差小,即構(gòu)建的代表性工況能準確地反映城市道路上汽車的實際行駛工況。

[1] 呂林,鎖國濤.不同路況下汽車行駛工況的測試與分析[J].武漢理工大學學報,2006,30(5)：781-784.

[2] 張建偉,李孟良,艾國和,等.車輛行駛工況與特征的研究[J].汽車工程,2005,27(2)：220-224.

[3] 杜青,楊延相,鄭偉,等.機動車實際道路排放特性及若干影響因素的研究[J].內(nèi)燃機學報,2002,20(4)：297-302.

[4] Lin J,Niemeier D A.An Exploratory Analysis Comparing a Stochastic Driving Cycle to California's Regulatory Cycle[J].Atmospheric Environment,2002,36(38)：5759-5770.

[5] Eisinger D S,Niemeier D A,Stoechenius T,et al.Collecting Driving Data to Support Mobile Source Emissions Estimation[J].Journal of Transportation Engineering,2006,132(11)：845-854.

[6] 姚志良,賀克斌,王岐東,等.IVE機動車排放模型應(yīng)用研究[J].環(huán)境科學,2006,27(10)：1928-1933.

[7] John K,Mitch C,Megan B,et al.MOVES2004 Software Design Reference Manual Draft[M].Ann Arbor：Environmental Protection Agency,2004.

[8] 陳超,鄒瀅.SPASS 15.0中文版常用功能與應(yīng)用實例精簡[M].北京：電子工業(yè)出版社,2009.

[9] 張濟民,王開文,池茂儒.基于AIC信息準則法的擺式列車傾擺伺服系統(tǒng)建模研究[J].中國鐵道科學,2003,24(6)：6-9.

[10] Zhen Dai,Niemeier D A,Eisinger D.Driving Cycles：a New Cycle-building Method that Better Represents Real-world Emissions.Task Order No.66[R].Berkeley：Department of Civil&Environmental Engineering University of California,2008.

[11] Lee T C,Judge G G,Zellner A.Estimating the Parameters of the Markov Probability Model from Aggregate Time Series Data[M].Amsterdam：North-Holland Publishing Company,1970.