楊洋 胡嫣然 袁振洲 王云鵬

（1. 北京航空航天大學(xué) 交通科學(xué)與工程學(xué)院，北京 100191；2. 北京航空航天大學(xué) 車路協(xié)同與安全控制北京市重點(diǎn)實(shí)驗(yàn)室，北京 100191；3. 北京交通大學(xué) 交通運(yùn)輸學(xué)院，北京 100044）

中共中央、國(guó)務(wù)院印發(fā)實(shí)施的《交通強(qiáng)國(guó)建設(shè)綱要》將“安全”和“高效”列入總體要求［1］，高速公路作為國(guó)家綜合立體交通網(wǎng)的重要組成部分，其運(yùn)行安全問題一直廣受關(guān)注。高速公路交通事故造成的死傷人數(shù)和財(cái)產(chǎn)損失往往遠(yuǎn)大于發(fā)生在其他類型道路中的交通事故［2］。此外，高速公路的交通事故還是造成節(jié)點(diǎn)癱瘓和擁堵蔓延的重要誘因［3］，據(jù)統(tǒng)計(jì)，高速公路中超過50%的非常規(guī)擁堵都是由交通事故引發(fā)［4］，嚴(yán)重影響行車效率。因此，挖掘高速公路交通事故影響機(jī)理，探究其時(shí)空傳播效應(yīng)，對(duì)于降低交通事故對(duì)行車效率的影響，并針對(duì)性地實(shí)施事故后的動(dòng)態(tài)交通控制具有重要的理論和現(xiàn)實(shí)意義。

針對(duì)高速公路事故時(shí)空影響傳播的研究，主要是采用數(shù)學(xué)模型、仿真實(shí)驗(yàn)等手段，描述事故后車輛的排隊(duì)和消散現(xiàn)象，量化事故造成延誤的時(shí)空范圍，并分析擁堵條件下的高速公路交通流運(yùn)行特性等。Chow［5］早年便采用沖擊波分析和排隊(duì)論兩種方法計(jì)算了交通事故發(fā)生后排隊(duì)車輛的消散時(shí)間和駕駛?cè)搜舆t時(shí)間。Morales［6］基于車輛守恒的基本概念分析了排隊(duì)和消散機(jī)制，并將這些分析結(jié)果應(yīng)用于兩種傳統(tǒng)的宏觀分析方法——確定性排隊(duì)分析和沖擊波分析。王建軍等［7］以流體力學(xué)車流波理論為基礎(chǔ)，分析道路交通事故下的沖擊波和對(duì)道路交通事故進(jìn)行干涉下的沖擊波，運(yùn)用速度—交通流密度的線性關(guān)系建立模型，以確定交通事故后的最大排隊(duì)距離、對(duì)事故進(jìn)行干涉下的最大排隊(duì)距離以及取消對(duì)事故干涉的時(shí)間點(diǎn)等參數(shù)。胡曉健等［8］首先確定了交通事故發(fā)生地周邊能分流的區(qū)域，在此基礎(chǔ)上進(jìn)行OD（起訖點(diǎn)）反推和交通流再分配，以事故干預(yù)后的區(qū)域內(nèi)行程時(shí)間為標(biāo)準(zhǔn)來確定道路交通事故的影響范圍，并通過實(shí)例驗(yàn)證了該方法的可行性，為事故后的交通疏導(dǎo)提供了一定理論指導(dǎo)。Wang 等［9］應(yīng)用改進(jìn)的確定性排隊(duì)理論量化不同交通事故類別產(chǎn)生的交通延誤，實(shí)驗(yàn)結(jié)果表明，該算法能夠較好地估計(jì)事故引起的延遲，并能捕獲事故持續(xù)時(shí)間內(nèi)高速公路交通流的演化。Chien 等［10］在考慮各種幾何條件和時(shí)變交通分布的情況下，利用CORSIM 算法和確定性排隊(duì)模型的概念，整合有限的仿真數(shù)據(jù)，提出了一種近似延遲的方法。Chung 等［11］開發(fā)了一種量化城市高速公路事故造成的延誤的方法，并確定了此類事故造成的總延誤的致因，該方法的基本思想是通過解決與時(shí)空區(qū)域拓?fù)湟恢碌腂IP（二進(jìn)制整數(shù)規(guī)劃）問題來估計(jì)由事故引起的擁堵區(qū)域最可能的時(shí)空范圍。Sun 等［12］提出了一種基于數(shù)據(jù)驅(qū)動(dòng)的事故影響評(píng)估方法，該方法可以根據(jù)處理后的數(shù)據(jù)識(shí)別每個(gè)交叉口的流量特征，然后通過異常值檢測(cè)量化事故影響，也可以通過異質(zhì)性模型來評(píng)估交通事故影響的持續(xù)時(shí)間；所開發(fā)的程序有助于捕捉城市交叉口或附近的事故影響持續(xù)時(shí)間，并識(shí)別出影響該持續(xù)時(shí)間的主要因素，包括事故特征、道路環(huán)境和附近交叉口的時(shí)間特征。Wang 等［13］開發(fā)了一種具有新穎約束的整數(shù)規(guī)劃模型來估計(jì)道路交通事故的時(shí)空影響；模型的輸入包括給定道路上的歷史速度，以及已知事故的位置和開始時(shí)間，模型的輸出為該事故影響的時(shí)空區(qū)域。該模型得到的結(jié)果與沖擊波的傳播理論是吻合的，其結(jié)果表明，該算法平均可減少95%~98%的計(jì)算時(shí)間。代萱［14］基于Van Aerde 模型，建立了交通流運(yùn)行的交通波方程，研究中考慮了交通波可能會(huì)發(fā)生突變這一情況，將事故影響區(qū)域分為了保護(hù)區(qū)、控制區(qū)和緩沖區(qū)，并分析了各區(qū)域的通行能力，之后利用VISSIM 仿真工具驗(yàn)證了所提方法的有效性和合理性。針對(duì)高速公路的實(shí)時(shí)事故風(fēng)險(xiǎn)，Yang 等［15］通過選取交通流基礎(chǔ)信息、交通流中車隊(duì)的變道行為、交通流中車輛的運(yùn)行狀態(tài)、交通流變量短時(shí)間內(nèi)的突變、車輛跟馳行為等反映交通流動(dòng)態(tài)特征的5 個(gè)維度指標(biāo)，構(gòu)建了風(fēng)險(xiǎn)識(shí)別模型，明確了不同交通狀態(tài)下的實(shí)時(shí)事故風(fēng)險(xiǎn)因子。

綜上，在高速公路事故時(shí)空影響傳播分析方面，目前雖取得了較為豐富的成果，但大部分既有研究是基于仿真結(jié)果且效率較低，適用性不理想；此外，在事故影響傳播的度量上，多采用交通波理論計(jì)算排隊(duì)長(zhǎng)度等指標(biāo)，角度較為單一，且對(duì)交通事故時(shí)空影響在上游傳播的具體量化工作中，鮮有較為合理的指標(biāo)被提出。有鑒于此，文中基于待研究路段的高速公路上微波雷達(dá)交通信息檢測(cè)器采集的交通流數(shù)據(jù)以及對(duì)應(yīng)路段的交通事故數(shù)據(jù)，引入速度變化率作為度量指標(biāo)來判斷某地點(diǎn)在某時(shí)刻是否受到交通事故的影響；應(yīng)用雙線性插值法對(duì)速度變化率在時(shí)空范圍進(jìn)行補(bǔ)全，得到了速度變化率輪廓圖；然后根據(jù)速度變化率閾值，基于Savitzky-Golay 濾波擬合法對(duì)交通事故影響時(shí)空區(qū)域的外圍輪廓點(diǎn)進(jìn)行擬合，從而得到了事故影響時(shí)空區(qū)域的外圍輪廓；最后計(jì)算和分析了不同速度變化率閾值下事故時(shí)空影響的各量化指標(biāo)。

1 數(shù)據(jù)準(zhǔn)備

1.1 研究區(qū)域

文中選取我國(guó)北京—哈爾濱高速公路（簡(jiǎn)稱京哈高速公路）的里程樁號(hào)k0-k39（四方橋至西集）全長(zhǎng)約39 km的路段作為研究路段，圖1給出了待研究路段的區(qū)位情況。研究路段內(nèi)可以收集到2013 年1月至2014年10月的交通事故數(shù)據(jù)（包含事故類型、嚴(yán)重程度、時(shí)間、地點(diǎn)）、交通流數(shù)據(jù)（包含1 min集計(jì)的海量高精度交通流量、速度、占有率以及85%速度等）。下面對(duì)兩組數(shù)據(jù)源特征進(jìn)行介紹。

圖1 文中研究區(qū)域Fig.1 Study area

1.2 交通流數(shù)據(jù)

待研究路段的上下游方向共布設(shè)20 組微波雷達(dá)檢測(cè)器，所采集的原始交通信息數(shù)據(jù)主要包括速度、流量、車道占有率等。上下游微波檢測(cè)器之間的最遠(yuǎn)距離達(dá)6.18 km，最短距離為0.8 km，平均距離約為1.9 km。微波檢測(cè)器主要采集的數(shù)據(jù)包括方向、各個(gè)車道的流量、車道占有率和速度等，采集時(shí)間間隔為1 min。以表1第1行數(shù)據(jù)為例，其含義為：在時(shí)間段2013-04-29 19：29：00—2013-04-29 19：29：59 內(nèi)，在方向2（出京方向），編號(hào)為523050003 的檢測(cè)器檢測(cè)到車道1 的流量為22 輛，速度為59 km/h，占有率為16%。

表1 交通流數(shù)據(jù)示例Table 1 Examples of traffic flow data

微波檢測(cè)器在實(shí)際工作過程中可能會(huì)因?yàn)橛挟愇镎趽酢⑼ㄓ?、惡劣天氣等問題而造成所采集數(shù)據(jù)的缺失或者異常，需要進(jìn)行數(shù)據(jù)清洗。以交通流三參數(shù)（流量、速度、占有率）的取值范圍以及它們之間的基本關(guān)系作為判斷依據(jù)，可篩選出異常數(shù)據(jù)［15］。表2給出了異常交通流數(shù)據(jù)的判別規(guī)則，如果交通流數(shù)據(jù)符合表中的一條或幾條規(guī)則，便將其做直接刪除處理。

表2 異常交通流數(shù)據(jù)的判別規(guī)則Table 2 Discriminant rules of invalid traffic flow data

1.3 交通事故數(shù)據(jù)

采集到研究路段內(nèi)2013年1月—2014年9月的交通事故數(shù)據(jù)共198條，數(shù)據(jù)字段包括事故發(fā)生時(shí)間、樁號(hào)、方向、事故類型、事故描述等。刪除由車輛起火、損壞等原因引發(fā)的單車事件，并對(duì)事故數(shù)據(jù)按照發(fā)生時(shí)間的先后順序進(jìn)行編號(hào)操作，以方便下一步的處理。經(jīng)過數(shù)據(jù)預(yù)處理后，共提取研究路段和時(shí)間段內(nèi)的交通事故共164起。交通事故發(fā)生的地點(diǎn)與交通流信息采集設(shè)備的空間關(guān)系如圖2所示。

圖2 交通事故與交通流檢測(cè)器的空間關(guān)系Fig.2 Spatial relationship between traffic crash and traffic flow detectors

2 模型構(gòu)建

2.1 相關(guān)理論方法

2.1.1 雙線性插值

在大部分實(shí)際問題中，y=f（x）這個(gè)函數(shù)經(jīng)常被用來描述某種潛在的數(shù)量關(guān)系，這種函數(shù)關(guān)系大部分是通過實(shí)驗(yàn)或者觀察測(cè)量所得。盡管在x∈［a，b］這個(gè)區(qū)間內(nèi)，f（x）是客觀存在的，但往往只能給出部分x對(duì)應(yīng)的f（x）值，不能得到［a，b］這個(gè)區(qū)間內(nèi)任意x對(duì)應(yīng)的函數(shù)值。并且，有時(shí)候雖然可以求得y=f（x）的表達(dá)式，卻鑒于計(jì)算復(fù)雜等因素，在使用時(shí)并不方便。在一些問題中，有時(shí)需要得到?jīng)]有給定的那部分x對(duì)應(yīng)的f（x）值，這時(shí)候就需要一個(gè)既和f（x）的特性相符合、計(jì)算又比較簡(jiǎn)單的函數(shù)P（x）來對(duì)f（x）進(jìn)行近似，求解所需要的f（x）值。一般情況下，P（x）為代數(shù)多項(xiàng)式或者分段代數(shù)多項(xiàng)式，同時(shí)還需要滿足P（xi） =f（xi）（i=1，2，3，…，n）。這樣求得的P（x）被稱為插值函數(shù)。文中所使用的插值法具體描述如下：

設(shè)函數(shù)y=f（x）在區(qū)間［a，b］內(nèi)有定義，且已知在點(diǎn)a≤x0

成立，那么，就稱P（x）為f（x）的插值函數(shù)，點(diǎn)x0，x1，…，xn為插值節(jié)點(diǎn)。線性插值即為P（x）是一次函數(shù)的情況，如圖3所示。

圖3 線性插值示意圖Fig.3 Schematic diagram of linear interpolation

設(shè) 函數(shù)y=f(x) 上兩 點(diǎn)（x0，y0）、（x1，y1），則[x0，x1]區(qū)間內(nèi)某一位置x對(duì)應(yīng)的y值為

在文中討論的問題中，只能獲得微波檢測(cè)器所在位置以1 min 為集計(jì)間隔的速率，相應(yīng)地，也就只能獲得相同位置和時(shí)間的速度變化率。要想獲得更具體時(shí)間和位置的速度變化率，就需要在時(shí)間和距離兩個(gè)維度上進(jìn)行插值。線性插值法廣泛應(yīng)用在信號(hào)處理、數(shù)字圖像和視頻處理等方面，包括一維線性插值法、二維線性插值法等。在一些需要多維的處理場(chǎng)景中，一維線性插值法往往無法滿足要求，此時(shí)可以采用多維線性插值法或非線性插值法。本研究中應(yīng)用的數(shù)據(jù)是一整條高速公路相鄰的多個(gè)交通信息檢測(cè)器所采集的交通流數(shù)據(jù)，具有很強(qiáng)的空間共線性，而非線性插值法無法處理共線性較強(qiáng)的數(shù)據(jù)。雙線性插值法作為數(shù)值分析中的一種插值算法，從數(shù)學(xué)意義而言，是具備兩個(gè)變量插值函數(shù)的線性插值擴(kuò)展，其核心思想是在兩個(gè)方向分別進(jìn)行一次線性插值，因此，文中采用二維線性插值法來實(shí)現(xiàn)事故前后交通流速度變化率的相應(yīng)指標(biāo)的構(gòu)建與計(jì)算。當(dāng)函數(shù)值由兩個(gè)變量確定時(shí)，需要分別在x和y兩個(gè)方向上進(jìn)行線性插值。在本文討論的問題中，x軸為時(shí)間，y軸為到事故發(fā)生地的距離，那么圖3就變換成了文中問題的基本原理，如圖4所示。

圖4 雙線性插值示意圖Fig.4 Schematic diagram of bilinear interpolation

假設(shè)已知函數(shù)f在Q11=（x1，y1）、Q12=（x1，y2）、Q21=（x2，y1）以及Q22=（x2，y2）這4 個(gè)點(diǎn)的值，如果要得到未知函數(shù)f在點(diǎn)P= （x，y）的值，那么，首先需要在x方向（即時(shí)間維度）進(jìn)行兩次線性插值，得到

再在y方向（即空間維度）進(jìn)行一次線性插值，得到

綜合起來就是雙線性插值最后的結(jié)果：

2.1.2 Savitzky-Golay濾波器

Savitzky-Golay 濾波器（S-G 濾波器）是一種在時(shí)域內(nèi)基于局域多項(xiàng)式最小二乘法擬合的濾波方法。它能減弱數(shù)據(jù)質(zhì)量因素對(duì)時(shí)序數(shù)據(jù)的不利影響，既可以對(duì)噪聲進(jìn)行濾除，又可以保持波形的寬度和形狀不發(fā)生變化［16］，這些特性使得S-G濾波器非常適用于對(duì)事故影響時(shí)空區(qū)域的外圍輪廓進(jìn)行擬合。

如果設(shè)置S-G 濾波器窗口的寬度為2m+1，說明在窗口中存在2m+1 個(gè)數(shù)據(jù)。設(shè)這一組數(shù)據(jù)為x［i］，i= -m，-m+1，…，0，…，m-1，m，i，且i的取值為整數(shù)，之后使用n階（n≤ 2m+1）多項(xiàng)式f（i）來對(duì)這組數(shù)據(jù)進(jìn)行擬合，公式如下：

擬合后的數(shù)據(jù)點(diǎn)與原數(shù)據(jù)點(diǎn)之間的殘差平方和可以用下式表示：

想要得到最好的擬合結(jié)果，殘差平方和應(yīng)為最小，即E對(duì)該多項(xiàng)式的系數(shù)bnk的偏導(dǎo)應(yīng)為0：

解得

當(dāng)已知擬合的單邊點(diǎn)數(shù)m、多項(xiàng)式的階次n以及待擬合的數(shù)據(jù)x[i]時(shí)，即可求出該多項(xiàng)式。對(duì)應(yīng)的窗口內(nèi)的中心點(diǎn)估計(jì)值通過擬合出來的多項(xiàng)式來求取，而對(duì)于后面的點(diǎn)，通過不斷地移動(dòng)窗口重復(fù)上述操作即可。利用Python 的scipy 包中的savgol_filter 函數(shù)可以實(shí)現(xiàn)Savitzky-Golay 濾波擬合法。該函數(shù)有兩個(gè)主要參數(shù)：window_length 和k。window_length 為窗口的長(zhǎng)度，取值為奇正整數(shù)；k為擬合窗口內(nèi)數(shù)據(jù)點(diǎn)多項(xiàng)式的階數(shù)，其值應(yīng)小于window_length的值。window_length的值越小，擬合后的曲線越貼近事故影響時(shí)空區(qū)域的外圍輪廓的實(shí)際曲線；window_length值越大，事故影響時(shí)空區(qū)域的外圍輪廓的曲線平滑效果越強(qiáng)。k值越大，擬合后的曲線越貼近事故影響時(shí)空區(qū)域的外圍輪廓的實(shí)際曲線；反之，k值越小，曲線越平滑。

2.2 事故時(shí)空影響傳播分析指標(biāo)的選取

高速公路交通事故的影響是多方面的，在交通流三參數(shù)（速度、占有率和流量）中，若事故前后高速公路上車輛都較少，則流量指標(biāo)不足以辨識(shí)事故前后的變化，那么可以明顯觀測(cè)到的交通流指標(biāo)變化主要是車輛速度和占有率。因?yàn)榭赡苡兄T多其他因素（例如不同路段限速、司機(jī)駕駛習(xí)慣等因）的影響和干擾，故不應(yīng)直接采用速度和占有率作為此類問題的度量指標(biāo)，而速度的二級(jí)衍生指標(biāo)能更好地消除這些因素帶來的影響。與交通流的歷史運(yùn)行特征進(jìn)行比較，可以更直觀的通過相應(yīng)指標(biāo)量化交通事故影響產(chǎn)生的效應(yīng)。在此，定義速度變化率為某地點(diǎn)在某時(shí)刻的速度相較于該時(shí)刻該地點(diǎn)的歷史速度的變化率，以速度變化率作為評(píng)估指標(biāo)來判斷某地點(diǎn)在某時(shí)刻是否受到交通事故的影響，計(jì)算方法如下：

式中，Δv（d，t）為地點(diǎn)d在t時(shí)刻的速度變化率，vp（d，t）為地點(diǎn)d在t時(shí)刻的過往速度均值，vn（d，t）為高速公路交通事故發(fā)生后地點(diǎn)d在t時(shí)刻的速度的值。速度變化率越大，說明速度變化的幅度越大。如果速度變化率超過某一閾值，就判定該地點(diǎn)在該時(shí)刻受到了交通事故的影響。

一般來說，高速公路交通事故發(fā)生地的下游不會(huì)受到交通事故的影響，速度較事故發(fā)生前基本沒有變化，所以只需要在事故發(fā)生地上游選擇交通信息檢測(cè)器即可。文中選擇距離交通事故發(fā)生地最近的上游4個(gè)檢測(cè)器，時(shí)間的選取根據(jù)實(shí)際情況來確定，用以計(jì)算事故發(fā)生后速度變化率的變化情況。

2.3 輪廓圖構(gòu)建

文中使用的交通流數(shù)據(jù)是以1 min 集計(jì)的，并且在地點(diǎn)上也只能得到4個(gè)檢測(cè)器位置的交通流數(shù)據(jù)，無法獲得更具體的時(shí)間和地點(diǎn)的速度變化率。基于此，應(yīng)用雙線性插值法實(shí)現(xiàn)速度變化率在時(shí)間和空間范圍內(nèi)的補(bǔ)全操作，以獲得更為具體的時(shí)間和空間的速度變化率。

使用雙線性插值法對(duì)速度變化率進(jìn)行補(bǔ)全后，可以得到速度變化率輪廓圖。速度變化率輪廓圖不能直接用以計(jì)算事故時(shí)空影響的量化指標(biāo)，需要選擇合適的速度變化率閾值q來判斷某地點(diǎn)在某時(shí)刻是否受到交通事故的影響。速度變化率大于閾值q意味著該地點(diǎn)在該時(shí)刻受到了交通事故的影響，速度變化率小于閾值q則意味著沒有受到事故影響。根據(jù)閾值q來提取受到交通事故影響的時(shí)空區(qū)域，得到交通事故影響時(shí)空區(qū)域圖后，提取輪廓圖外圍的各個(gè)點(diǎn)，然后使用Savitzky-Golay 濾波擬合法對(duì)這些點(diǎn)進(jìn)行擬合，可得交通事故影響時(shí)空區(qū)域的外圍輪廓，根據(jù)事故影響時(shí)空區(qū)域的外圍輪廓圖來計(jì)算各個(gè)事故時(shí)空影響評(píng)價(jià)指標(biāo)。

2.4 事故時(shí)空影響量化指標(biāo)

文中選取3類指標(biāo)對(duì)高速公路交通事故在上游的時(shí)空影響傳播進(jìn)行量化，分別是時(shí)間相關(guān)指標(biāo)、空間相關(guān)指標(biāo)以及速度指標(biāo)。圖5給出了各指標(biāo)的計(jì)算示意。這3類指標(biāo)都是在事故影響時(shí)空區(qū)域的外圍輪廓圖的基礎(chǔ)上進(jìn)行計(jì)算的［17］。下面對(duì)這3類指標(biāo)的計(jì)算方法進(jìn)行介紹。

圖5 各量化指標(biāo)的計(jì)算示意圖Fig.5 Schematic diagram of the calculation of each quantita?tive index

2.4.1 時(shí)間相關(guān)指標(biāo)

事故影響開始時(shí)間：文中定義事故影響開始時(shí)間為t0，即從t0時(shí)刻起，速度變化率開始超過閾值q，并在一段時(shí)間內(nèi)都處于超過閾值q的狀態(tài)。

事故影響結(jié)束時(shí)間：文中定義事故影響結(jié)束時(shí)間為t1，即從t0時(shí)刻起到t1時(shí)刻，速度變化率都保持在閾值q以上；t1時(shí)刻后，速度變化率開始低于閾值q，并且之后都保持在不超過閾值q的狀態(tài)。

事故影響持續(xù)時(shí)間：文中定義事故影響持續(xù)時(shí)間為Δt，計(jì)算方式為事故結(jié)束時(shí)間t1與事故影響開始時(shí)間t0的差值，公式如下：

2.4.2 空間相關(guān)指標(biāo)

事故影響最近距離：文中定義事故影響最近距離為l0，即在事故影響開始的t0時(shí)刻到事故影響結(jié)束的t1時(shí)刻這一時(shí)間范圍內(nèi)，與事故發(fā)生地相距最近的受到事故影響的地點(diǎn)到事故發(fā)生地的距離。

事故影響最遠(yuǎn)距離：文中定義事故影響最遠(yuǎn)距離為l1，即在事故影響開始的t0時(shí)刻到事故影響結(jié)束的t1時(shí)刻這一時(shí)間范圍內(nèi)，與事故發(fā)生地相距最遠(yuǎn)的受到事故影響的地點(diǎn)到事故發(fā)生地的距離。

事故影響空間范圍：文中定義事故影響空間范圍為Δl，計(jì)算方式為事故影響最遠(yuǎn)距離l1與事故影響最近距離l0的差值，即

2.4.3 速度指標(biāo)

文中定義事故影響傳播速度為v1，即在交通事故影響向上游傳播的過程中，單位時(shí)間內(nèi)事故影響距離的變化值。計(jì)算公式為

式中：Δt1為在交通事故向上游傳播的過程中的某一時(shí)間間隔；Δl1為在Δt1這一時(shí)間間隔內(nèi)，交通事故向上游傳播的距離。

定義事影響傳播速度為v2，即在交通事故影響消散的過程中，單位時(shí)間內(nèi)事故影響距離的變化值：

式中：Δt2為在交通事故影響消散的過程中的某一時(shí)間間隔；Δl2為在Δt2這一時(shí)間間隔內(nèi)，交通事故影響距離的變化值。

2.5 分析流程設(shè)計(jì)

將前面的分析匯總，可得到高速公路事故時(shí)空影響動(dòng)態(tài)效應(yīng)傳播分析流程，主要分為5個(gè)步驟。

步驟1選擇高速公路交通事故上游距離事故發(fā)生地最近的4個(gè)交通流檢測(cè)器，時(shí)間范圍根據(jù)實(shí)際事故情況選擇，觀察、分析事故前后交通流檢測(cè)器的速度變化。

步驟2計(jì)算選定的時(shí)間范圍內(nèi)各交通流檢測(cè)器的速度變化率。

步驟3利用雙線性插值法對(duì)速度變化率進(jìn)行填充，得到更具體時(shí)間和空間范圍內(nèi)的速度變化率，繪制速度變化率輪廓圖。

步驟4選擇合適的速度變化率閾值q，根據(jù)閾值q來繪制受到交通事故影響的時(shí)空區(qū)域，并對(duì)其外圍輪廓點(diǎn)進(jìn)行提取；利用Savitzky-Golay 濾波擬合法對(duì)受到交通事故影響的時(shí)空區(qū)域外圍輪廓點(diǎn)進(jìn)行擬合，得到事故影響時(shí)空區(qū)域的外圍輪廓圖。

步驟5根據(jù)事故影響時(shí)空區(qū)域的外圍輪廓圖來計(jì)算各事故時(shí)空影響量化指標(biāo)，包括時(shí)間相關(guān)指標(biāo)、空間相關(guān)指標(biāo)和速度指標(biāo)，并對(duì)結(jié)果進(jìn)行分析。

3 算例分析

以京哈高速公路出京方向里程數(shù)為19.3 km 的張家灣附近發(fā)生的一起交通事故作為實(shí)例分析對(duì)象。該事故的發(fā)生時(shí)間是2013 年4 月18 日18 時(shí)02 分，事故起因?yàn)橐惠v大貨和一輛小客車發(fā)生剮蹭，無人員傷亡，有路損和遺灑，此遺灑影響進(jìn)京方向車輛，車輛行駛緩慢。距離事故發(fā)生地點(diǎn)最近的4個(gè)上游檢測(cè)器編號(hào)分別為523050015、523050014、523050013 和523050012，其對(duì)應(yīng)的距離事故點(diǎn)的里程數(shù)分別為0.300、1.316、2.500、3.700 km，如表3所示。

表3 事故上游檢測(cè)器的位置信息Table 3 Location information of the upstream detectors

3.1 事故前后速度變化率分析

圖6 給出了事故發(fā)生前3.5 h 和事故發(fā)生后4.5 h時(shí)，4個(gè)上游檢測(cè)器以及距離事故地點(diǎn)最近的1 個(gè)下游檢測(cè)器速度值的變化?？梢钥闯觯鹿拾l(fā)生后，離交通事故地點(diǎn)最近的檢測(cè)器1最先檢測(cè)到速度的急劇下降，檢測(cè)器2的速度值隨后也出現(xiàn)下降趨勢(shì)，檢測(cè)器3 和4 處的速度值并沒有出現(xiàn)較為明顯的變化趨勢(shì)。同時(shí)也可以看出，交通事故前后，距離事故地點(diǎn)最近的下游檢測(cè)器的速度沒有明顯變化，說明交通事故對(duì)下游車輛的速度沒有顯著影響。

圖6 事故前后的速度變化Fig.6 Speed change before and after the crash

由于計(jì)算速度變化率所需的歷史速度均值數(shù)據(jù)量太過龐大，提取所有歷史數(shù)據(jù)求取速度均值顯然是不可行的，為此，隨機(jī)抽取15 份相同檢測(cè)器在相同時(shí)刻的歷史速度數(shù)據(jù)來計(jì)算速度均值。例如，要計(jì)算檢測(cè)器1 位置在18：00 的歷史速度均值，就首先隨機(jī)抽取15 個(gè)日期，假設(shè)抽取的日期為2013 年4 月1 日、2013 年1 月3 日、2013 年2 月9日、……，提取檢測(cè)器1位置在這些日期18：00的速度，并計(jì)算其均值，該均值即為檢測(cè)器1位置在18：00 的歷史速度均值。根據(jù)歷史速度均值即可計(jì)算出事故前后的速度變化率，如圖7所示。

圖7 事故前后速度變化率的折線圖Fig.7 Curve chart of speed change rate before and after the crash

3.2 速度變化率輪廓圖的繪制

以時(shí)間為橫軸、距事故發(fā)生地的距離為縱軸，并以顏色的深淺來表示速度變化率的大小，可以繪制出事故前后速度變化率的時(shí)空?qǐng)D，如圖8所示。

圖8 事故前后速度變化率的時(shí)空?qǐng)DFig.8 Spatio-temporal diagram of speed change rate before and after the crash

利用雙線性插值對(duì)速度變化率進(jìn)行時(shí)間和空間范圍內(nèi)的填充。在縱向即空間范圍內(nèi)以1 m 為間隔、在橫向即時(shí)間范圍內(nèi)以10 s為間隔對(duì)速度變化率進(jìn)行填充，可以得到速度變化率輪廓圖，如圖9所示。通過速度變化率輪廓圖可以得到檢測(cè)器1到檢測(cè)器4 的距離范圍和14：35—22：35 時(shí)間范圍內(nèi)任意距離、任意以10 s 為間隔的時(shí)刻的速度變化率。

圖9 速度變化率輪廓圖Fig.9 Contour diagram of speed change rate

需要補(bǔ)充的是，由于文中最近只能獲得距離事故地點(diǎn)300 m 的檢測(cè)器1 的速度變化情況，而無法獲得距離事故地點(diǎn)0~300 m 這一距離范圍內(nèi)的速度變化率（這是插值法的局限性所在，也是此類交通信息采集手段在該類問題中的缺陷），所以文中事故影響的最近距離為檢測(cè)器1的位置，而非事故發(fā)生地。

3.3 事故影響時(shí)空區(qū)域外圍輪廓圖的提取

速度變化率輪廓圖不能直接用以計(jì)算事故時(shí)空影響量化指標(biāo)，需選擇合適的速度變化率閾值q，根據(jù)閾值q來提取受到交通事故影響的時(shí)空區(qū)域，并判斷某地點(diǎn)在某時(shí)刻是否受到交通事故的影響。由圖7可以看出，在交通事故發(fā)生前，速度變化率并非一直保持為零，多有大于零的情況發(fā)生，但總體上講，速度變化率都是小于20%的。此外，其他相關(guān)研究基于不同數(shù)據(jù)集建模［17］并首選20%作為速度變化率閾值時(shí)，也取得了指標(biāo)值量化值差異較為顯著的結(jié)果，所以在此將20%作為速度變化率閾值的起點(diǎn)。文獻(xiàn)［17］在討論此類問題時(shí)將速度變化率閾值終點(diǎn)設(shè)置為0.4，而圖9 的速度變化率輪廓圖顯示，40%所處的顏色為輪廓圖的邊緣位置，所以文中可以將40%作為速度變化率閾值的終點(diǎn)，并依據(jù)此類研究的慣例將10%設(shè)為速度變化率閾值的步長(zhǎng)，以便于分析［17］。綜上，取20%、30%和40%作為速度變化率閾值的取值，討論不同閾值下的交通事故時(shí)空影響。下文以閾值取20%為例對(duì)分析過程進(jìn)行具體闡述。

速度變化率閾值取20%，即如果某地點(diǎn)在某時(shí)刻的速度變化率大于20%，則該地點(diǎn)在該時(shí)刻被判斷為受到了交通事故的影響；如果小于20%，則被判斷為沒有受到交通事故的影響。提取受到交通事故影響的時(shí)空區(qū)域，如圖10所示。

圖10 交通事故影響的時(shí)空區(qū)域Fig.10 Spatio-temporal region affected by traffic crashes

對(duì)交通事故影響時(shí)空區(qū)域的外圍各個(gè)點(diǎn)進(jìn)行提取，可以得到圖11（a）。借助Python 中scipy 包的savgol_filter 函數(shù)可以實(shí)現(xiàn)Savitzky-Golay 濾波擬合法，對(duì)受到交通事故影響的時(shí)空區(qū)域的外圍點(diǎn)進(jìn)行擬合，其中window_length 取71，k值取3，從而可以得到圖11（b）。

圖11 擬合前后的交通事故影響時(shí)空區(qū)域外圍情況Fig.11 Peripheral and the outer contour of the spatio-temporal region affected by the crash before and after fitting

3.4 結(jié)果分析與討論

根據(jù)以上過程，可以得到不同速度變化率閾值（20%、30%、40%）下事故影響時(shí)空區(qū)域的外圍輪廓圖，在此基礎(chǔ)上可以得到事故時(shí)空影響的各個(gè)量化指標(biāo)，包括事故影響開始時(shí)間、事故影響持續(xù)時(shí)間、事故影響最近距離、事故影響最遠(yuǎn)距離和事故影響空間范圍，結(jié)果如表4所示。

表4 不同閾值下的事故時(shí)空影響量化指標(biāo)結(jié)果Table 4 Results of quantification index of spatio-temporal im?pact of crashes under different thresholds

圖12（a）給出了速度變化率的閾值取20%時(shí)的結(jié)果。如果事故影響傳播速度為正，說明事故影響在向上游傳播；如果事故影響傳播速度為負(fù)，說明事故影響正在消散。在18：05：10時(shí)刻，也就是事故發(fā)生3.17 min后，事故影響開始傳播到交通事故地點(diǎn)上游第1個(gè)微波雷達(dá)信息檢測(cè)器，之后事故影響的傳播速度開始減弱至接近0.2 m/s。從18：11：50起，事故影響的傳播速度有微弱回升，并開始出現(xiàn)波動(dòng)，但一直處于大于0 m/s 的狀態(tài)，說明事故影響一直在向上游進(jìn)行傳播。在18：43：20—18：47：30這個(gè)時(shí)間段內(nèi)，事故影響的傳播速度下跌至0 m/s以下，事故影響出現(xiàn)短暫消散傾向。之后，事故影響傳播速度在［-1，1］m/s 之間上下波動(dòng)。21：10：10為集結(jié)波和消散波的交匯點(diǎn)，在此之后，也就是事故發(fā)生后188.17 min，事故影響傳播速度開始保持在小于0 m/s 的狀態(tài)，即事故影響開始消散；直到21：34：00，事故影響基本結(jié)束。

圖12（b）描述了速度變化率的閾值取30%時(shí)的結(jié)果。如果事故影響傳播速度為正，說明事故影響在向上游傳播；如果事故影響傳播速度為負(fù)，說明事故影響正在消散。在18：07：30 時(shí)刻，也就是事故發(fā)生5.50 min后，事故影響開始傳播到交通事故地點(diǎn)上游第1個(gè)微波雷達(dá)信息檢測(cè)器，之后事故影響的傳播速度開始減弱至接近0 m/s。從18：12：30起，事故影響的傳播速度有微弱回升，并穩(wěn)定在0.27 m/s 左右，直到18：31：50。在18：42：30—18：45：00這個(gè)時(shí)間段內(nèi)，事故影響的傳播速度下跌至0 m/s 以下，事故影響出現(xiàn)短暫消散傾向。在這之后，事故影響傳播速度在［-1，1］m/s 之間上下波動(dòng)。21：10：20為集結(jié)波和消散波的交匯點(diǎn)，在此之后，也就是事故發(fā)生188.33 min后，事故影響傳播速度開始保持在小于0 m/s 的狀態(tài)，即事故影響開始消散；直到21：28：10，事故影響基本結(jié)束。

圖12 不同速度變化率閾值下的事故影響傳播速度變化Fig.12 Change of crash impact propagation speed at various speed change rate thresholds

圖12 （c）為速度變化率的閾值取40%時(shí)的結(jié)果。如果事故影響傳播速度為正，說明事故影響在向上游傳播；如果事故影響傳播速度為負(fù)，說明事故影響正在消散。在18：08：50 時(shí)刻，也就是事故發(fā)生6.83 min后，事故影響開始傳播到交通事故地點(diǎn)上游第1個(gè)微波雷達(dá)信息檢測(cè)器，之后事故影響的傳播速度開始下降，在18：12：00—18：13：40、18：24：40—18：26：10、18：32：50—18：34：20 以及18：42：30—18：44：20 這4 個(gè)時(shí)間段內(nèi)，事故影響的傳播速度下跌至0 m/s以下，事故影響出現(xiàn)短暫消散傾向。之后，事故影響傳播速度在［-0.6，0.6］ m/s之間上下波動(dòng)。21：10：20 為集結(jié)波和消散波的交匯點(diǎn)，在此之后，也就是事故發(fā)生188.33 min 后，事故影響傳播速度開始保持在小于0 m/s 的狀態(tài)，即事故影響開始消散；至21：23：30，影響基本結(jié)束。

由表4可以看出，速度變化率的閾值越小，事故影響開始的時(shí)間就越早，事故影響結(jié)束的時(shí)間就越晚，事故影響持續(xù)時(shí)間就越長(zhǎng)，事故影響的距離也越長(zhǎng)。反之，速度變化率閾值越大，事故影響開始的時(shí)間就越晚，事故影響結(jié)束的時(shí)間就越早，事故影響持續(xù)時(shí)間就越短，事故影響的距離也越短。從事故影響傳播速度變化圖來看，不同的速度變化率閾值下，事故影響的傳播速度在時(shí)間維度的發(fā)展態(tài)勢(shì)各有不同，但其產(chǎn)生的集結(jié)波和消散波所匹配得到的交匯點(diǎn)所處的時(shí)間相差并不大。這樣的結(jié)果驗(yàn)證了文中所提出的方法的可靠性。

4 結(jié)論

文中通過對(duì)高速公路交通事故時(shí)空影響動(dòng)態(tài)效應(yīng)的傳播分析，得到如下主要結(jié)論。

1）速度變化率的閾值越小，高速公路交通事故影響開始的時(shí)間就越早，事故影響結(jié)束的時(shí)間就越晚，事故影響的距離就越長(zhǎng)；反之，速度變化率閾值越大，交通事故影響開始的時(shí)間就越晚，事故影響結(jié)束的時(shí)間就越早，事故影響的距離就越短。

2）不同的速度變化率閾值下，高速公路交通事故影響的傳播速度隨時(shí)間的變化趨勢(shì)存在差異，但集結(jié)波和消散波的交匯點(diǎn)所在的時(shí)間相差無幾，這進(jìn)一步驗(yàn)證了文中所提出的方法的可靠性。

3）文中所構(gòu)建的方法計(jì)算操作便捷、結(jié)果分辨率高，在高速公路安全管理工作中，可基于該方法明確交通事故傳播的時(shí)空特征；在交通事故發(fā)生后，可通過主動(dòng)調(diào)控并設(shè)置合適的交通引導(dǎo)信息，來降低二次事故發(fā)生風(fēng)險(xiǎn)并減少擁堵時(shí)間。

文中應(yīng)用了微波雷達(dá)交通信息檢測(cè)器采集的交通流數(shù)據(jù)，其動(dòng)態(tài)交通特征并不全面。后續(xù)研究中，可通過視頻檢測(cè)器等獲得車輛的變道行為和跟馳行為等參數(shù)，或者基于真實(shí)駕駛數(shù)據(jù)采集車輛的加速度、車頭時(shí)距等參數(shù)，完善建模輸入的信息。