基于歷史預警準確率的時空重排掃描最大掃描半徑優(yōu)化方法

張 亞 楠，龍 華，邵 玉 斌，杜 慶 治，陳 騰 飛

(昆明理工大學信息工程與自動化學院，云南 昆明 650500)

0 引言

隨著大數(shù)據(jù)時代的到來，通過時空數(shù)據(jù)異常探測可發(fā)現(xiàn)諸多重要信息。例如：識別路況中的時空異常數(shù)據(jù)，有助于檢測造成交通擁堵的交通事件[1]；搜尋疾病[2]、犯罪[3]、火災[4]、極端高溫[5]等發(fā)生的熱點地區(qū)，有助于分析事件發(fā)生的規(guī)律并對未來事件提出預警。為利用時空異常數(shù)據(jù)中蘊含的重要信息并挖掘其中的價值，越來越多的學者開始探索時空異常數(shù)據(jù)的探測方法[6-8]。

異常探測可分為事物異常探測、空間異常探測、時空點事件異常探測、時空序列異常探測、時空軌跡異常探測5類。其中，時空點事件包括離群事件和熱點事件：前者指存在于時空域內(nèi)的孤立事件點以及少量事件的聚集；后者指顯著程度較大的局部聚集[9]。Kulldorff于1997年首次提出掃描統(tǒng)計方法[10]并對時空熱點事件進行異常探測，2001年在僅考慮空間屬性的掃描模型中加入時間屬性，提出時空掃描統(tǒng)計方法[11]，2005年進一步提出無需人口數(shù)據(jù)，僅根據(jù)區(qū)域病例數(shù)即可進行時空掃描分析的時空重排掃描統(tǒng)計方法[12]。該方法使用圓柱體掃描窗口對研究區(qū)域進行尺寸限定，但圓柱體掃描窗口對于不規(guī)則大型實際數(shù)據(jù)集的探測存在局限性[13]。因此，相關學者對圓柱體掃描窗口的關鍵參數(shù)之一——掃描形狀進行了優(yōu)化。例如：Duczmal等提出非圓形簇方法[14]，以提升檢測非常不規(guī)則形狀簇的能力；Takahashi等通過將相鄰子區(qū)域組合，生成形狀不規(guī)則的窗口進行掃描，提出一種靈活時空掃描統(tǒng)計方法[15]，并在算法中加入集群尺寸(一般為整個研究區(qū)域面積的10%～15%)控制，可進行較小尺寸異常點的檢測；Kulldorff等提出橢圓形掃描窗口[16]，對潛在聚集區(qū)的形狀進行限定，防止識別出過度不規(guī)則的聚集區(qū)；萬幼等提出一種改進的不規(guī)則形狀時空異常聚類模式挖掘方法[17]，基于時空鄰近單元格構建時空鄰接矩陣，對蟻群最優(yōu)化掃描統(tǒng)計方法進行改進，使其適應時空區(qū)域掃描，有效識別了時空范圍內(nèi)的不規(guī)則形狀異常聚類。然而，目前對圓柱體掃描窗口的另一關鍵參數(shù)——最大掃描半徑的優(yōu)化研究較少。在時空重排掃描方法中，多以50%研究區(qū)面積所對應的圓形半徑作為最大掃描半徑[12]，掃描半徑隨研究區(qū)擴大而增加，但最大掃描半徑過大將導致計算資源浪費。另外，由于未有效區(qū)分和篩選實際異常點，會影響預警的準確率。為此，楊威等提出基于歷史命中率的時空重排掃描方法[18]選取最大搜索半徑，但該方法逐一計算不同最大掃描半徑下的預警命中率，耗時較長，且仍未有效區(qū)分和篩選實際異常點，影響合適的掃描半徑篩選，從而影響預警準確率。鑒于此，本文提出一種基于歷史預警準確率的時空重排掃描最大掃描半徑優(yōu)化方法，在保證預警準確率的同時，能快速選取較小的最大掃描半徑。

1 研究方法與數(shù)據(jù)

1.1 時空重排掃描方法

時空重排掃描方法的基本原理為：假設事件發(fā)生的概率在時空范圍內(nèi)服從泊松分布，在研究區(qū)域內(nèi)劃分掃描區(qū)域并將該區(qū)域內(nèi)事件發(fā)生概率與掃描區(qū)域外事件發(fā)生概率之比作為掃描統(tǒng)計量；不斷擴大掃描范圍并改變位置中心，尋找出整個研究區(qū)域中掃描統(tǒng)計量值較大的區(qū)域。為消除掃描結果的隨機性影響，對尋找出的監(jiān)測點用蒙特卡洛模擬方法，篩選保留可信度較大的地點作為事件發(fā)生的預警點[12]。對研究區(qū)域以街道或其他行政區(qū)分配ID，每個ID對應一個監(jiān)測點；時間節(jié)點(t=1,2,…,T)可選取天、月、年，本文以天為單位。假設研究區(qū)域某時段內(nèi)某時空事件的數(shù)量和為C，其計算公式為[12]：

(1)

式中：Cot表示監(jiān)測點o(o=1,2,…,O)在時間點t監(jiān)測到的事件總數(shù)。

設Co表示監(jiān)測點o在研究時間T內(nèi)發(fā)生的事件總數(shù)，Ct表示所有監(jiān)測點在時間點t發(fā)生的事件總數(shù)，則監(jiān)測點o在時間點t上期望事件數(shù)量μot的計算公式為[12]：

(2)

時空重排掃描過程可模擬為一個圓柱體M的移動過程，設圓柱體M的底面圓心為掃描區(qū)域中一個監(jiān)測點OM，掃描半徑為RM，M的高度為研究的時間范圍。每次掃描指定圓心OM，遍歷所有o=1,2,…,O，o∈M，選出事件地點間距小于RM的監(jiān)測點，并根據(jù)期望事件數(shù)構成期望矩陣μ=[μot],則圓柱體M內(nèi)實際發(fā)生事件數(shù)量CM和期望值μM為[12]：

(3)

(4)

當CM遠小于研究時間范圍內(nèi)發(fā)生的事件總數(shù)C時，可認為CM近似服從泊松分布，其均值為μM[19]。采用泊松模型的掃描統(tǒng)計量對事件聚集區(qū)進行識別和探測，通過判斷事件發(fā)生數(shù)量服從泊松分布的程度，識別事件發(fā)生的獨立性程度。對圓柱體M用廣義似然比GLRM檢驗泊松分布[12，20]：

(5)

本文研究內(nèi)容為熱點事件聚集的監(jiān)測，即只考慮CM>μM的情況。為計算方便，使用對數(shù)化處理后的形式，即對數(shù)廣義似然比(LGLRM)(式(6))。LGLRM越大，說明柱體M內(nèi)事件聚集性越強。通過蒙特卡洛模擬方法對聚集性較高的監(jiān)測點進行顯著性檢驗，計算在對數(shù)似然比下的概率估計P值，最后選擇P值較小的監(jiān)測點作為發(fā)生異常事件的預警點。

(6)

式中：u為指示函數(shù)，在CM>μM情況下，指示函數(shù)為階躍函數(shù)，u=1。

1.2 實際異常點判定及歷史預警準確率

(7)

(8)

利用該方法得到的預警結果中包含的假異常點更少，歷史預警命中率較高，從而可更準確地預警出類似事件發(fā)生可能性更大的監(jiān)測點。本文將歷史預警準確率(α)定義為：在給定空間范圍內(nèi)，對歷史數(shù)據(jù)通過時空重排掃描方法預警到發(fā)生異常事件的監(jiān)測點數(shù)量(β1)與實際發(fā)生異常事件的監(jiān)測點數(shù)量(β2)之比(式(9))，其中β1≤β2，β2>0，β1/β2∈[0,1]；將預警準確率(λ)定義為：在給定空間范圍內(nèi)，對研究數(shù)據(jù)通過時空重排掃描方法預警到發(fā)生事件的監(jiān)測點數(shù)量(γ1)與實際發(fā)生事件的監(jiān)測點數(shù)量(γ2)之比(式(10))，其中γ1≤γ2，γ2>0，γ1/γ2∈[0,1]。

α=β1/β2×100%

(9)

λ=γ1/γ2×100%

(10)

1.3 基于歷史預警準確率的時空重排掃描最大掃描半徑二分選取方法

基于歷史預警準確率的時空重排掃描方法最大掃描半徑的優(yōu)化，就是在搜索區(qū)間內(nèi)尋找使得歷史預警準確率最高的最小半徑。本研究采用一維搜索方法尋找最佳半徑。精確一維搜索常用于求解非線性函數(shù)極值點[21，22]，對函數(shù)的連續(xù)性、可微性沒有嚴格要求，只要求選定的插入點有對應的函數(shù)值即可，因此普適性較好。精確一維搜索算法中的二分法通過取值試探的方式，求解原函數(shù)的導函數(shù)，即通過求解非線性方程的根獲得最優(yōu)解。假定在搜索區(qū)間[a,b]內(nèi)取中點c(c=(a+b)/2)，同理在區(qū)間[a,c]、[c,b]內(nèi)分別取中點d、e(d

本文以歷史預警準確率為目標函數(shù)的原函數(shù)，在選定的掃描半徑區(qū)間內(nèi)是一個非連續(xù)性函數(shù)，不存在導函數(shù)。為得到最優(yōu)半徑或最優(yōu)半徑所在區(qū)間，借鑒二分法思想，不求解目標函數(shù)的導函數(shù)，而采用試探選點的方法,即:確定搜索起始區(qū)間后，選擇一個搜索區(qū)間的中點作為試探點，計算試探點相應的函數(shù)值并進行比較，以確定新的搜索區(qū)間；不斷重復該過程，將區(qū)間縮小至給定搜索區(qū)間的精度范圍，若搜索區(qū)間達到設定的區(qū)間間隔精度，停止二分搜索，逼近最優(yōu)值。綜合不同歷史時段預警結果，選擇半徑長度最小、歷史預警命中率最高的掃描半徑作為最優(yōu)掃描半徑。與文獻[18]確定最大掃描半徑的方法相比，該方法在保證歷史預警準確率的前提下，試探選點比順序選點選取的次數(shù)更少，能快速選取較小的最大掃描半徑，再利用該半徑進行前瞻性的區(qū)域時空重排掃描，可減少因最大掃描半徑過大引起的計算資源損耗。具體步驟如圖1所示。

圖1 最大掃描半徑二分選取流程Fig.1 Flow chart of dichotomy selection of maximum scanning radius

1.4 實驗數(shù)據(jù)

相關研究[24-26]證明，火災事件可以使用時空掃描方法進行時空異常探測。為便于實驗結果的驗證，本文從舊金山地區(qū)數(shù)據(jù)協(xié)調(diào)網(wǎng)站(https://datasf.org/opendata/)提供的“Fire Department Calls for Service”數(shù)據(jù)集中提取2018-2020年的火災事件數(shù)據(jù)進行實驗。由于公共安全事件預警通常為短期預警，因此，本文選取1周作為事件監(jiān)測與預警的時間閾值，以增強事件間的相關性[27]；同時通過耗時驗證方法有效性，即前瞻性時空重排掃描統(tǒng)計分析的運行時間。為減少實驗結果的偶然性和隨機性，本文設計了3組實驗對結果進行驗證;為驗證選擇不同最大掃描半徑的探索效果，實驗中引入q統(tǒng)計量作為分層異質(zhì)性的探測方法[28]。通過對舊金山地區(qū)2018-2020年的火災事件進行預處理，獲取可進行時空掃描的數(shù)據(jù)。對研究的40個監(jiān)測點進行編號(表1)，并篩選統(tǒng)計出監(jiān)測點相應的火災事故數(shù)據(jù)(表2)。

表1 監(jiān)測點對應編號Table 1 Corresponding numbers of monitoring points

表2 監(jiān)測點火災事件統(tǒng)計Table 2 Statistics of fire incidents at monitoring points

2 實驗與結果分析

采用回顧性時空重排掃描分析方法對40個監(jiān)測點進行實驗，選取2019年4個不同時間段(3月1-14日、6月1-14日、9月1-14日、12月1-14日)的實驗數(shù)據(jù)，作為預警對比數(shù)據(jù)。為避免局部異常，基于2018年相同研究月份的數(shù)據(jù)，繪制相應監(jiān)測點的火災事件分布箱線圖(圖2)。

圖2 2018年3月、6月、9月、12月監(jiān)測點火災事件分布Fig.2 Fire incident distribution at monitoring points in March,June,September and December 2018

2018年3月、6月、9月、12月，分別在一天內(nèi)發(fā)生大于2起、1起、0起、1起火災事件(圖2中粗橫線對應數(shù)值)，即異常行為。本研究分析時間為一周(d=7)，利用式(8)計算出各監(jiān)測點的異常事件數(shù)閾值分別為14、7、0、7。根據(jù)閾值大小以及2019年3月8-14日、6月8-14日、9月8-14日、12月8-14日各監(jiān)測點實際發(fā)生火災的情況(圖3中虛線對應刻度值為觀測點研究時間段內(nèi)發(fā)生異常事件數(shù)的閾值)，篩選出研究時間段內(nèi)發(fā)生火災事件的實際異常監(jiān)測點(表3)。因2019年3月8-14日各監(jiān)測點均無火災事件發(fā)生，故實驗中省略3月的數(shù)據(jù)。針對其余3個時間段的數(shù)據(jù)，參照本文方法，以[0 km,10 km]為二分法的搜索起始區(qū)間，以0.1 km為搜索區(qū)間的最小間隔精度(即當搜索區(qū)間的間隔小于等于100 m時停止搜索)，通過比較試探選取最大掃描半徑，利用式(9)計算歷史預警準確率，連續(xù)二分縮小搜索區(qū)間，得到最佳的最大掃描半徑為1.25 km(表4)。

表3 2019年6月、9月、12月實際異常監(jiān)測點統(tǒng)計Table 3 Statistics of actual abnormal monitoring points in June,September and December 2019

圖3 2019年3月8-14日、6月8-14日、9月8-14日、12月8-14日實際異常監(jiān)測點Fig.3 Actual abnormal monitoring points from March 8th to 14th,June 8th to 14th,September 8th to 14th,and December 8th to 14th,2019

表4 2019年6月、9月、12月基于本文方法的掃描結果歷史預警準確率Table 4 Historical warning accuracy of the scanning results based on the proposed method in this paper in June,September and December 2019

地理現(xiàn)象普遍具有空間分異性。分異及因子探測器可以探測因變量(研究區(qū)域火災事件發(fā)生的預警準確率)的空間分異性，探測影響因子(最大掃描半徑)對因變量的空間分異性解釋程度(本文用q統(tǒng)計量表示[28])。對本文方法在不同最大掃描半徑下的空間聚集情況進行分析(圖4)，可以看出，第1組實驗數(shù)據(jù)(2019年6月)的q統(tǒng)計量在最大掃描半徑為1.25 km時最大，且接近1；第2組實驗數(shù)據(jù)(2019年9月)的q統(tǒng)計量隨最大掃描半徑的增大而減小，在最大掃描半徑為1.25 km時，q統(tǒng)計量不再發(fā)生變化；第3組實驗數(shù)據(jù)(2019年12月)的q統(tǒng)計量隨最大掃描半徑的增大而增加，在最大掃描半徑為1.25 km時，q統(tǒng)計量不再發(fā)生變化。綜上，在最大掃描半徑為1.25 km時，空間分異性顯著，影響因子對研究區(qū)域火災事件發(fā)生的預警準確率解釋度較好。

圖4 3組實驗數(shù)據(jù)在不同最大掃描半徑下的空間異質(zhì)性情況Fig.4 Spatial heterogeneity of three groups of experimental data under different maximum scanning radii

本文借助測試集(2020年1月8-14日數(shù)據(jù))驗證各方法的預測效果與泛化能力。如表5所示，本文方法與文獻[12]、文獻[18]方法的RMSE、MAE與MAPE相同，表明3種方法在泛化能力、預測效果方面一致，但本文方法選取的最大掃描半徑最小、耗時最短，優(yōu)于另外兩種方法。

表5 3種方法泛化能力、預測效果及耗時比較Table 5 Comparison of generalization ability,prediction effect and time consumption of the three methods

為減少實驗結果的偶然性和隨機性，選取2020年1月1-7日、2月1-7日、3月1-7日3組數(shù)據(jù)進行對比實驗，并將2020年1月8-14日、2月8-14日、3月8-14日的數(shù)據(jù)作為驗證數(shù)據(jù)。用3種最大掃描半徑進行時空重排掃描并統(tǒng)計掃描耗時(表6)，發(fā)現(xiàn)本文方法在保證預警準確率的前提下，選擇1.25 km作為最大掃描半徑，耗時最短，效率最高。

表6 2020年1月、2月、3月時空重排方法及其優(yōu)化方法時空重排掃描結果預警準確率及耗時對比Table 6 Comparison of early warning accuracy and time consumption of the scanning results of spatiotemporal rearrangement scan statistic method and its optimization method in January,February and March 2020

3 結論與展望

時空重排掃描統(tǒng)計方法是時空事件異常探測的常用方法，為擴展其普適性，本文提出一種基于歷史預警準確率的時空重排掃描最大掃描半徑優(yōu)化方法：考慮歷史數(shù)據(jù)集的同期平均值對實際異常點的影響，采用二分法對歷史時空數(shù)據(jù)集進行回顧性時空重排掃描統(tǒng)計分析；選擇預警準確率最高的最大掃描半徑作為前瞻性掃描統(tǒng)計分析的最大掃描半徑。該方法在保證預警準確率的前提下，縮短了尋找合適最大掃描半徑的計算時間，從而優(yōu)化了時空重排掃描統(tǒng)計方法的性能；在半徑選擇過程中，雖然會占用一些計算資源，但考慮到同一觀測點的長期觀測研究，可為未來預警節(jié)省更多計算資源。該方法仍存在不足之處：在篩選實際異常點時，使用歷史數(shù)據(jù)的同期均值作為判斷閾值，掃描預警可能會遺漏發(fā)生事件數(shù)較少的部分異常點；在選擇合適的最大掃描半徑時，需對研究區(qū)域歷史數(shù)據(jù)進行處理，對于歷史事件數(shù)量較少或缺失的監(jiān)測點，會降低其普適性。今后將繼續(xù)對上述問題進行改進。