基于WBS2的一種多維數(shù)據(jù)多變點檢驗方法及其應(yīng)用

毛佳慧，施三支

（長春理工大學(xué) 理學(xué)院，長春 130022）

變點是指在該時間點上，樣本的分布或者數(shù)字特征在突然發(fā)生變化。研究變點問題可以判斷過程中某參數(shù)發(fā)生變化的時刻并有效控制該參數(shù)，也能夠分析系統(tǒng)的穩(wěn)定性，從外部控制變量出發(fā)檢驗或預(yù)測形態(tài)發(fā)生的變化。變點檢驗就是估計變點的數(shù)量和位置，該研究現(xiàn)已廣泛用于工業(yè)質(zhì)量控制[1]、醫(yī)學(xué)診斷[2]、交通流[3]和網(wǎng)絡(luò)安全[4]等許多領(lǐng)域。檢驗方法也從參數(shù)檢驗發(fā)展到非參數(shù)檢驗，檢驗對象也從一維擴展到多維甚至是高維數(shù)據(jù)。關(guān)于一維數(shù)據(jù)的多變點問題，許多學(xué)者給出了檢驗效果非常好的方法，如2012年Killick[5]考慮大型數(shù)據(jù)集的多變點問題，提出PELT法來找到變點可能的數(shù)量和位置的成本函數(shù)最小值，從而得到具有計算成本的變點的最佳數(shù)量和位置，其計算效率是O()n。Rigaill[6]于2015年提出的修剪動態(tài)規(guī)劃將最佳情況下的一些基于懲罰項的變點估計方法的計算耗時加速到線性時間，并提供快速實現(xiàn)。2014年Frick[7]針對指數(shù)族回歸中的變點提出了控制FWER的SMUCE法。該方法首先最小化α級的一個多尺度檢驗的接受域上的變點數(shù)量，從而估計未知步長函數(shù)，再通過構(gòu)建漸進未知步長函數(shù)和變點的置信集，得到估計變點位置的指數(shù)界限。Li和 Munk[8]在 2016年提出了基于 FDR的相關(guān)方法。Fryzlewicz和 Subba Rao[9]于 2014 年以及 Cho 和 Fryzlewicz[10]于 2012 年利用二元分割對單變量時間序列數(shù)據(jù)進行分割，Cho和Fryzle?wicz[11]于2016年對多變量甚至高維時間序列數(shù)據(jù)進行分割。2007年Venkatraman和Olshen[12]提出的圓形二元分割，2014 年 Fryzlewicz[13]提 出 的WBS和 2018年 Baranowski等人[14]提出的最窄閾值法旨在提高二元分割的性能。2019年Fryzle?wicz[15]在克服了WBS算法缺點的基礎(chǔ)上提出了WBS2（Wild Binary Segmentation2），不同于 WBS中所有子分段都是預(yù)先分割好的，WBS2算法是數(shù)據(jù)自適應(yīng)的，每一組子分段的位置都是由先前檢驗到的變點的位置決定的，并使用SDLL進行模型選擇。關(guān)于高維數(shù)據(jù)變點檢驗問題的可用工具較少，并且大多不具有普適性，2019年Fryzlewicz針對英國32個城區(qū)的房價指數(shù)，利用主成分分析將一個32維、長度為284的數(shù)據(jù)降維成兩組一維數(shù)據(jù)，再分別進行變點檢驗。

通過研究交通流多維數(shù)據(jù)的特點，給出了一種基于線性投影的多維數(shù)據(jù)變點檢驗方法，能夠同時檢驗多條線或者多個站點的突變點，并對其變點位置從發(fā)生時間進行估計。通過生成階梯狀模擬數(shù)據(jù)，利用不同的變點檢驗方法進行對比分析，發(fā)現(xiàn)WBS2.SDLL能得到較好的檢驗結(jié)果，由于檢測的交通數(shù)據(jù)具有高度相關(guān)性，先采用主成分分析和線性投影將多維數(shù)據(jù)降為一維數(shù)據(jù)后，再利用WBS2.SDLL進行變點檢驗，檢驗結(jié)果同每個站點單獨進行多變點檢驗相比，具有快速準確的效果。再高效地找出交通高峰時間段，并推斷其發(fā)生的原因。

1 數(shù)據(jù)降維方法

首先利用主成分分析將分量相關(guān)的原隨機向量通過正交變換轉(zhuǎn)化成分量不相關(guān)的新隨機變量。具體步驟如下：

假設(shè)進行主成分分析的指標變量有m個：x1,x2,…,xm，共有n個評價對象，第i個評價對象的第j個指標的取值為aij。將各指標aij轉(zhuǎn)換成標準化 指 標

其中i=1,2,…,n;j=1,2,…,m；且

為標準化指標變量。再計算相關(guān)系數(shù)矩陣R=(rij)m×m：

式中，rii=1；rij=rji；rij是第i個指標與第j個指標的相關(guān)系數(shù)。再計算相關(guān)系數(shù)矩陣R的特征值λ1≥λ2≥…≥λm≥0，及對應(yīng)的特征向量u1,u2,…,um，其中uj=(u1j,u2j,…,umj)T，由特征向量組成m個新的指標變量：

式中，y1是第1主成分；y2是第2主成分，…，ym是第m主成分。再選擇p(p≤m)個主成分，計算綜合評價值。首先計算特征值λj(j=1,2,…,m)的信息貢獻率和累積貢獻率。主成分yj的信息貢獻率為：

αp為主成分y1,y2,…,yp的累積貢獻率：

當(dāng)αp接近于1(αp=0.85,0.90,0.95 )時，則選擇前p個指標變量y1,y2,…,yp作為p個主成分，代替原來m個指標變量。文章提取出累積貢獻率達到95%的前三個主成分，再把每個主成分的載荷占三個主成分載荷的和的比例值作為這三個主成分的投影系數(shù)，從而能夠通過結(jié)合主成分分析和線性投影的方法將多維數(shù)據(jù)降成一維數(shù)據(jù)。

2 WBS2.SDLL

考慮模型：

其中，ft是一個確定的一維分段持續(xù)信號，它的變點數(shù)量N和位置η1,…,ηN未知。隨機序列是相互獨立且服從均值為0且方差為σ2的正態(tài)分布，序列是有界的，即對于t=1,…,T，有|ft|0時，δT=δT，間隔長度，滿足，故對于一個足夠大的常數(shù)C，δT和B由要求。Fryzlewicz給出了用于識別候選變點位置的主要定位統(tǒng)計量是累積和，即對于數(shù)據(jù)(Xs,…,Xe) ：

其中，s≤b≤e，且n=e-s+1。

Fryzlewicz[15]在 2019 年提出了一種遞歸算法：WBS2。它在第一次遍歷數(shù)據(jù)時，抽取少量M個子區(qū)間樣本，起點和終點隨機選擇，并均勻且獨立地從集合{1,…,T} 中替換，WBS2將M個CUSUM的絕對值的argmax記為第一個候選變點，然后利用此候選變點將集合{1,…,T} 劃分成兩個，并再次遞歸地在候選變點的左側(cè)和右側(cè)抽取M個子區(qū)間樣本，依此類推，在短子域上所有可能子區(qū)間樣本將小于M，在這種情況下，WBS2將繪制每個這樣的子域的所有子區(qū)間樣本。

對于檢驗的置信水平為90%和95%的情況，將兩個算法分別稱為WBS2.90和WBS2.95。對降維后的一維數(shù)據(jù)，利用找到(fs,…,fe)中最有可能的一個變點位置。WBS2僅在子域長度為1時停止。完成該步驟后，WBS2按照CUSUM絕對值大小對候選變點進行降序排列。WBS2方案在計算時是數(shù)據(jù)自適應(yīng)的，每個下一批M子區(qū)間的子域都是由先前檢驗到的候選變點的位置確定，由于該過程在{1,…,T}不斷變短的子域上操作，使得WBS2的計算速度加快，且其只需要設(shè)置M這一個參數(shù)值。

接下來將對模擬數(shù)據(jù)進行對比實驗，再利用WBS2對降維后的真實數(shù)據(jù)進行變點檢驗，相關(guān)算法流程圖如圖1所示。

圖1 算法流程圖

3 模擬實驗與實證分析

先利用多變點檢驗法PELT、SMUCE和WBS作為對比對象對階梯狀模擬數(shù)據(jù)進行檢驗，并從中選出表現(xiàn)最好的方法，再對交通行人流數(shù)據(jù)進行變點檢驗。

3.1 模擬計算

使用兩組階梯狀模擬數(shù)據(jù)，其中一組是1 000個帶有10個不同均值相同方差的服從正態(tài)分布的模擬數(shù)據(jù)，其均值分別是 1，3，5，7.5，8，6，5，5.5，6，4.5，方差取 1，其數(shù)據(jù)量均為 100，模擬次數(shù)分別為100，200，500和1 000次。第二組模擬數(shù)據(jù)與第一組的唯一區(qū)別是，其數(shù)據(jù)量是從50，140，120，75，40，160，125，80，60，150中隨機無放回地抽取。模擬次數(shù)分別是500和1 000。

令N是數(shù)據(jù)集的真實變點數(shù)量，?是使用算法估計變點數(shù)量，則表示變點估計數(shù)量的絕對平均值，它與估計的變點位置的許多精度測量值強烈相關(guān)。令?為一個分段常數(shù)函數(shù)，其每對連續(xù)估計變點之間的值是這對變點界定的區(qū)間內(nèi)的數(shù)據(jù)的平均值。令為模擬中估計ft的均方誤差，作為變點位置估計誤差的度量。選用作為變點估計的準確率的評價標準，T表示模擬次數(shù)，模擬結(jié)果如表1、表2和表3所示。

表1 多維數(shù)據(jù)變點檢驗?zāi)M結(jié)果1

表2 多維數(shù)據(jù)變點檢驗?zāi)M結(jié)果2

表3 各組數(shù)據(jù)量不等時數(shù)據(jù)變點檢驗?zāi)M結(jié)果

表1和表2在第一組模擬數(shù)據(jù)的基礎(chǔ)上得到上述結(jié)果，而表3使用的是第二組模擬數(shù)據(jù)，僅考慮模擬次數(shù)分別為500和1 000的情況。由表中結(jié)果可知，利用WBS2.95做變點檢驗得到的變點數(shù)量和變點位置的偏差都是最小的，而隨著模擬次數(shù)的增加誤差逐漸越來越小。第二組數(shù)據(jù)的檢驗誤差相較于第一組的要大，說明隨機抽取間隔數(shù)使得其中某幾組的數(shù)據(jù)量偏小，會直接影響檢驗的準確率。

為了更直觀地展示檢驗結(jié)果和所用算法的優(yōu)劣性，以隨機產(chǎn)生的其中一組模擬數(shù)據(jù)為例，將數(shù)據(jù)和利用不同算法檢驗得到的變點個數(shù)和位置繪制在圖2中。

圖2 不同算法得到的變點個數(shù)和位置

圖2從上至下使用的變點檢驗方法分別是WBS2.95，WBS2.90，WBS，PELT和 SMUCE 法。實線表示其中一個含有10組服從不同均值相同方差的正態(tài)分布，且數(shù)據(jù)量均為100的模擬數(shù)據(jù)集，豎直的虛線表示算法檢驗出的變點位置。

圖2中使用的模擬數(shù)據(jù)集的真實變點個數(shù)為 9，變點位置分別是 100，200，300，400，500，600，700，800和900。由圖2中的豎直虛線的個數(shù)和位置可知使用的5個方法中只有WBS2.95法檢驗得到的變點個數(shù)為9，其他方法得到的變點個數(shù)均為8，盡管利用五個方法檢驗得到的變點位置與實際變點位置均有一定的偏差，但是WBS2.95法的偏差相對較小。通過圖2與表1和表2的對比結(jié)果，采用WBS2.95法對降維后的地鐵進出閘機口行人流數(shù)據(jù)進行變點檢驗。

3.2 實證分析

選取2015年4月1日上海市地鐵一號線上25個站點（富錦路，友誼西路，寶安公路，共富新村，呼蘭路，通河新村，共康路，彭浦新村，汶水路，上海馬戲城，延長路中山北路，上?；疖囌荆瑵h中路，新闡路，人民廣場，黃陂南路，陜西南路，常熟路，衡山路徐家匯，上海體育館，漕寶路，上海南站，錦江樂園）的閘機進出口的行人流數(shù)據(jù)，每五分鐘作為一個計數(shù)節(jié)點。為保證各維度數(shù)據(jù)量一致，且通過檢驗發(fā)現(xiàn)每天6：30之前和22：00之后地鐵閘機進出口的人數(shù)較少，不影響之后的檢驗結(jié)果，故將行人流數(shù)據(jù)的時間段截取為6：20-22：00，從而得到一個25維長度為188的數(shù)據(jù)集。通過運用主成分分析和線性投影將多維的數(shù)據(jù)降成1維數(shù)據(jù)后，再利用WBS2.SDLL對地鐵閘機口進站和出站的人數(shù)進行變點檢驗，得到行人進出的高峰期，通過算法檢驗得到高峰時間段的行人流均值如表4所示。

表4 多個閘機進出口高峰時間段

由表4可看出在時間段7：00-9：00和17：30-18：00內(nèi)，進站閘機口人流量達到最高峰，在時間段 8：05-9：30和 17：40-18：30內(nèi)，出站閘機口人流量達到最高峰。結(jié)合實際與人們通勤時間相符合。同時，為驗證上述結(jié)果的可信度，接下來將對單個站點閘機進出口行人流進行變點檢驗，選取一號線中的一個換乘站：漢中路站（可換乘12和13線）作為表4中變點檢驗結(jié)果的對比驗證對象，結(jié)果如表5所示。

表5 漢中路站閘機進出口高峰時間段

由表5可看出在時間段8：00-8：30和17：20-18：00內(nèi)，漢口站進站閘機口人流量達到最高峰，在時間段8：30-9：00內(nèi)，該出站閘機口人流量達到最高峰。上述時間段均被涵蓋于多個站點檢驗得出的高峰時間段內(nèi)，說明通過結(jié)合主成分分析和線性投影來將數(shù)據(jù)進行降維，再對降維后的數(shù)據(jù)進行變點檢驗的方法對地鐵行人流數(shù)據(jù)是有效的，檢驗結(jié)果也是可信的。

通過對上海市一號線25個地鐵站的閘機進出口行人數(shù)量的變點檢驗分析研究，根據(jù)算法可以計算出行人流發(fā)生突變的次數(shù)和突變發(fā)生的具體時間，從而給交通管理部門對特定時間段的進出閘機口數(shù)量的合理開放提供相關(guān)數(shù)據(jù)支持。

4 結(jié)論

首先對大量階梯狀的模擬數(shù)據(jù)利用不同的的多變點檢驗算法進行變點估計，得到表現(xiàn)較好的WBS2.95。再利用主成分分析和線性投影將一個25維的上海市一號線地鐵站閘機進口的行人流數(shù)據(jù)有效地降成一維數(shù)據(jù)，再利用WBS2.95進行變點檢驗，結(jié)果表明7：00-8：30和17：20-19：00為進站高峰期，其中 7：46-8：20人流量最大，8：05-9：00和 17：45-20：45為出站高峰期，其中18：15-18：55人流量最大。

將WBS2.SDLL方法運用到行人交通樞紐中，同時考慮多個站點行人流的變點問題，為深入研究多個站點或多條線上的行人流突變問題提供一種新方法，但是WBS2.SDLL只適用于一維數(shù)據(jù)，之后還將嘗試研究一種不需要事先進行數(shù)據(jù)降維就能夠直接用于檢驗多維數(shù)據(jù)的方法，為多維甚至高維數(shù)據(jù)的變點檢驗問題提供參考。