楊秋菊，任杰，向晗

陜西師范大學物理與信息技術學院,西安 710061

0 引言

亞暴是一種發(fā)生在地球夜半球太空區(qū)域的巨大能量瞬間釋放現(xiàn)象(Akasofu,1964),亞暴發(fā)生過程伴隨著極光形態(tài)和亮度的劇烈變化,研究人員將這種極光活動稱為極光亞暴.極光亞暴的發(fā)生伴隨著由電離層電流引起的強烈磁干擾(Akasofu,2017),會對人類活動造成很大的影響,例如引起地球同步軌道上的衛(wèi)星充電(李柳元等,2006)、高緯度地區(qū)無線電通訊中斷(Wen et al.,2012;Chernyshov et al.,2020)和供電系統(tǒng)故障(Boteler,2019)等.目前,極光亞暴的產生機制仍存在爭議(Schindler,1974;Angelopoulos et al.,2008;Maimaiti et al.,2019),行星際磁場由北轉向南是否是亞暴發(fā)生的原因仍有待確認.從海量的觀測數(shù)據(jù)中自動檢測亞暴事件并確定其發(fā)生位置有助于統(tǒng)計分析亞暴的產生機制(Ieda et al.,2019).

目前,檢測極光亞暴事件的途徑可以分為三種.第一種途徑是從海量極光圖像中人工挑選亞暴事件,如Liou(2010)基于Polar衛(wèi)星的全域極光圖像手動標記極光亞暴起始時刻及其發(fā)生位置(磁地方時-地磁緯度),Ieda等(2019)結合全域極光圖像和全天空極光圖像確定亞暴起始時刻,但人工標注費時費力且在圖像數(shù)據(jù)量非常大時容易出現(xiàn)漏標、誤標等問題(He et al.,2014).第二種途徑是通過各類物理參數(shù)自動檢測亞暴事件,如Sutcliffe(1997)基于Pi2脈沖利用神經網絡的方法自動檢測亞暴事件,Maimaiti等(2019)結合太陽風速度、質子數(shù)密度和行星際磁場基于深度學習的方法自動檢測亞暴.雖然該途徑克服了人工挑選的局限性,但亞暴事件真實發(fā)生時間和所測量的物理指標之間存在時延(Meng and Liou,2004;Ieda et al.,2019),導致檢測結果存在一定的誤差.為了解決以上兩種問題,人們提出了第三種途徑,基于衛(wèi)星紫外成像儀(Ultraviolet Imager,UVI)圖像利用機器學習方法自動檢測亞暴事件.如楊秋菊等(2013)基于模糊C均值聚類通過分析亮斑的變化情況確定亞暴的起始時刻,該方法雖然實現(xiàn)簡單并取得了較高的召回率,但準確率有待提升且需要人工設置多個閾值得到檢測結果;結合亞暴事件發(fā)生區(qū)域的地磁先驗信息,利用SCSLD(Shape-Constrained Sparse and Low-Rank Decomposition)算法通過粗檢測、序列運動分析、細檢測實現(xiàn)亞暴事件檢測,該方法在多年的極光觀測上獲得了很好的檢測效果(Yang et al.,2016),但實現(xiàn)過程中需要人工多次參與分析算法結果,操作復雜而且運行效率較低,實際應用難度大.

基于深度學習的視頻時序行為檢測研究(Zhao et al.,2017;Lin et al.,2018)為亞暴自動檢測提供了新的選擇,但與現(xiàn)有時序行為檢測不同的是亞暴事件檢測重點關注極光點亮的時刻(onset),準確檢測亞暴事件需要構建密集的特征序列;此外,衛(wèi)星觀測數(shù)據(jù)時間分辨率較低(0.0056幀/s VS常規(guī)視頻24幀/s),極光亞暴持續(xù)的視頻幀長度較短,模型需要包含輸入圖像序列的全部時序信息.基于此,本文利用深度學習技術提出了一種端到端的亞暴檢測網絡(Substorm Onset Detection Network,SODN),該方法在取得亞暴檢測高準確率的同時,極大提高了亞暴檢測速度.SODN由兩部分組成:(1)亞暴特征提取模塊;(2)亞暴起始檢測模塊.亞暴特征提取模塊基于雙流網絡(Simonyan and Zisserman,2014)構建,包含兩個并行的卷積神經網絡:空間網絡和時序網絡.其中,空間網絡用于提取亞暴空間形態(tài)特征,時序網絡用于提取亞暴時序運動特征,本文將上述兩種特征級聯(lián)獲得亞暴特征表示.亞暴起始檢測模塊利用三個一維時序卷積(Lea et al.,2016)融合上述兩種特征,輸出亞暴起始的概率序列.與第三種途徑不同的是本文方法利用雙流網絡自動提取亞暴時-空特征,不需要手動設計亞暴特征表示.此外,利用三個一維時序卷積層構建亞暴起始檢測模塊獲取亞暴起始的邊界信息,無需繁瑣的篩選分析過程直接得出檢測結果,極大提高了亞暴檢測效率.

本文的結構安排如下:第1節(jié)簡要介紹極光亞暴事件特點和實驗數(shù)據(jù)來源,第2節(jié)詳細介紹SODN框架,實驗結果和分析在第3節(jié)給出,第4節(jié)為結論和討論.

1 數(shù)據(jù)集與數(shù)據(jù)預處理

1.1 UVI圖像

本文實驗數(shù)據(jù)由Polar衛(wèi)星攜帶的UVI拍攝(Torr et al.,1995).該衛(wèi)星于1996年2月24日進入高橢圓形軌道,遠地點高度為9個地球半徑,近地點高度1.5個地球半徑,軌道傾角86°,周期約17.5 h.Polar衛(wèi)星于2008年4月完成相關計劃,在此期間收集了大量的UVI極光圖像,首次詳細觀察了地磁亞暴發(fā)生期間極光的演變過程.Polar衛(wèi)星的軌道平面在觀測期間緩慢向南移動,本文選用Polar衛(wèi)星位于北半球期間1996年12月,1997年1月、2月、12月和1998年12月采集的UVI圖像作為數(shù)據(jù)集,其中共包含371個亞暴事件;從每月中隨機選取5天作為測試數(shù)據(jù)集,共包含亞暴事件66個;其余為訓練數(shù)據(jù)集,共包含亞暴事件305個.

UVI采集了多個波段的極光圖像(Torr et al.,1995),借鑒文獻(Kullen and Karlsson,2004)的處理方式本文只使用LBHL(Lyman-Birge-Hopfield long)濾波通道下的UVI圖像.該波段下UVI圖像的空間分辨率為30 km,圖像覆蓋范圍約4104000 km2,時間分辨率為3 min,大小為228×200;UVI圖像信息除了像素亮度外,還包含了每個像素對應的MLT、MLAT等信息(楊秋菊等,2013).本文采用Liou(2010)提供的人工標注亞暴事件列表,將該標注作為挑選訓練數(shù)據(jù)集和衡量亞暴檢測結果的基準.

1.2 數(shù)據(jù)預處理

Polar衛(wèi)星采集UVI圖像時位置在不斷變化,UVI圖像的同一位置在不同圖像中對應著不同的磁地方時和地磁緯度(楊秋菊等,2013),不利于從圖像序列中統(tǒng)計亞暴發(fā)生時亮斑的位置變化.本文使用和我們前期工作(楊秋菊等,2013)相同的方法將UVI圖像轉換至磁地方時-地磁緯度(MLT-MLAT)坐標下,轉換后圖像中的相同像素位置對應相同的磁地方時和地磁緯度,具體步驟為:

(1)根據(jù)亞暴發(fā)生的時間和位置,重點關注夜側和高緯區(qū)域,MLT-MLAT網格圖像限定為磁地方時18—24 MLT及00—06 MLT,地磁緯度50°—90°MLAT;

(2)設定網格大小,轉換后的網格圖像大小為200×120((90°—50°MLAT)/0.2°=200行,12 h/0.1 h=120列);

(3)對網格中每個坐標點,找到其在UVI圖像中的位置,計算UVI圖像中落入每個網格的像素的均值,該均值作為每個網格點的像素強度.

由MLT-MLAT網格圖像的定義知,各網格圖像相同位置對應相同的MLT和MLAT值,極光亞暴過程中亮斑區(qū)域的形態(tài)變化在MLT-MLAT網格圖像中主要表現(xiàn)為水平方向和豎直方向的形態(tài)變化(如圖1所示),極光形態(tài)隨時間的演變規(guī)律可以直接從MLT-MLAT圖像序列中計算.

圖1為一組轉換后的MLT-MLAT網格圖像示例，該圖像展示了一個典型的極光亞暴事件發(fā)生過程,根據(jù)該過程中極光的形態(tài)特征可以將亞暴事件分為三個階段(Rostoker et al.,1980;Liang et al.,2018):

圖1 1996年12月04日21∶51∶41—22∶28∶10UT采集的UVI圖像轉為MLT-MLAT網格圖像的結果，MLT范圍為18—06，MLAT范圍為50°—90°.圖中粗體標注的時間為亞暴膨脹階段，其中22∶00∶34UT為亞暴事件的起始時刻，斜體標注的時間為亞暴消退階段Fig.1 The MLT-MLAT grid images of the UVI images captured on December 04,1996 between 21∶51∶41—22∶28∶10UT.The range of MLT is 18—06,and the range of MLAT is 50°—90°.The time marked in bold is the substorm expansion phase,where 22∶00∶34UT is the onset of the substorm,and the time marked in italics is the substorm recovery phase

(1)增長階段:地面成像儀可以觀測到極光弧向赤道向移動,衛(wèi)星成像儀觀測到極光卵區(qū)域仍保持靜止.

(2)膨脹階段:以極光卵局部區(qū)域突然點亮為開始標志,點亮區(qū)域迅速向極向和東西向膨脹,同時伴隨著亮度的增強,最后到達一個極向最大邊界.

(3)消退階段:該階段發(fā)生在點亮區(qū)域膨脹至極向最大位置后,表現(xiàn)為點亮區(qū)域亮度逐漸減弱、膨脹區(qū)域消退,最后恢復到亞暴之前的水平.

其中,增長階段的持續(xù)時間通常為30～60 min,膨脹階段和消退階段共持續(xù)1～2 h(Ebihara,2019).膨脹期內,從UVI圖像上觀測亞暴的兩個關鍵特征(Kepko and McPherron,2001;Frey,2004)是:(1)極光卵部分區(qū)域極光亮度迅速增強;(2)亮度增強區(qū)域向極向和赤道方向膨脹,尤其在磁地方時(MLT)21—03時和地磁緯度(MLAT)65°—75°之間.只有空間形態(tài)變化和時序變化同時滿足上述兩個特征才被稱為亞暴,只滿足第一個特征可能是偽暴(Nakamura et al.,1994;Aikio et al.,1999;Kullen and Karlsson,2004)或圖像中的噪聲點.基于此,SODN使用空間網絡和時序網絡分別關注亞暴形態(tài)信息和時序運動信息,最后通過時序卷積融合兩種信息獲得檢測結果.

2 SODN模型

本文提出的SODN模型如圖2所示.首先將轉換后的MLT-MLAT圖像序列分為若干個連續(xù)等長的片段,每個片段生成一個單元,一個單元包含片段內的單幀圖像和片段內的堆疊光流圖(x方向和y方向).然后利用特征提取模塊分別對每個單元提取亞暴空間特征和時序特征,并級聯(lián)兩個特征作為單元級特征向量.最后按時間順序拼接每個單元的特征向量作為亞暴起始檢測模塊的輸入,輸出為整個UVI圖像序列的亞暴起始概率向量,向量中的每個值為其對應單元作為亞暴起始的置信度分數(shù).

2.1 亞暴特征提取模塊

為了提取亞暴時-空特征,本文采用雙流網絡構建亞暴特征提取模塊,該網絡由空間網絡和時序網絡組成.如圖2所示,本文將單幀MLT-MLAT圖像輸入空間網絡提取亞暴空間形態(tài)特征,將堆疊的光流圖輸入時序網絡提取亞暴時序運動特征.光流計算采用TV-L1算法(Zach et al.,2007),空間網絡和時序網絡都采用具有多個歸一化層的BN-Inception(Ioffe and Szegedy,2015)作為其骨干網絡.該網絡由多個Inception模塊組成,在Inception網絡的基礎上使用多個小尺度卷積代替尺度較大的卷積降低模型參數(shù)量,并通過增加BN層使得模型更快收斂,使該模型在準確性和效率之間實現(xiàn)了較好的平衡(Wang et al.,2016).本文在BN-Inception的最后一個全連接層之前加入一個維度為50的全連接層,用于輸出空間特征向量和時序特征向量.

圖2 SODN模型示意圖Fig.2 The diagram of SODN

2.2 亞暴起始檢測模塊

圖3 亞暴起始時刻檢測模塊輸出的概率序列示意圖.選擇輸出概率值大于0.5的時刻作為候選亞暴起始點,在候選起始點中選擇概率峰值點作為檢測結果Fig.3 The diagram of the probability sequence outputted by the substorm onset detection module.The moments when the output probability is greater than 0.5 are selected as the candidate substorm onset points,and we choose the time corresponding to the maximum probability as the final detection result

一個時序卷積層可以表示為：Conv1d(cf,ck,Act)，其中Conv1d表示一維卷積操作，cf,ck,Act分別表示一維卷積的卷積核數(shù)量、卷積核尺寸和卷積層的激活函數(shù).亞暴起始檢測器由三個連續(xù)的時序卷積層組成，定義為:

Conv1d(256,3,Relu)→Conv1d(256,3,Relu)

→Conv1d(1,3,Sigmoid),(1)

3 實驗結果及分析

3.1 SODN和SCSLD方法的結果對比

為了充分評估SODN與SCSLD的檢測性能,本文選用準確率(precision)、召回率(recall)、F1值和每秒幀數(shù)(Frame Per Second,FPS)來衡量兩者的檢測精度和檢測效率,

(2)

(3)

(4)

其中，TP(true positive)表示檢測的亞暴事件在人工標記的亞暴事件列表中;FP(false positive)表示檢測的亞暴事件未出現(xiàn)在人工標記的亞暴事件列表中;FN(false negative)表示出現(xiàn)在人工標記的亞暴事件列表中的事件未被檢測到.

SODN與SCSLD的對比結果如表1所示.SCSLD方法在檢測亞暴過程中,人工多次參與了結果分析和篩選,所以準確率比SODN高1.83%.但從檢測效率方面來看，SODN處理一幀圖像僅需0.003 s，檢測速度高達393幀/s；更重要的是，將UVI圖像輸入SODN后能直接輸出檢測結果，即SODN實現(xiàn)的是端到端的亞暴自動檢測，在實際應用中非常簡便.而相比之下，SCSLD方法整個檢測流程復雜，需要經過粗檢測、時序運動分析、細檢測三步完成檢測，且在時序運動分析中需要對圖像進行兩次運動特征提取，僅其中一次運動特征提取就需要0.15 s(6.93幀/s)，完成檢測一幀的時間很長，從而使其檢測速度很低.

表1 SODN與SCSLD檢測結果Table 1 The detection results of SODN and SCSLD

3.2 SODN與人工標記的結果對比

SODN與人工標記的結果對比如表2所示，SODN從測試數(shù)據(jù)集共檢測出64個亞暴事件，其中與人工標記(Liou,2010)一致的56個，多檢亞暴事件8個,漏檢亞暴事件10個.分析SODN檢測結果中與人工標記不一致的案例(多檢和漏檢),我們根據(jù)其原因不同將這些案例大致歸為以下三類.

表2 SODN與人工標記檢測結果Table 2 The detection results of SODN and manual approach

(1)本文僅考慮了LBHL波段的觀測數(shù)據(jù)，丟失了亞暴起始關鍵信息

如圖4所示,人工標注列表中1996年12月11日09∶53∶53—09∶54∶36UT為亞暴事件起始時段,但SODN漏檢了該事件.這是因為人工標注亞暴事件時參考了多個波段的UVI極光圖像,而本文僅考慮了LBHL波段,對于圖4所示的亞暴事件來說,其演化速度較快,LBHL波段未采集上亞暴膨脹期極光亮斑從點亮到膨脹的關鍵信息,導致該事件從我們的數(shù)據(jù)集上看更像一個偽暴事件.類似的情況在1997年12月03日16∶16∶08—16∶16∶55UT也有出現(xiàn).

圖4 1996年12月11日09∶53∶36—10∶05∶52UT采集的UVI圖像，MLT范圍為18—06，MLAT范圍為50°—90°Fig.4 The UVI images captured on December 11,1996 between 09∶53∶36—10∶05∶52UT.The range of MLT is 18—06,and the range of MLAT is 50°—90°

(2)人工標注中漏標了部分真實存在的亞暴事件

圖5所示為1996年12月14日17∶46∶48—19∶32∶17UT采集的UVI圖像序列.分析該時段UVI圖像序列,用粗體標注的多幀圖像符合亞暴點亮和膨脹相的特征,是一個典型的亞暴事件,但在人工標注列表中該事件未被標注.該天13∶26∶44UT也存在類似的情況,而這些事件均被SODN檢測為亞暴事件.

圖5 1996年12月14日17∶46∶48—18∶14∶24UT采集的UVI圖像序列，MLT范圍為18—06，MLAT范圍為50°—90°Fig.5 The UVI images captured on December 14,1996 between 17∶46∶48—18∶14∶24UT.The range of MLT is 18—06,and the range of MLAT is 50°—90°

(3)較為復雜的亞暴事件

圖6所示為1998年12月14日14∶28∶39—14∶29∶16UT采集的UVI圖像序列,人工標注該時段為亞暴起始時段.分析該時段UVI圖像序列,14∶29∶36UT對應UVI圖像中白色圓圈標記一塊亮斑,其上方存在兩塊細長條亮斑且亮度相近,隨后該亮斑出現(xiàn)微弱的亮度增強和膨脹現(xiàn)象.上述情況增加了該亞暴事件的檢測難度且訓練數(shù)據(jù)集中類似的亞暴事件很少,導致了SODN漏檢了該事件.

圖6 1998年12月14日14∶29∶36—14∶31∶26UT采集圖像序列，MLT范圍為18—06，MLAT范圍為50°—90°Fig.6 The UVI images captured on December 14,1998 between 14∶29∶36—14∶31∶26UT.The range of MLT is 18—06,and the range of MLAT is 50°—90°

3.3 空間定位統(tǒng)計

亞暴發(fā)生的位置信息對亞暴研究同樣具有重要意義.Liou(2010)利用區(qū)域生長算法確定亞暴起始時刻圖像中的亮斑區(qū)域,再計算亮斑的質心位置作為亞暴發(fā)生的位置(MLT-MLAT坐標)進行統(tǒng)計分析.該方法的缺點是每處理一張圖像都需要采用人機交互的方式確認區(qū)域生長的種子點.為了解決這一問題,本文利用Grad-CAM(Selvaraju et al.,2017)方法定位亞暴發(fā)生位置.該方法利用網絡輸出對指定的卷積層輸出計算梯度得到一個二維矩陣,矩陣中的每個數(shù)值表示該區(qū)域對網絡決策分類的重要程度,將該二維矩陣渲染為偽彩圖與輸入圖像疊加生成類激活圖(Class Activation Map，CAM)，從圖中可以直觀看出網絡更為關注的區(qū)域.本文選取SODN空間網絡的卷積層輸出結果計算CAM圖,利用該圖自動定位亞暴起始時刻亮斑區(qū)域,將定位的亮斑區(qū)域作為亞暴發(fā)生位置進行統(tǒng)計分析.雖然該定位區(qū)域可能無法與整個亮斑完全重合,但該區(qū)域是整個亮斑區(qū)域的子集,其統(tǒng)計結果仍可揭示亞暴發(fā)生位置的分布規(guī)律.定位過程如圖7所示,具體步驟如下:

圖7 空間網絡(BNInception V2)中Inception5生成CAM圖.從左到右的每一列依次表示MLT-MLAT圖像、生成的CAM圖、通過CAM圖粗略定位亮斑位置、截取出的圖像亮斑部分和最終定位到的亮斑區(qū)域.(a)圖中的UVI圖像發(fā)生于1998年12月04日12∶12∶49UT,(b)圖中的UVI圖像發(fā)生于1998年12月11日08∶13∶45UTFig.7 CAM map generated by Inception5 module in spatial network (BNInception V2).Each column from left to right represents the MLT-MLAT image,the generated CAM map,the rough spot location map through the CAM map,the extracted bright spot,and the final detected spot locations respectively.The UVI image in (a)was captured at 12∶12∶49UT on December 04,1998,and the UVI image in (b)was captured at 08∶13∶45UT on December 11,1998

(1)利用SODN時序檢測結果確定亞暴事件起始幀;

(2)將起始幀輸入空間網絡,獲得空間網絡最后一個Inception模塊的CAM圖;

(3)通過對CAM圖設置閾值(設置多組候選值篩選后得到)獲得圖中權重較大的區(qū)域，該區(qū)域大致定位亮斑所在位置;

(4)對定位出的區(qū)域繪制其灰度直方圖,由于亮斑區(qū)域相對其他區(qū)域占比較大且亮度高,在直方圖上呈現(xiàn)為一個峰值，將小于峰值對應的灰度值范圍灰度置零，剩下的非零區(qū)域就是亮斑所在位置.

利用上述方法對所有檢測出的亞暴事件統(tǒng)計分析其發(fā)生位置.圖8給出了亞暴發(fā)生位置在MLT-MLAT坐標下的分布情況.圖中亞暴發(fā)生位置集中在20—02 MLT、60°—70°MLAT之間、該結果與以往亞暴研究的物理結論高度一致(Liou,2010).

圖8 亞暴起始時刻亮斑在MLT-MLAT網格圖中的分布位置及發(fā)生頻率(%)Fig.8 The distribution and occurrence frequency (%)of the substorm onset locations in MLT-MLAT grid chart

4 結論

針對現(xiàn)有亞暴自動檢測方法的不足,本文提出一種基于深度學習的亞暴事件時-空自動檢測模型SODN.該模型利用兩個并行的卷積神經網絡分別提取亞暴的空間形態(tài)特征和時序運動特征,然后將兩種特征級聯(lián)生成亞暴特征表示.該特征可以有效表示亞暴從增長相到恢復相整個過程中極光卵在時空上的變化情況,解決了傳統(tǒng)方法需要分步提取亞暴特征的問題,簡化了特征提取流程，降低了模型復雜度,使得該模型在實際應用中擁有良好的部署性和擴展性.同時為了兼顧亞暴檢測的準確性和效率,SODN利用三個緊湊的時序卷積層融合多個時序范圍的亞暴特征,直接輸出亞暴起始的概率序列,克服了傳統(tǒng)自動檢測方法需要人工參與分析、手動設置閾值和檢測效率較低的問題.在多年的Polar衛(wèi)星觀測數(shù)據(jù)上SODN實現(xiàn)了端到端的亞暴事件自動檢測,檢測速度高達393幀/s;而且利用該方法對亞暴發(fā)生位置的統(tǒng)計結果與已有物理結論非常吻合,證明了SODN可以應用于大規(guī)模亞暴時-空檢測,具有良好的實用價值.

亞暴特征提取的優(yōu)劣會直接影響亞暴檢測結果,后續(xù)工作將嘗試采用最新的自監(jiān)督光流模型(Sun et al.,2018)改進光流效果,進一步提升檢測準確性.此外,數(shù)據(jù)集的質量直接影響著網絡的性能,本文僅考慮單個波段的數(shù)據(jù)，時間分辨率較低且部分數(shù)據(jù)存在標注錯誤,因此SODN在檢測性能上仍有較大的提升空間.后續(xù)工作將結合衛(wèi)星多個波段采集的UVI圖像提升數(shù)據(jù)集的時間分辨率,并將該模型應用于其他衛(wèi)星(如IMAGE衛(wèi)星)觀測數(shù)據(jù),進一步驗證模型的泛化性能.

致謝本文UVI圖像數(shù)據(jù)由美國國家空間科學和技術中心(National Space Science and Technology Center,NSSTC)提供,亞暴事件的人工標記由Liou提供.在此感謝上述機構和專家對本文工作的幫助.