基于探測過程建模的探地雷達多目標識別

高 翔 姬光榮 姬婷婷 王 群

(1.中國海洋大學信息科學與工程學院電子系,山東青島266100;2.國家深海基地管理中心,山東青島266061;3.中國人民解放軍總裝工程兵研究一所,江蘇無錫214035)

1.引 言

目前,全世界平均每年有2萬余人因戰后的未爆地雷致死或致殘,地雷的探測與排除技術在國際上始終受到關注[1]。在現代地雷中,反步兵地雷等淺層目標的體積小、金屬成分少,傳統的金屬探測器無能為力,而探地雷達(GPR)利用高頻電磁波的發射與接收實現地下無損探測,具有更好的普適性,因此,得到了廣泛應用[2]。

已有大量文獻對GPR數據的雜波抑制與地下目標定位進行了討論[3],本文在其基礎上專注于利用目標區域的回波信號對埋藏目標進行識別。GPR信號的目標識別包含兩方面關鍵內容：對回波信號進行有效的特征提取以及選取適合的判別算法對地下目標進行識別。在特征提取中,維格納-威利(WVD)分布及其擴展[4-5]、S變換[6]、小波變換[7]以及功率譜估計[8]等非平穩隨機信號處理方法被許多學者所關注,線性判別分析與主元分析等方法也有學者進行研究[7,9]。另一類思路針對B-scan圖像中目標形成的雙曲線,采用各種圖像處理手段描述目標[10-14]。相比較而言,時頻分析方法更貼近信號的本質特性,有利于進行分類識別。在目標識別中,神經網絡[7,15]較好的并行計算和學習能力以及支持向量機(SVM)[16-17]較強的小樣本推廣能力使它們在GPR信號分析中已有成功的應用。此外,隱馬爾可夫模型(HMM)也是報道較多的方法之一。目前國際上絕大部分HMM應用于GPR信號的討論來自于Paul D.Gader領導的研究小組[11-14],但他們的HMM檢測算法都以回波信號中雙曲線目標的形狀特征為基礎,并未涉及信號中豐富的時頻聯合信息。

在探地雷達掃描過程中,天線沿水平方向移動并采集數據,在某特定位置收到的一維時域信號反映了此點下方鉛直方向的地下介質分布狀況,一般將此一維時域信號稱為A掃描或A-scan;由若干A-scan并行組成的位置—深度二維數據稱為B掃描或B-scan;由若干平行的B-scan組合形成的三維數據稱為C-scan.A-scan是探地雷達回波信號的基本組成單元,而時頻分析正是研究A-scan這類時域非平穩隨機信號的優良工具。時頻分析靈活性強、能夠定量描述信號在時域和頻域的能量分布特征。考慮到時頻分析結果的二維特點,本文結合圖像的紋理特征分析,提出采用紋理描述算子提取時頻分析結果的特征,以便反映出A-scan序列的變化情況。在實際探測中,探地雷達天線與目標在空間上存在著“由遠及近再及遠”的變化過程,這種空間變化通過一系列A-scan組成的序列的變化得到表現。掃描過程中相鄰的A-scan之間具有較強的相關性,因此,可將掃描過程得到的A-scan序列作為一個整體來建模。HMM的優勢正在于其對變化過程進行建模的能力,即它能夠定量描述對象在不同狀態間轉換并表現出來的過程,這一點是神經網絡與SVM等僅利用特征點集進行分類的識別算法無法做到的。本文正是以HMM的這一特點為基礎,利用A-scan的時頻紋理描述算子構造特征量以反映探測回波的變化過程,結合探測實際采用無跨越單向連續HMM對變化過程建模,實現地下目標的識別。仿真實驗證實了基于HMM算法的分類性能優于SVM.

2.基于時頻分析紋理描述的特征提取

在時頻分析中,由于短時傅里葉和小波變換等線性方法無法描述信號的瞬時功率譜密度,因此,本文采用雙線性方法,以使結果更加直觀和合理;對于包含時間窗的Cohen類分布,時間窗無疑限制了應用的靈活性,因此,本文采用不含時間窗的WVD分布對A-scan信號進行時頻分析。同時,由于本文關注的是時頻圖像紋理間的區別而非信號內各分量間的干擾,因此,WVD分布的交叉項問題并不影響此處的應用。信號x(t)的WVD分布如式(1)所示

式中,z(t)是x(t)的解析信號。鑒于不同位置A-scan的WVD分布圖像僅在某些特定時窗與頻寬內顯著變化,在得到WVD分布圖像后采用文獻[4]所述的3σ局部化算法提高計算效率。

一維時域信號經WVD分布可得到二維時域-頻域聯合分布圖像,需在此基礎上提取其特征。本文結合圖像處理的有關內容,采用圖像紋理分析方法進行特征提取,而紋理分析中選擇合適的紋理描述算子構建特征向量是問題的關鍵。紋理描述算子的選取原則是使彼此具有盡可能小的相關性,以便全面反映信號特征。本文一方面借鑒了灰度共生矩陣(GLCM)[18]與紋理特征編碼方法(TFCM)[19],在得到紋理特征編號共生矩陣(TFNCM)后,選取式(2)～(5)所述的4個描述算子;另一方面部分地吸收了Savelyev等人[4]的方法,采用WVD分布圖像的第一奇異值和奇異向量作為描述算子,如式(6)～(7)所述。

1)能量指數或均勻度(EN)

2)四階簇變形指數(CP)

3)二階差分矩(SODM)

4)相關性指數(CO)

式中：pΔ(i,j|d,θ)表示當毗鄰像素灰度差異閾值為 Δ時,方向角為θ、距離為d的TFNCM 矩陣中的(i,j)元素;μi和 μj表示給定(d,θ)下 TFNCM 矩陣中的行方向均值與列方向均值。d與θ的組合表示八個與中心像素的毗鄰像素,紋理描述算子的最終結果取八個情況的均值即可。

5)第一奇異向量指數(FSVeI)

6)第一奇異值指數(FSVaI)

式(6)～(7)中：K與L分別為時頻分析圖像的時寬及頻寬 ;σ1,σ2,u1,2與 ν1,2分別為最大和次大奇異值及其對應的奇異向量。

3.基于HMM的目標識別

HMM起源于20世紀60年代后期,是一種描述隨機過程統計特性的參數化方法。HMM的顯著特點是雙重隨機性：馬爾可夫鏈用來描述隨機狀態轉移,而一般隨機過程用來描述某狀態的隨機觀察序列,兩類隨機問題均采用概率密度函數進行描述。模型中的狀態轉移過程不可見,固稱之為“隱”馬爾可夫模型。HMM 用以下參數表征：λ=(A,B,π),其中A為狀態轉移概率矩陣,B為觀察序列概率矩陣,π為模型初始狀態。

若狀態總數為N,狀態序列Q={qt},則A=[aij]且aij=P(qt=j|qt-1=i)其中i,j=1,…,N.在實際探測中,探測天線與地下目標之間存在著“由遠及近再由近及遠”的變化過程,如圖1所示。考慮到盡可能使模型簡化、減小計算量,本文采用無跨越單向HMM建模。根據上述距離變化的物理意義,可取隱狀態個數N=3,狀態集合為{“遠離”,“接近”,“遠離”},狀態間彼此互不相同且只能單向切換,狀態轉換如圖2所示。關于隱狀態個數N的確定,實驗部分將進一步討論。

本文的可觀測符號為連續值的時頻紋理特征向量,因此使用連續HMM建模。在探測范圍內,設觀察序列為正態分布且概率密度函數各不相同,對于任意位置的A-scan,其觀察序列概率可采用高斯混合模型(GMM)描述,采用一個或若干個高斯函數分量的線性相加來表征總概率。記觀察序列概率矩陣B=[bj(ot)],其中

式中：M為GMM中高斯分量的個數;ot為前述的6維特征向量;μj(m)為第j個狀態下第m個高斯分量的均值;為相應的協方差;ωj(m)為混合權值,且若每個觀察序列的長度為T,則有為6維高斯函數

模型初始狀態可取“遠離狀態”;A的各行元素均勻取值即可,且滿足;對于矩陣B,將訓練集觀察序列中的特征向量根據三類狀態進行3-均值聚類,計算每個聚類的均值、方差及歸一化權重,得到初始的高斯混合模型。

在訓練過程中,首先將數據歸一化,再采用Baum-Welch算法估計模型參數。考慮到實際計算中的下溢問題,加入比例因子項。記α為前向概率,β為后向概率,則

比例因子由下式計算

則前向概率計算公式改寫為

后向概率為

由此,狀態轉移概率矩陣A的重估公式為

式中：

式中,k表示第k個觀察序列,1≤k≤K.對矩陣B的重估轉化為對混合高斯概率密度函數的重估,三個參數的重估公式分別為

另外,在利用Viterbi算法選擇最佳狀態序列時,為了避免使用比例因子,采用對數形式對概率進行計算[20]。

在識別環節,在數據歸一化后,對于每個訓練得到的HMM采用Viterbi算法[20]計算其產生待識別序列的概率,取最大者作為識別結果。

4.實驗結果與分析

本文所用數據來自比利時皇家軍事學院超寬帶探地雷達實驗室[21],埋地目標為PMN地雷、VS/50地雷、PRB409地雷、石塊4類物體。數據共包含18組C-scan,分別在不同目標、不同土壤介質類型與不同埋藏深度的狀況下獲得,每組C-scan包含50×50組A-scan。本文對包含埋藏目標的17組C-scan進行訓練與識別。對于每一組C-scan,本文取目標周圍平行于掃描方向的26個B-scan(第15道至第40道)作為訓練與測試數據,如圖3所示,圖中深灰色為埋地目標,淺灰色區域為本實驗所用數據區,箭頭所指為掃描方向。

圖3 實驗數據空間分布示意圖

在基于HMM的識別過程中,對于掃描沿徑B-scan中的每一條A-scan,首先計算其WVD分布,再對結果進行時頻紋理特征提取,形成6維特征向量,將每一道B-scan轉化為長度為50的特征向量觀測序列,以26個觀測序列為基本單位對其相應類別進行訓練并識別。在得到特征向量后,本文采用文獻[17]所述基于SVM的識別方法作為對比,在交叉驗證選取參數后,先采用基于單通道的SVM識別方法進行計算,再將掃描沿徑上的每一道B-scan作為多通道輸入,以提高準確率。

考慮到訓練數據與測試數據的選取對識別結果影響較大,本文針對訓練充分和訓練不充分兩種情況進行實驗。實驗方案Ⅰ(訓練不充分情況)：將第15～30道B-scan作為訓練集,共16×17=272組觀測序列,分為17類,分別對SVM與HMM兩個模型進行訓練并測試;將第31～40道共10×7=170組觀測序列作為測試集進行測試。此方案中訓練集與測試集所用數據相對目標位置完全不同,兩組數據包含不同信息,可看作訓練不充分情況下的測試。實驗方案Ⅱ(訓練集充分情況)：將第15～40道中的奇數道作為訓練集,偶數道作為測試集,兩組數據各有13×17=221組觀測序列。此方案中訓練集與測試集所涉及數據位置相似且完全覆蓋目標,包含信息相似,可看作訓練充分情況下的測試。

按照上述兩個實驗方案,可得表1和表2所示結果,可以看出：訓練不充分方案的識別率明顯低于訓練充分方案;在利用訓練集生成模型后,模型對訓練集的識別率也必然高于測試集。

表1顯示了SVM與隱含狀態個數為3的HMM的識別結果對比,其中每個狀態GMM中的高斯分量個數M的取值不同。從表 1可知：1)HMM與SVM分類思想的不同帶來了識別性能上的差異。一般而言,SVM是一種判決模型即分類器,而HMM是一種生成模型即統計模型。在分類對象是無序集的情況下,通常分類器的識別率優于生成模型。但對于探地雷達回波序列這類前后具有相關性的有序信號而言,側重對序列變化過程建模的HMM更為有效。HMM將若干A-scan組成的序列作為一個整體,整體的各部分之間具有相關性,因此,序列這一整體降低了其中每個A-scan錯分的可能性;而SVM僅限于對A-scan組成的無序集合分類,割裂了原序列A-scan之間的相互關系,即便采用多通道統計方法進行修正,其性能較HMM也略遜一籌。2)在隱含狀態數N相同的情況下,GMM中的高斯分量個數越多,對模型的描述越精細,識別效果越好。

表1 SVM與HMM(隱狀態數N=3)識別結果對比

在HMM的實際應用中,隱狀態個數N與GMM中的高斯分量個數M的確定始終是討論的熱點問題。本文設隱狀態個數N=3,物理意義明確,但其合理性也需實驗驗證。表2即驗證了在GMM分量個數確定的情況下不同隱狀態個數對HMM識別結果的影響。在方案Ⅰ中的不充分訓練實驗中,隱狀態個數為3時識別效果最好,由此可以得出結論：當訓練集所提供的信息并不充分時,隱狀態個數為3的HMM具有最好的推廣能力。在大多數實際情況下,訓練集數據所提供的信息往往不可能充分,因此,前述假設在實際應用中具有合理性。當然,當訓練集信息較為充分時,模型本身的精細程度會隨著隱狀態個數的增加而增加,其識別性能也相應提升。

還需指出,上述實驗將2類土質、4類目標、3個不同深度的情況細化為17個類別,對此17類分別建立模型并進行識別,模型庫細致且全面,涵蓋了諸如對同一深度、同一土壤中的不同的目標進行識別等各種子問題。在實際應用中,可針對不同需要建立條件更為寬松的模型庫。例如,可將上述實驗中同一目標在不同深度、不同土質的數據所訓練得到的模型合并為單一模型庫,則此模型庫在埋設深度與埋設土質方面更具普適性而對埋設目標更具專一性。

在本文所述算法中,隱馬爾科夫模型的模型訓練過程計算量最大但可離線完成,模型的在線識別并不影響實時運算速度。對于計算A-scan信號的時頻聯合分布圖像及其紋理描述算子,可在前期加入適當的A-scan信號截斷及時頻圖像主要區域分割以減少計算量。另外,所采用時頻分析方法的抗噪能力很大程度上決定了本文算法的抗噪性能,可在實際應用中對時頻分析方法加以選擇。

表2 不同隱狀態個數HMM的識別結果對比(觀察序列概率中GMM分量個數M=1)

5.結 論

本文以A-scan為基本對象,提出采用EN、CP、SODM等六個紋理描述算子對時頻分析得到的時頻聯合圖像進行描述,直接反映A-scan信號的變化特征。在識別環節,考慮到實際掃描過程中天線與目標之間“由遠及近再及遠”的變化過程,采用連續HMM對不同目標的探測過程建模,實現了基于掃描過程的多目標識別。實驗對比分析表明：這種基于過程的目標識別方法比對無序特征集合進行分類的SVM方法具有更好的效果。

[1] FURUTA K,ISHIKAWA J.Anti-personnel Landmine Detection for Humanitarian Demining[M].1st Edition.London：Springer-Verlag London Limited,2009.

[2] DANIELS D J.Ground Penetrating Radar[M].2nd Edition.London：Institution of Electrical Engineers,2004.

[3] 吳仁彪,劉家學,張 蓓.探地雷達地雜波抑制方法研究進展[J].信號處理,2005,21(4A)：510-513.WU Renbiao,LIU Jiaxue,ZHANG Bei.Survey on GPR ground bounce removal[J].Signal Processing,2005,21(4A)：510-513.(in Chinese)

[4] SAVELYEV G T,KEMPEN VAN L,Sahli H,et al.Investigation of time-frequency features for GPR landmine discrimination[J].IEEE Trans.Geoscience and Remote Sensing,2007,45(1)：118-129.doi：10.1109/TGRS.2006.885077.

[5] WONG D C,NGUYEN L H,GAUNAURD G C,et al.Time-frequency analysis for radar classification of land-mineimages[J].Journal of Electronic Imaging,2007,16(3)：033014.doi：10.1117/1.2776354.

[6] GUHA S,KRUSE S and WANG P.Joint time-frequency analysis of GPR data over layered sequences[J].The Leading Edge,2008,27(11)：1454-1460.doi：10.1190/1.3011017

[7] 周輝林,田 茂,熊俊志.探地雷達回波信號的特征提取和分類[J].電波科學學報,2006,21(4)：586-589.ZHOU Huilin,IAN Mao,IONG Junzhi.Feature extraction and classification of ground penetrating radar echo signals[J].Chinese Journal of Radio Science,2006,21(4)：586-589.(in Chinese)

[8] 李建勛,鄭軍庭.超寬帶探地雷達自動目標識別研究[J].電波科學學報,2005,20(4)：535-540.LI Jianxun,ZHENG Junting.Automatic target recognition for ultra-wideband ground-penetrating radar[J].Chinese Journal of Radio Science,2005,20(4)：535-540.(in Chinese)

[9] 胡進峰,周正歐.基于核方法的探地雷達目標特征提取方法[J].電波科學學報,2005,20(5)：671-674.HU Jinfeng,ZHOU Zhengou.GPR target feature extraction based on kernel method[J].Chinese Journal of Radio Science,2005,20(5)：671-674.(in Chinese)

[10] JOONKI P,CHEOLHA P L,MONGIA.Imageprocessing-based mine detection techniques：A review[J].Subsurface Sensing Technologies and Applications,2002,3(3)：153-202.

[11] GADER P D,MYYSTKOWSKI M,ZHAO Y.Landmine detection with ground penetrating radar using hidden markov models[J].IEEE Trans.Geoscience and Remote Sensing,2001,39(6)：1231-1244.

[12] FRIGUI H,HO K C,GADER PD.Rea-ltime landmine detection with ground penetrating radar using discriminative and adaptive hidden markov models[J].EURASIP Journal on Applied Signal Processing,2005,2005(12)：1867-1885.doi：10.1155/ASP.2005.1867.

[13] ZHANG X P,GADER P D,HICHEM F.Feature learning for a hidden markov model approach to landmine detection[C]//Proc.SPIE,6553,2007,pp：655327.doi：10.1117/12.722593.

[14] FRIGUI H,MISSAOUI O,GADER PD.Landmine detection using discrete hidden markov models with gabor features[C]//Proc.SPIE,6553,2007,pp：65532A.doi：10.1117/12.722241.

[15] 王 群,何云龍,王春和,等.基于神經網絡的探地雷達探雷研究[J].電波科學學報,2001,16(3)：398-403.WANG Qun,HE Yunlong,WANG Chunhe,et al.Study of mine detection with GPR based on neural networks[J].Chinese Journal of Radio Science,2001,16(3)：398-403.(in Chinese)

[16] ZHANG J,LIU Q,NATH B.Landmine featureextraction and classification of GPR data based on SVM method[C]//ISNN 2004,LNCS 3173,2004,636-641.doi：10.1007/978-3-540-28647-9ˉ 104

[17] 胡進峰,周正歐,孔令講.探地雷達多目標識別方法的研究[J].電子與信息學報,2006,28(1)：26-30.HU Jinfeng,ZHOU Zhenou,KONG Lingjiang.Researches on GPR multi-object recognition[J].Journal of Electronics&Information Technology,2006,28(1)：26-30.(in Chinese)

[18] HARALICK R M,SHANMUGAM K,DINSTEIN I.Textural features for image classification[J].IEEE Trans.Systems,Man and Cybernetics,1973,SMC-3(6)：610-621

[19] HORNG M-H.Texture feature coding method for texture classification [J].Optical Engineering,2003,42(1)：228-238.doi：10.1117/1.1527932.

[20] 王炳錫,屈 丹,彭 煊,等.實用語音識別基礎[M].北京：國防工業出版社,2005：188-190.

[21] Royal Military Academy of Belgian.HUDEM Project,UWB GPRdata[DB/OL].Ispra：DG Joint Research Centre,1999-12-16[2010.7].http：//apl-da-tabase.jrc.it/Home/Radar/Non-JRC/RMAˉ UWBˉGPR/.