龍 華,王 美,楊 威,邵玉斌
(昆明理工大學(xué) 信息工程與自動化學(xué)院,昆明 650500)
E-mail:longhua@kmust.edu.cn
空間和時空掃描模型在預(yù)測預(yù)警模型中一直備受關(guān)注,自kulldorff[1,2]提出相應(yīng)模型后,被廣泛用于流行病[3-5]、生態(tài)學(xué)[6],犯罪學(xué)[7,8]等諸多領(lǐng)域中.在時空掃描模型中,每個掃描窗口為一個集群對象,模型定義底面為掃描覆蓋區(qū)域,高度為時間閾值的圓柱體為掃描窗口,根據(jù)窗口內(nèi)外統(tǒng)計的發(fā)生與未發(fā)生事件得到似然函數(shù),引入蒙特卡洛方法基于似然函數(shù)檢測得到集群評估統(tǒng)計的顯著性[9],在最大掃描范圍和最大掃描時間閾值內(nèi),不斷迭代空間和時間得到各圓柱體窗口的顯著性,從而獲取到較高風(fēng)險聚集區(qū)域.對時空掃描模型課題的研究中,決定掃描窗口的掃描半徑、掃描時間閾值和掃描窗口形狀,決定掃描統(tǒng)計量的測試統(tǒng)計量假設(shè)分布以及次級掃描區(qū)域即除了最大似然比對應(yīng)區(qū)域以外的掃描區(qū)域處理等均是討論的熱點問題[10].
掃描窗口問題的討論中,Patil和Taillie[11]創(chuàng)新性的采用上層設(shè)置的概念減小窗口大小使得能檢測任意形狀的簇,提升集群檢測能力,Duczmal和Assun??o[12]則采用模擬退火方法,基于圖的算法連通子圖空間找到似然比局部最大的區(qū)域,Kulldorff[13]和Christiansen[14]討論了使用橢圓形窗口替代經(jīng)典模型中圓形掃描窗口,并證明在特定數(shù)據(jù)集中有一定的優(yōu)勢;測試統(tǒng)計量假設(shè)分布的討論中,較有代表性的是2014年Can?ado[15]提出了零膨脹掃描統(tǒng)計量(zero-inflated Poisson ZIP),在該模型中當數(shù)據(jù)集包含大量零時,為了提高掃描結(jié)果中的空間精確度和空間召回率,將零分為符合泊松分布的抽樣零和獨立隨機變量得到的結(jié)構(gòu)零,對該模型的有效性檢驗可查閱Allévius 和H?hle[16]于2019做的模擬研究實驗,此外還有測試統(tǒng)計量符合超幾何分布的情況等.在前述方法模型中均會涉及一個問題,即掃描過程會產(chǎn)生成千上萬相互重疊的重要集群,這類集群會模糊較風(fēng)險區(qū)域的判斷,傳統(tǒng)時空掃描以及前面提到的改進模型均按Satscan層次結(jié)構(gòu)將其直接刪除,只報告不相交不重復(fù)的集群,但很明顯該方法會損失很多對信息量更有貢獻的集群[17].
為了充分利用潛在集群的信息,Gangnon[18]使用加權(quán)平均似然比檢驗統(tǒng)計量,將重疊集群用似然比加權(quán)合并,該方法只適合集群相對較少的情況,集群較多時會合并成較大而無實際報告意義的集群,李小洲和王勁峰[19]用優(yōu)化選擇格網(wǎng)點間隔方法和多重排序算法,以達到減少候選聚集區(qū)域的遺漏,并在較短的時間內(nèi)刪除所有的重復(fù)候選聚集的目的,該方法刪除的是完全重復(fù)的集群,仍然留下大量交集較大的集群模糊較風(fēng)險區(qū)域的判斷,同樣不能達到很好處理次級區(qū)域的效果,Gangnon與李小洲各自所提方法總的來說都存在不能量化分析實際效果的問題,直到2016年Han和Junhee[17]提出使用基尼系數(shù)--一種量化分析次級區(qū)域的方法來確定最優(yōu)的集群,如當一個大的集群里面包含幾個小的集群時,比較大的集群根據(jù)觀察數(shù)據(jù)和預(yù)期數(shù)據(jù)得到的基尼系數(shù)與其包含的小集群同樣方法計算得到的基尼系數(shù),選擇報告基尼系數(shù)較大的區(qū)域,減小信息量的損失,從所舉例子可知其使用對象導(dǎo)致Han和Junhee所提方法有一定的局限性,Han和Junhee所提方法是在保持傳統(tǒng)方法使用層次結(jié)構(gòu)刪除掃描重復(fù)窗口的情況下,再次使用基尼系數(shù)做判斷,所以分析的過程中會存在很多問題,如文獻[19]說的他們只考慮了緊湊的簇所以顯示結(jié)果良好,但實際使用中并不理想,如文獻[20]提出的基于基尼系數(shù)分析時在即使單個較大的群集確實有一定意義的情況下,往往會導(dǎo)致報告的是多個較小的群集.
為了充分利用潛在集群的信息,同時不出現(xiàn)Han和Junhee誤報告的問題,我們提出一種基于信息量的時空深度掃描模型IN-scan model,引入信息量I(p,G),p為統(tǒng)計值,G為基尼系數(shù),從信息量的角度出發(fā)量化每一個掃描集群,在不做刪除重復(fù)掃描區(qū)域,充分利用有效信息的情況下,報告出較有意義的集群.
時空掃描模型,主要是找出聚集性較高,即置信度(1-p)較大的集群,模型中似然函數(shù)值反應(yīng)一個窗口為聚集域的可能性,即最可能的聚集域是最大化對數(shù)似然比對應(yīng)的掃描窗口,此時測試統(tǒng)計量T=maxZLog(LR(z)).假設(shè)我們研究區(qū)域是Z,z表示某個掃描窗口,其似然比表示如下:
(1)
其中,cz和μz分別表示掃描窗口z內(nèi)的觀察案例和預(yù)期案例,C=∑z∈Zcz和N=∑z∈Zμz分別表示掃描區(qū)域Z內(nèi)的總觀察案例和總預(yù)期案例,在模型中如果只對較高可能為聚集域的掃描分析感興趣則使用如公式1中的指示函數(shù)I(cz>μz),若考慮的是較不可能為聚集域的使用I(cz<μz),若二者均考慮則刪除指示函數(shù),在本文的研究中僅考慮第一種情況.
得到各掃描域的似然比值,還需要對其進一步分析屬于非隨機的置信度(1-p),求解p值目前有兩個方法一個是測試統(tǒng)計量T看做近似服從于極值分布[21]做計算,一個是采用蒙特卡羅(Monte Carlo)假設(shè)檢驗方法求解,因為T真實服從的分布還有待研究,故文中選后者方法.使用Monte Carlo計算得到p=rank(LLR)/(M+1),LLR為真實數(shù)據(jù)集掃描域計算得到的似然比,M個隨機數(shù)據(jù)集是根據(jù)真實數(shù)據(jù)集采用重排算法得到的,M個隨機數(shù)據(jù)集對應(yīng)計算得到M個LLR′,將真實的LLR與M個LLR′放在一起由大到小排序返回排序值rank(LLR),由此計算得到p值.p值越小對應(yīng)聚集域?qū)儆诜请S機的置信度越大.后面提到的pi表示第i個掃描窗口對應(yīng)的統(tǒng)計值.
LR-scan模型是使用了對數(shù)似然比為測試統(tǒng)計量并引入了Junhee提出的基尼系數(shù)判別方法的最新時空掃描模型,該模型首先是根據(jù)置信度得到最大的聚集域?qū)χ貜?fù)域進行一次刪除,然后根據(jù)基尼系數(shù)對重復(fù)域進行再一次篩選.該方法依然存在信息浪費與報導(dǎo)有偏差問題,故本文提出IN-scan模型.
LR-scan模型已于SaTScan中更新,SaTScan是用于時空掃描統(tǒng)計分析的一個開源軟件,其集成了很多時空掃描的方法.當前,LR-scan 模型是最新添加的方法,因此文中也將LR-scan方法稱為最新SaTScan方法.
時空掃描方法屬于預(yù)警模型,集群置信度越大越容易報警,即其區(qū)域內(nèi)發(fā)生事件的概率越大.在IN-scan 模型中,我們將每個集群的置信度視為集群內(nèi)事件發(fā)生的概率,進而計算事件均值和基尼系數(shù),然后根據(jù)求得的p和G計算信息量.圖1是我們提出模型的結(jié)構(gòu),首先與LR-scan 模型一致以對數(shù)似然比為測試統(tǒng)計量進行掃描計算,但IN-scan不對掃描結(jié)果進行刪除,而是對其求均值,以進行二次以信息量為測試統(tǒng)計量的掃描計算,然后依次報導(dǎo)較風(fēng)險的聚集域.本節(jié)將給出各個參數(shù)的計算方法.
圖1 基于信息量的時空深度掃描模型結(jié)構(gòu)Fig.1 IN-scan model structure
模型特點為:
1)只考慮較高可能發(fā)生性事件,故始終有地點i觀察案件數(shù)oi大于等于預(yù)測案件數(shù)μi;
2)使用均值數(shù)據(jù),n′i和μ′i;
3)于L中二次掃描,測試統(tǒng)計量使用信息量I(p,G).
假設(shè)圖2中3個集群的掃描條件均一樣,且掃描時間也一致,此時結(jié)合表1記錄以地點A為例有:nA1=nA2=nA3=nA,μA1=μA2=μA3=μA,加入事件發(fā)生的概率,計算地點A數(shù)據(jù)的均值:
(2)
地點B和C與A分析一致.這樣計算完成后相當于把cluster1、cluster2、cluster3的信息映射到掃描區(qū)域L中,二次掃描時只用分析L中的數(shù)據(jù)即可,見方法介紹模塊.
圖2 均值分析-聚集域分布圖Fig.2 Mean analysis-aggregation area distribution map
表1 均值分析-聚集域記錄Table 1 Mean analysis-aggregation area recording
基尼系數(shù)[22]是洛倫茲曲線的一個度量值,洛倫茲曲線[23]主要用于經(jīng)濟領(lǐng)域,如圖3所示基尼系數(shù)G=S1/(S1+S2),G越大說明收入分配越不平等.Han和Junhee首次將G用于時空掃描模型中,其橫縱軸分別表示觀察到案例的累加
圖3 聚集域基尼系數(shù)顯示圖Fig.3 Aggregation area of Gini coefficient
百分比和預(yù)期案例的累積百分比,G=0時,即說明該集群無顯著性,而G越大,對應(yīng)掃描域的置信度(1-p)越大,Han使用G是為了報告出更有意義的集群,在我們的研究模型中,研究G亦是為了報告出能提供更多信息量I的集群.
在對IN-scan model模型的研分析中,由于我們使用的是均值數(shù)據(jù),地點i觀察案件數(shù)n′i與預(yù)測案件數(shù)μ′i能綜合反映i點的風(fēng)險度,此時地點i與地點j信息相對獨立,所以在掃描區(qū)域L中進行二次掃描時,測試統(tǒng)計值信息量I(p,G)中的G我們使用掃描域內(nèi)的n′從小到大排序后計算得到,見圖3,假設(shè)掃描域為table1中的cluster3,掃描時間為1天,此時若n′c>n′B>n′A,見圖4中t=1的G顯示圖,圖中陰影部分占下三角的面積百分比即為G值,t=2和t=3類似分析,從圖3中可看出就基尼系數(shù)來說,cluster3掃描時間為1天時G較大.在IN-scan model中,結(jié)合p與G值計算各條件下的集群的信息量,報告I最大的集群.
圖4 t=1的G顯示圖Fig.4 G display of t=1
在考慮全部潛在集群的分析中,計算基尼系數(shù)時,使用事件均值計算可減小誤差,如圖5是基于Satscan樣本數(shù)據(jù)集紐約市醫(yī)院記錄的發(fā)燒案例,在Satscan和我們最終模擬實驗中知道圖中兩個掃描區(qū)域是以171為中心的窗口較異常優(yōu)先報告,當直接使用觀測值計算基尼系數(shù)時,見表2,169對應(yīng)的基尼系數(shù)明顯高于171對應(yīng)的值,在傳統(tǒng)模型中由于使用層次結(jié)構(gòu)的方法提前刪除了169對應(yīng)窗口,所以不會報告,說明在不做刪除重復(fù)掃描區(qū)域的情況下,直接使用觀測結(jié)果和基尼系數(shù)做判斷會產(chǎn)生誤差,基于此問題,我們充分利用有效信息,對各掃描點做均值后再次掃描計算基尼系數(shù),比較出基尼系數(shù)較大的窗口,表3是使用提出方法實驗的結(jié)果,可看出能正確比較出171對應(yīng)窗口.考慮到基尼系數(shù)和顯著性統(tǒng)計值,均能反應(yīng)信息量的大小且相互影響,所以我們將引入信息量I(p,G)分析最終結(jié)果.
圖5 相交窗口分析Fig.5 Intersection of windows with center 169 and 171
在很多領(lǐng)域信息量都是作為選優(yōu)條件的評估指標,如地質(zhì)學(xué)[24]和工程學(xué)[25],尤其在預(yù)警領(lǐng)域,很多時候會由于信息量較少導(dǎo)致大量的假預(yù)警報告,所以更應(yīng)該考慮信息量因素,信息量表達式為I(x)=-log2p(x),其中p(x)表示事件x發(fā)生的概率,p(x)越小信息量I越大.對于我們模型來說I∝(1/p),I∝G,所以有:
表2 Han和Junhee方法重疊圓分析結(jié)果Table 2 Analysis result of the intersecting windows by Han and Junhee
表3 使用均值方法重疊圓分析結(jié)果Table 3 Analysis result of the intersecting windows with mean value
備注:表中center:掃描中心點,r (Km):掃描半徑(單位是千米),ids:掃描窗口包含的地點集合,observed:掃描窗口內(nèi)總的觀察案例數(shù),expected:掃描窗口內(nèi)總的期望案例數(shù),p:顯著性統(tǒng)計值,G:基尼系數(shù).
I(p,G)=-logp-log(1-G)p∈(0,1],G∈[0,1)
(3)
顯著性統(tǒng)計值p越小,基尼系數(shù)G越大,信息量I越大.集群I越大,越優(yōu)先報警.
為了評估基于信息量的時空深度掃描模型的空間精度,結(jié)合機器學(xué)習(xí)中的評估方法[26]以及2012年Neil[27]提出的空間精確度SP和空間召回率SR,文中使用F-Score作為模型評估指標.此處若Z*表示檢測到的聚集空間區(qū)域,ZT表示真實爆發(fā)區(qū)域,則有:
(4)
|Z*|和|ZT|分別是Z*和ZT包含的地點數(shù),若向集群中多添加地點將有利于召回率同時精確率會有所下降,若從檢測的區(qū)域移除地點會造成以召回率為代價提高精確率結(jié)果,故需要同時權(quán)衡SP和SR.在模型中我們不希望|Z*|很大,導(dǎo)致SP很小,但SR很大的情況,SP和SR同樣重要,故評估參數(shù)F-Score表達式如下:
(5)
時空掃描統(tǒng)計方法應(yīng)用于多領(lǐng)域,2008年Vadrevu[28]和Tuia應(yīng)用時空掃描方法監(jiān)測火災(zāi)風(fēng)險區(qū)域取得一定成效,本次實驗中我們將采用舊金山的火災(zāi)數(shù)據(jù)展開討論.實驗數(shù)據(jù)集來自舊金山地區(qū)數(shù)據(jù)協(xié)調(diào)網(wǎng)站DataSF1提供的“Fire Department Service”火災(zāi)記錄公開數(shù)據(jù)集,數(shù)據(jù)集中包含舊金山18年和19年上半年共四萬余條火災(zāi)記錄數(shù)據(jù),且網(wǎng)站在實時更新中,考慮實驗的回顧性驗證所以實驗數(shù)據(jù)集隨機提取2018年1月和3月的火災(zāi)事件進行研究.圖6是數(shù)據(jù)集中舊金山火災(zāi)發(fā)生的地點標記圖.
圖6 舊金山火災(zāi)觀測地點圖Fig.6 San Francisco fire observed sites
本節(jié)中用IN-scan表示我們提出的基于信息量的時空深度掃描模型,LR-scan表示現(xiàn)有最新時空掃描模型,實驗中使用IN-scan與LR-scan的掃描結(jié)果比較,說明新模型的有效性.首先選取時間2018年1月8號-14號的數(shù)據(jù)進行分析,掃描半徑設(shè)置為r=5Km,最大掃描時間t=3day,首先基于似然比測試統(tǒng)計量不刪除重復(fù)窗口進行掃描,并使用蒙特卡羅計算掃描得到的每個掃描窗口的置信度,表4中輸出前5個掃描結(jié)果.最大掃描時間為3天,故實驗結(jié)果中包括掃描時間分別為1天、2天和3天的情況,實驗中會分開計算,這里以分析t=1為例,如表4中第一個模塊中最大似然比對應(yīng)包含掃描區(qū)域為12,獲取以12為圓心的所有同心圓包含的地點集合D12,之后獲取以D12包含的地點為掃描中心的所有窗口,然后按照事件均值求解方法算出各地點的均值,計算結(jié)果見表5,表中obaerves和expected是地點的觀測值和預(yù)期值,obaerves-mean和expected-mean是地點的觀測值均值和預(yù)期值均值,表中有觀測值不一定有觀測期望值,如locations 2 觀測值為2,但觀測期望值為0,是因為包含有地點2的全部窗口置信度都為0導(dǎo)致,也有觀測值很小,但觀測期望值很大,如locations 12,因為地點12自身的觀測值很大,包含有地點12的窗口置信度多數(shù)偏高,故觀測期望值很大.基于期望值求取各窗口的基尼系數(shù),使用基尼系數(shù)和置信度隨之求取各窗口的信息量,表6即為表4中第一個掃描窗口計算得到的信息量I,在計算結(jié)果中留下最大信息量對應(yīng)窗口,同時取剩下的第一個集群重復(fù)分析取出信息量最大的窗口,如此循環(huán)直到I=0停止.
表4 1.8-1.14,r=5Km前五個掃描結(jié)果Table 4 Top five scan results from 1.8-1.14,r=5K
表5 1.8-1.14 t=1 均值表Table 5 Mean values from 1.8-1.14,t=1
表6 地點12,t=1信息量計算Table 6 Information calculation of location 12,t=1
實驗選取2018年1月8號-14號與2018年3月8號-14號做掃描對比實驗,首先確定真實爆發(fā)區(qū)域ZT.實驗最終比較的是模型對異常地點檢測的準確性和回歸率,因為我們考慮的是短期爆發(fā)事件,故文中爆發(fā)定義為掃描時間后對應(yīng)的第一周內(nèi)發(fā)生事件,見圖7可得1.15-1.21爆發(fā)地點集[1,2,3,7,8,10,11,12,13,14,16,20,23,24,25,26,27,29,32,34],3.15-3.21爆發(fā)地點集為[1,2,3,5,6,7,10,12,13,15,17,20,21,25,33,39].對比結(jié)果見表7-表9.
表7 1.8-1.14,r=5 km掃描結(jié)果比較Table 7 F-Score comparison from 1.8-1.14,r=5 km
表8 3.8-3.14,r=5 km掃描結(jié)果比較Table 8 F-Score comparison from 3.8-3.14,r=5 km
表9 3.8-3.14,r=3 km掃描結(jié)果比較Table 9 F-Score comparison from 3.8-3.14,r=3 km
表7掃描時間1月8號-14號,掃描半徑為5km,運用評估模塊的計算方法,得到相比LR-scan,IN-scan的F-Score提高14%;為了驗證方法的有效性使用時間段3月8號-14號,掃描半徑為5km進行掃描,見表8相比LR-scan,IN-scan的F-Score提高18%;實驗3選取同一時間段3月8號-14號,但不同掃描半徑(3km)進行掃描,見表9相比LR-scan,IN-scan的F-Score提高12.6%.總體比較,IN-scan模型較現(xiàn)時空掃描模型F-Score性能評估值提升10%以上.
圖7 1.15-21與3.15-21火災(zāi)爆發(fā)統(tǒng)計圖Fig.7 Fire outbreak statistics comparisonbetween 1.15-21 and 3.15-21
從實驗結(jié)果來看,現(xiàn)存在的時空掃描方法中對掃描過程產(chǎn)生的大量重復(fù)域直接進行刪除確實存在一定的缺陷,或是留下了置信度大但信息量很小的聚集域,或是剔除了置信度小但信息量大的聚集域,即使向LR-scan模型加入了基尼系數(shù)進行第二次判斷也不能避免此類問題,一定程度上會加重該問題,因為LR-scan模型中進行基尼系數(shù)計算時,往往會報告多個較小的群集,而該類集群往往是置信度大但信息量較小.可見我們提出的IN-scan模型在綜合處理權(quán)衡聚集域的置信度與信息量上有一定優(yōu)勢.
在時空掃描方法中,每次掃描都會產(chǎn)生成千上萬大量有重復(fù)的集群,現(xiàn)有研究方法會按Satscan層次結(jié)構(gòu)刪除此類集群,或者用似然比加權(quán)集群分析.每個集群都會對應(yīng)一個置信度,都能得到直觀反應(yīng)他們異常度的基尼系數(shù)值,故不應(yīng)該直接刪除或合并,所以為了充分利用掃描域提供的信息,量化分析選取報警區(qū)域,我們提出了一種基于信息量的時空深度掃描模型,深度是因為進行了二次掃描,基于信息量是因為我們利用每個聚集域的觀測值均值和預(yù)期值均值計算得到基尼系數(shù)G和基于似然比掃描的得到的統(tǒng)計值p得到信息量IN展開比較報告較有意義的集群,經(jīng)過實驗對比,確實有一定有效性.
另外本文介紹的IN-scan模型具有很強的穩(wěn)定性.現(xiàn)有模型中當存在多個集群時需要使用適用于多個集群的估計方法進行進一步分析,而IN-scan模型可以對重復(fù)的多個群集進行分析,直接計算相交群集中的事件所包含案例.就目前的研究來說,我們僅考慮單個事件類型,而忽略了事件間的內(nèi)部結(jié)構(gòu)差異造成的影響等,即考慮多元變量問題,這是我們下一個主要研究的方向.