基于高斯噪聲的HLSM穩健性多變點檢驗

朱小雪，施三支，銀炳皓

（長春理工大學 數學與統計學院，長春 130022）

變點檢驗是時間序列分析的重要組成部分，時間序列數據的均值、方差或模型系數在未知的時刻突然發生了改變，可能標志著數據生成的重大變化，因此變點檢驗必不可少。變點檢驗問題第一次被提出可追溯到20世紀50年代［1］，現在已經被拓展應用于股票、生物、基因數據、網絡流量分析、交通等領域，并且已經有著廣泛的研究，變點檢驗方面已經有著大量文獻，如Kim（2019）［2］、Cheng（2015）［3］、Liu（2020）［4］、Kurt（2018）［5］、王紅玉（2018）［6］等。

對于變點檢驗，在早期工業質量控制上，Page（1954）［1］基于 CUSUM 檢驗了觀測值的累積和超過閾值的位置。Scott（1974）［7］提出一種二進制分割（Binary Segmentation）算法用于近似多變點檢驗，用成本函數最小化將序列分割成互不相交的段。Killick（2012）［8］提出 PELT算法，將成本函數最小化來求解變點的方法加以拓展，對于潛在的變點進行了修剪。Maidstone（2018）［9］提出一種FPOP算法，通過成本函數最小化求解序列中多變點的精確值。

穩健性是變點檢驗的一個重要方面，近幾年也有文章研究變點的穩健性。如Fearnhead（2019）［10］提出一個動態算法能最優分割序列，選擇對異常點不敏感的損失函數，解決當序列中存在異常時變點被高估的問題。Dehling（2019）［11］基于兩樣本的Hodges-Lehmann檢驗證明當噪聲是重尾分布的情況下有著好的檢驗效果。許多變點檢驗算法，在高斯噪聲異常值存在的情況下檢驗均值變化缺乏魯棒性。基于似然比檢驗的方法（Yau，2015）［12］、基于懲罰似然的方法（Killick，2012）［8］，這類方法對重尾噪聲表現出一定的魯棒性，但在實際應用中，檢驗出的變點數會高于實際變點。

本文提出了一種基于高斯噪聲Hodges-Lehmann Scan方法（HLSM），可在全局序列中推斷多變點的問題，通過窗口滑動簡化成在局部窗口中檢驗單變點的問題。通過數值模擬實驗與PELT和BS方法進行了比較，結果表明HLSM方法與帶有準則的HLSM方法具有較高的穩健性。

1 HLSM方法

1.1 模型與基本假設

假設觀測序列{Xt}t=1,...,n={X1，…，Xn}，將自回歸（AR）過程分為m+1段。對于j=1，…，m，第j個變點τj表示自回歸過程中第j段突變到第j+1段的位置，假設τ0=0和τm+1=n，并且τ0<τ1< … <τm<τm+1。序列的第j段表示為：

其中，μt是每段序列中未知的常數；Yt,j是均值為零的自回歸過程，并滿足：

其中，i=1,…,m，假設每一段都是獨立的時間序列，將多個假設檢驗問題視為檢驗兩個相鄰段 (Xτi-1+1,Xτi-1+2,…,Xτi)和 (Xτi+1,Xτi+2,…,Xτi+1)的μi是否相等。本文采用了窗口滑動方法結合Hodges-Lehmann檢驗來解決假設檢驗問題（4），稱為Hodges-Lehmann Scan方法。

1.2 Hodges-Lehmann Scan

在1.1節中，未定義每段階數pi取值，假設每段的階數最大值為整數pmax。通過Hodges-Lehmann檢驗的掃描統計量估計出序列中潛在的變點。假設窗口半徑h=h(n)取決于樣本量n。定義在t處的窗口和窗口中的觀測值分別為：

為了在窗口中檢驗變點，選擇Dehling（2019）［11］中采用的Hodges-Lehmann檢驗統計量：

通過Mn(t)檢驗統計量，在窗口中能夠得到一組Hodges-Lehmann檢驗統計量的值(Mn(h)，Mn(h+1)，…，Mn(n-h))。若t處是變點，則Mn(t)值是在局部窗口中的最大值。因此通過HLSM得出一組潛在的變點值，定義如下：

1.3 模型選擇

時間序列中常用最小化準則函數估計多變點，在這里選擇 Davis（2006）［13］提出的最小描述長度（minimum description length；MDL）準則。

給定一組樣本z=(z1,…,zn)，定義準似然函數為：

給出最小描述長度的定義為：

1.4 HLSM算法步驟

在算法1中給出了HLSM算法的具體步驟：

步驟（一）：輸入一系列觀測值：X1,…,Xn，Xi∈R，設置窗口半徑為h=dlog(n)2；

2 模擬實驗

在本節中，比較了PELT、BS（Binary Segmentation）和HLSM的兩種方法，即HLSM和加最小化描述長度準則（MDL）的HLSM，模擬數據由式（1）生成，給出了在同方差與異方差的高斯噪聲下，無異常值和有異常值等不同場景下參數變化的實驗結果。

2.1 模擬基礎設置

為了評估變點估計的性能，根據Killick（2012）［8］中提出的真陽率（True positive rate，TPR）和假陽率（False positive rate，FPR）評估每個方法檢驗變點的能力。對于每個數值實驗，定義：

其中，l表示算法檢測到的總變點數；l1表示l中檢測正確的變點數，真實的變點是τi∈Τ，若方法檢驗到的變點τ?i是距離真實變點τi的10個點的范圍內，同樣認為檢驗到的該變點為真實點。

2.2 模擬實驗

對于均值結構的變化，設置漂移的均值μi=Ui(a,a+0.5)，服從獨立同分布的均勻分布，設置a=1.5是低信號水平，a=2.5是高信號水平。

圖1給出了低信號場景下（a=1.5），不同算法和不同樣本量下估計的真陽率和假陽率（TPR和FPR）。圖1中的場景A，從左圖的TPR折線圖中看出，除了樣本量在500時，PELT較大以外，其他樣本量都是HLSM方法更好，結合右圖FPR折線圖，當樣本量為500時，PELT的FPR是最大的，能判斷PELT算法有較大的TPR是因為高估了變點值。在場景B中，對于額外的方差變化，從TPR折線圖中看出，HLSM也很穩健，另外兩種經典的算法PELT和BS，在FPR折線圖顯示出，不同樣本量之下，FPR都較高，由此得知這兩種方法都會高估變點。對于加了異常值的場景C下，HLSM方法的TPR在各樣本量下，保持最高，能夠精確估計變點位置，此外應用了模型準則的HLSM方法的FPR在場景A、B、C下都為最優。綜合看，HLSM的兩種方法有著不錯的變點檢驗效果，兩種方法的側重優點不同，其中HLSM方法有較高的真陽率（TPR），而加了模型準則的HLSM方法則有著更小的假陽率（FPR）。但是相比于PELT和BS方法而言，本文提出的兩種方法在真陽率有著較好檢驗效果的同時，能夠保持更小的FPR。

圖1 低信號水平下：真陽率與假陽率折線圖

圖2給出的是在高水平信號（a=2.5）下，樣本量為250的場景C中，四種方法實驗模擬得出的一次結果，其中左邊虛線是真實變點的位置，右邊虛線是方法估計出的變點位置。從圖中觀測出，PELT和BS方法誤判異常值的位置為變點，導致高估變點的數目。但HLSM的兩種方法在有著異常值的場景下，忽略異常值的干擾，不會高估變點，對比PELT和BS方法有著較高的穩健性。

圖2 變點估計位置圖

表1與表2分別是在場景A下，樣本量為1 000時，四種方法測得的真陽率與假陽率。從表中看出，在低水平信號下，HLSM方法的真陽率最高，說明HLSM測得的變點最準確，PELT與BS算法檢驗效果相差無幾，而就假陽率來看，加了準則的HLSM方法最小，說明該方法最穩健，測得的變點幾乎都是正確的。在高水平信號下，BS和PELT真陽率較高，但是這兩種算法的假陽率也高，說明會把不是變點的值誤判成變點。對比低水平信號與高水平信號，變化越大，假陽率的值都有所上升，并且除了BS方法以外，其他三種方法的假陽率也有所降低，說明變化越大，越容易檢測出變點值。

表1 場景A下n=1 000時四種方法的真陽率/%

表2 場景A下n=1 000時四種方法的假陽率/%

表3與表4是在場景B下，樣本量為1 000時四種方法的真陽率與假陽率的數值模擬實驗結果。對于加入額外方差變化的場景，在低水平信號下，HLSM有著最高的TPR，且假陽率也較低。在高水平信號下，PELT方法的TPR最高，但是假陽率也較高，HLSM的兩種方法的真陽率較高，且假陽率很低，綜合表現最好，PELT和BS方法雖然有著較高的真陽率，但是假陽率同樣很高。

表3 場景B下n=1 000時四種方法的真陽率/%

表4 場景B下n=1 000時四種方法的假陽率/%

表5和表6給出的分別是樣本量為1 000時，在加有異常值的場景C中，四種方法真陽率和假陽率的值。從表中能觀察到，同一種方法下，高水平比低水平有著更高的TPR，說明變點較大更易檢測到變點，且除了PELT方法，其他三種方法的FPR也降低。表5中，HLSM有著最高的TPR，對比PELT和BS能更準確地檢測出變點位置。在表6中，加了準則的HLSM方法的FPR最小。結合TPR和FPR分析，HLSM方法檢驗效果更好，BS和PELT方法會高估變點。

表5 場景C下n=1 000時四種方法的真陽率/%

表6 場景C下n=1 000時四種方法的假陽率/%

圖3比較了在高水平下（a=2.5），三種場景中不同的方法得到的TPR和FPR折線圖，場景A中，左TPR圖中觀察出，除樣本量為1 000下，HLSM比經典的PELT和BS方法有著更高的TPR，在樣本量為1 000下，結合FPR觀測，PELT的FPR較大，推測PELT的FPR較大是由于高估變點導致。且PELT和BS方法得到每個樣本量下FPR的值變化較大，而HLSM兩種方法則較平穩。場景B中，同樣也是在樣本量為1 000時，PELT的FPR最大，但是有著較高的FPR，表示在有額外的方法變化時，PELT和BS算法會高估變點。場景C是加入了異常值的情況，結合TPR和FPR觀測，HLSM提出的兩種方法檢測表現最好，有較高的TPR且FPR較低，不僅能精確地估計變點位置，還不會高估變點值，HLSM的兩種方法在異常值存在的情況下，變點檢驗有著穩健性。在三種場景下，綜合TPR和FPR來看，還是HLSM的兩種方法有著較優的檢驗效果。從TPR來看，在高水平下，對比低水平信號和高水平信號的TPR看出，變化越大，就更容易檢驗出變點值。

圖3 高信號水平下：真陽率與假陽率折線圖

3 應用

實際數據來自百度地圖開放平臺2020年8月11日（星期二）長春市生態大街的交通流量數據，數據結構為每五分鐘記錄一次該路段的車輛速度。將提出的兩種HLSM方法應用于交通流數據。

圖4（a）為HLSM檢驗出的變點位置與數量，數量為5，位置分別在41、80、144、180、232。圖4（b）為加MDL準則的HLSM檢驗的變點，數量為3，位置分別在80、144、180。可以看出兩種方法都檢驗出變點80、144和180，對應的時間點為7：45、13：00和16：00，而HLSM 方法則還認為 4：25和21：30也為變點。從檢驗結果可以看出交通變化的特征，車輛速度在7：45的前后，速度差值較大，很可能因為8月11日是工作日，人們上班早高峰發生了交通擁堵。在7：45和13：00這段時間的車輛速度一直較低，早晨人們上班時間不同，無論是出行旅游還是上班一般都是傾向于上午的時間段，所以上午出行很可能會遇上擁堵。

圖4 變點檢驗位置圖

在16：00以后，基本為人們下班的時間，下班時間不同，所以能從圖中看出，車輛速度較低的狀態也是持續了一段時間。對比于上午的早高峰，下午的晚高峰造成的擁堵時間更短。而在21：30之后，基本沒有大量的車輛出行，路況良好，車輛也就暢行無阻。

以上的實例分析，證明了所提出的方法在實際應用中的合理性和較強的實用性。

4 結論

本文提出一種基于窗口滑動和Hodges-Lehm‐ann檢驗的方法——Hodges-Lehmann Scan，識別高斯噪聲下的變點。通過數值模擬實驗表明，在同方差與異方差的高斯噪聲下，認為所提出的新方法能夠精確地估計出變點的數量和位置，并且對于有異常值的場景，對比PELT與BS方法顯示更好的穩健性。對于在實際問題中出現異常值干擾的情況，HLSM避免一般方法會高估變點的現象。將提出的方法應用于交通數據的實例分析，說明該方法的實用性。