利用磁梯度張量識別地下小型目標(biāo)體

鄭建擁 范紅波 張 琪 李志寧

(①陸軍工程大學(xué)石家莊校區(qū)，河北石家莊 050003； ②中國人民解放軍94019部隊，新疆和田 848000 )

0 引言

利用地磁場異常探測地下或水下物體，一直是勘探地球物理界的研究熱點。目前地磁技術(shù)研究的主要對象是地下或水下大型異常體(如潛艇)，而對諸如管道、未爆彈藥等小型鐵磁性物體，始終缺乏成熟的技術(shù)進行探測。當(dāng)前，磁梯度張量方法受到廣泛關(guān)注，成為探測小型磁性目標(biāo)體形狀、大小、姿態(tài)等信息的新手段。

孫刃[1]較系統(tǒng)地分析總結(jié)了利用磁異常正反演識別目標(biāo)的方法，并將其用于水下沉船、管線、未爆軍火等目標(biāo)的探測；吳國超[2]提出利用張量矩陣特征值和6個矩陣元素對大型地質(zhì)體進行邊界識別，分辨率較高，試驗效果很好；王林飛等[3]利用磁梯度張量研究大型磁源體的位置與邊界；李金朋等[4]研究了地下小目標(biāo)體的反演及三維成像；尹剛等[5]利用張量不變量和矩陣特征值定位磁目標(biāo)體；朱慧慧等[6]建立了基于磁偶極子構(gòu)造法幾何建模的地下管線磁異常正演模型，分析了多種磁化因素對鐵管磁異常正演結(jié)果的影響，為高精度地下管線磁探測奠定了理論基礎(chǔ)。上述方法主要針對磁性目標(biāo)的二維邊界識別、三維姿態(tài)反演、磁性目標(biāo)定位等，是局部磁異常識別方法，而且無法避免磁測數(shù)據(jù)的正、反演，運算量大，受測量數(shù)據(jù)精度的限制較大，因而不能清晰地表現(xiàn)磁測目標(biāo)體。謝永茂[7]提出了建立樣本庫，通過模板匹配識別地下磁性體，提供了一種新的磁異常體識別思路。

本文首先建立地下小型磁性目標(biāo)體不同形狀、不同方位條件下模型的正演數(shù)據(jù)庫，分析了9個磁梯度張量矩陣屬性參量，論證其作為分類屬性的可行性； 利用這9個矩陣屬性參量融合識別磁源目標(biāo)體的形狀，可在一定程度上彌補磁測數(shù)據(jù)精度的不足； 再將這9個屬性參量用于支持向量機(Support Vector Machine，SVM)對目標(biāo)進行分類； 然后利用量子粒子群算法優(yōu)化的支持向量機對仿真模型數(shù)據(jù)分類，推斷目標(biāo)體形狀，取得了很好的識別效果； 最后設(shè)計實驗驗證其實用性。

1 磁梯度張量矩陣及屬性參量

1.1 磁梯度張量矩陣

地磁總場包括地磁場和異常場，其中地磁場在空間上的變化很小，在小范圍內(nèi)可以看做一個穩(wěn)定場，所以磁梯度張量可以認為是磁異常矢量在x、y、z三個方向上的變化率。梯度張量元素受地磁場傾向、偏角影響較小，通過反演能夠更精確地描述場源體的磁化方向和幾何形態(tài)[8]。

地磁場是矢量場，由矢量場三分量推導(dǎo)梯度張量矩陣G表達式為

(1)

式中U為磁標(biāo)勢，矩陣中的9個元素即為磁場矢量B在空間三個方向分量(Bx、By、Bz)的變化率。磁法勘探中地磁場及鐵磁物質(zhì)產(chǎn)生的異常場可看作無源的靜磁場[5]。由式(1)可知梯度張量矩陣是對稱的，且旋度和散度均為0，因此梯度張量矩陣G的9個元素中，只有5個(Bxx、Byy、Bxy、Bxz、Byz)相互獨立[9]。

1.2 特征提取

作為一種機器學(xué)習(xí)方法，支持向量機需對分類對象進行特征提取，通過對特征進行分類來識別目標(biāo)。每一個特征都與分類對象本身相關(guān)，且不同對象的特征之間都應(yīng)該保持離散。特征的數(shù)量越多，分類越準(zhǔn)確。為保證信息不冗余，不同的特征之間必須互不相關(guān)，或相關(guān)度很低，這些特征即可作為支持向量機分類的屬性參量。從樣本庫中提取不同屬性參量的數(shù)據(jù)構(gòu)造特征向量，作為支持向量機的訓(xùn)練和測試樣本數(shù)據(jù)。

選擇磁總場模量(TMI)、矩陣特征值的最大值(ME)、兩個矩陣不變量(I1、I2)作為分類的屬性參量。

TMI可由磁場三分量計算

(2)

設(shè)A為n階矩陣，λ為常數(shù)，X為n維非零向量。若AX=λX，則λ稱為矩陣A的特征值。因為式(1)是3階張量矩陣，因此最多只有3個特征值。其中最大特征值對磁目標(biāo)體邊界和形狀都很敏感，識別能力較強。

矩陣不變量能夠表征一個矩陣的特征，經(jīng)由梯度張量矩陣計算得到的張量不變量不需要額外的處理就可以很好地描述磁場源。該矩陣的不變量蘊含著很多磁目標(biāo)體的特征信息，可以作為地下磁目標(biāo)體識別的分類屬性參量。式(1)中矩陣的特征值即是基本不變量[10]，矩陣中的其他不變量可由矩陣特征值表示。

定義λ1、λ2、λ3為張量矩陣G的特征值，則其余兩個矩陣不變量為

(3)

除以上4個屬性參量之外，還提取了式(1)中的Bxx、Bxy、Bxz、Byy、Byz等5個獨立元素作為支持向量機分類的屬性參量，總共使用這9個屬性參量識別目標(biāo)。

2 基于量子粒子群優(yōu)化的支持向量機

支持向量機可用于模式分類和非線性回歸，具有結(jié)構(gòu)風(fēng)險最小化的優(yōu)勢，具有很好的泛化性能。支持向量機的算法公式為

(4)

加入松弛因子ξ(i)，用于衡量實際指示值Li與支持向量機輸出之間的距離，支持向量機優(yōu)化的目標(biāo)函數(shù)為[11]

Li(wTXi+b)≥1，i=1,2,…,l

(5)

式中：w即是對目標(biāo)進行分類的超平面的法向量；l為分類時特征向量的維數(shù)；c是懲罰因子。

影響SVM的因素非常多，如懲罰因子c以及核函數(shù)參數(shù)σ，對于小型樣本分類精準(zhǔn)度的影響最大[12]。因此，在應(yīng)用SVM做分類預(yù)測時，這兩個參數(shù)的選取尤為關(guān)鍵。

一般通過網(wǎng)格劃分來尋找參數(shù)c和σ的最佳值。雖然采用網(wǎng)格搜索能夠獲得最高的分類準(zhǔn)確率，即全局最優(yōu)解，但如果要在更大范圍內(nèi)尋找參數(shù)c和σ的最佳值會很費時。本文應(yīng)用一種通過尋找最優(yōu)解來求解最佳參數(shù)的算法，即量子粒子群優(yōu)化算法(Quantum Particle Swarm Optimization，QPSO)[11]。

粒子群優(yōu)化算法是智能計算領(lǐng)域中，除蟻群算法之外的另一種基于群體智能的優(yōu)化算法。該方法通過粒子在解空間的搜索來尋找最優(yōu)解。相比基因等其他算法，粒子群優(yōu)化算法無需交叉、變異、選擇等流程，正確率比較高，全局最優(yōu)值收斂速度快，參數(shù)尋優(yōu)耗時短。

Sun等[13]認為粒子具有量子性，從量子力學(xué)的角度提出了一種新的粒子群優(yōu)化算法模型，并經(jīng)實驗證明其具有更好的收斂性。粒子群優(yōu)化算法[14]具體過程如下。

在有M個粒子進行尋優(yōu)的J維搜索空間中，假設(shè)第i個粒子的位置為Xi=(Xi,1，Xi,2，…，Xi,J)，這里t為進化代數(shù)； 其每代最優(yōu)位置為Pi=(pi,1，pi,2，…，pi,J)，整個量子粒子群目前所搜索到的全局最優(yōu)解為Pg=(pg,1，pg,2，…，pg,J)。通過200次迭代，選取全局最優(yōu)的c和σ。

粒子i的位置為

(6)

最后得到進化t代之后的pi,j，即可得到最優(yōu)參數(shù)c和σ。基于量子粒子群算法優(yōu)化后的c和σ，可建立最優(yōu)化QPSO-SVM分類模型。

3 仿真分析

3.1 仿真數(shù)據(jù)庫的建立

為了進行仿真驗證，首先構(gòu)建一個垂直磁化的磁性體進行正演。在觀測面上建立以觀測原點為中心的61×61個網(wǎng)格節(jié)點的數(shù)據(jù)平面。數(shù)據(jù)平面上每個網(wǎng)格節(jié)點都包含該位置的梯度張量信息，即一個由式(1)計算的3階梯度張量矩陣。

設(shè)置8組長方體磁性體正演模型，模型中心位置均處于坐標(biāo)原點，埋深均為5m，磁化率相同，姿態(tài)都為水平，但他們的形狀和磁化方向不同，即構(gòu)建8種不同尺寸的長方體模型。每個模型的磁化偏角范圍為0°～180°，間隔設(shè)置為10°(0°與180°重合)，這樣得到共180種不同磁化方向組合，總共8×9×180個包含61×61個網(wǎng)格點的數(shù)據(jù)平面。分別計算平面上61×61個點的TMI、3個不變量和5個張量矩陣元素，這樣每個平面都可以衍生出9個61維的獨立屬性矩陣。具體的數(shù)據(jù)結(jié)構(gòu)參見圖1。

圖1為9個屬性參量的等值線圖，每幅圖都代表一個屬性參量的屬性矩陣。圖中每個像素點的數(shù)值就是矩陣中相應(yīng)位置的值。

3.2 特征提取

為論證9個屬性參量作為支持向量機分類的屬性參量的可行性，以長、寬、高均為10m、埋深為15m的正方體模型為例進行正演。設(shè)磁傾角為90°，磁偏角為0°，磁化率為0.4。

由圖1可見，TMI的等值線分布(圖1a)比較均勻，但磁化體邊角特征不明顯；磁性體矩陣的最大特征值ME(圖1b)和I1(圖1c)相較于其他屬性參量，能更詳細和清晰地描述磁性體的邊界和形狀特征；I2(圖1d)對邊角特征更敏感，能凸顯淺部源； 5個矩陣獨立元素(圖1e～圖1d)中，Bxx描述y軸方向的邊界情況，出現(xiàn)了三個極值點，Bxy可描述磁性體的邊角，Bxz可描述y軸方向的邊界情況，出現(xiàn)了兩個極值點，Byy可描述x軸方向的邊界情況，出現(xiàn)了三個極值點，Byz主要描述x軸方向的邊界，出現(xiàn)了兩個極值點[15]。可見這5個矩陣獨立元素也能較清晰地表現(xiàn)磁性體的部分特征。可以看出，這5個屬性參量都與磁性體的形狀有著較高的相關(guān)性，可從不同角度反映磁性體的形狀特征。

圖2為所有樣本的9個屬性參量的分維屬性圖，可據(jù)此分析這9個獨立屬性參量對不同形狀和傾角的敏感性。圖中每個點線代表不同的磁偏角和磁傾角。對于傾角的敏感性，可根據(jù)每條點線的傾斜角度進行分析，因為每個點都代表不同的傾角，而這些點的排列是有規(guī)律的，因此對傾角變化敏感性強的屬性參量就是可用于研究磁化體形狀和姿態(tài)特征的分類。

圖1 9個屬性參量的等值線圖

不失一般性，可以認為這9個屬性參量與磁源的形狀、磁化方向等信息有著很高的相關(guān)性。由圖2可見，不同形狀時正演結(jié)果變化明顯，如基于屬性參量TMI(圖2a)與屬性參量I2(圖2d)可以很明顯地區(qū)分出8種形態(tài)；而屬性參量I1(圖2c)和Bxx(圖2e)對于磁化角度的敏感性隨著磁偏角的變化明顯不同，分別出現(xiàn)遞增和遞減兩種變化趨勢。從圖2還可以看出，分類屬性值與偏角呈明顯的線性關(guān)系。屬性參量Bxx與Byy、Bxz與Byz形態(tài)雖然相似，但數(shù)值大小相差兩個數(shù)量級。由式(1)可知，Bxx與Byy、Bxz與Byz相互獨立，因此它們之間具有離散性，即各屬性參量之間的數(shù)值相差很大，遠遠大于屬性參量內(nèi)部的數(shù)值差，不同參量的數(shù)值在數(shù)域上明顯分別開。因此，可以參與分類。

上述分析表明，這9個屬性參量對磁性體有良好的辨識能力，能夠作為支持向量機分類的屬性參量。每個屬性參量的樣本數(shù)據(jù)都可以構(gòu)造一個特征向量，作為支持向量機分類的輸入向量。

圖2 9個屬性參量的分維屬性

從樣本庫中調(diào)取數(shù)據(jù)，計算每個獨立屬性矩陣的最大特征值，分別對應(yīng)9個特征向量中的某一個元素。對于這8種形狀的磁化體，每種磁化體共有180個屬性參量平面，即可構(gòu)造出8組9×180的特征矩陣。以此作為PSO-SVM的訓(xùn)練樣本，并隨機抽選每種形狀下各90組數(shù)據(jù)為訓(xùn)練樣本，另外90組數(shù)據(jù)可作為測試樣本。

3.3 識別測試

在支持向量機中， 8種長方體形狀分別對應(yīng)的輸出代碼為整數(shù)1～8，即類別標(biāo)簽設(shè)置為1～8，分類結(jié)果通過數(shù)字體現(xiàn)。

在仿真數(shù)據(jù)輸入前，必須對數(shù)據(jù)進行歸一化預(yù)處理

(7)

式中p、q∈R。歸一化的目的是將原始數(shù)據(jù)規(guī)一化到[0，1]范圍。結(jié)果輸出后再進行反歸一化，最后經(jīng)過處理獲得1～8的分類標(biāo)簽。

設(shè)置量子粒子群優(yōu)化的初始搜索參數(shù)為：最大進化數(shù)量為200，種群最大數(shù)量為20。經(jīng)過200次進化尋優(yōu)，可以得到參數(shù)c和σ的最優(yōu)解分別為77.4375和6.2357。

適應(yīng)度用于評價粒子尋找最優(yōu)解的能力強弱，適應(yīng)度越大粒子尋優(yōu)能力越強，反之適應(yīng)度越小。適應(yīng)度度量函數(shù)因具體情況而異。適應(yīng)度的度量公式為

(8)

圖3為最終適應(yīng)度曲線，其中最佳適應(yīng)性曲線為每代粒子中的最佳適應(yīng)度，平均適應(yīng)度曲線為每代粒子的平均適應(yīng)度，其最佳適應(yīng)度最終為87.36%。

圖3 最終適應(yīng)度曲線

利用得到的c和σ的最優(yōu)解構(gòu)建基于QPSO-SVM的支持向量機。輸入720個訓(xùn)練樣本，對支持向量機進行訓(xùn)練，最后用測試樣本進行測試，結(jié)果見圖4，識別正確率約為93%。

圖4 測試分類與實際分類對比圖

圖5 實驗系統(tǒng)設(shè)計

4 實驗驗證

為驗證本文方法的可行性，采用4個磁通門傳感器構(gòu)建了平面十字形磁梯度張量系統(tǒng)。圖5為實驗所使用的實驗臺，圖6為實測結(jié)果。

實驗采用直徑分別為10、20cm、對應(yīng)長度為100、70cm的兩種鐵圓柱和一個半徑為20cm、厚度為10cm的鐵圓板三種樣品作為識別對象。分類標(biāo)簽分別設(shè)置為1、2、3。設(shè)置水平埋深為2m。搭建邊長(x，y軸)為2m的水平正方形實驗測量系統(tǒng)，x和y軸方向每隔0.1m設(shè)置一測量點，測量三種樣品在同一空間位置時的正方形實驗面中20×20個測量點的磁場數(shù)據(jù)，因此得到一個20×20的矩陣。對三種樣品的4種姿態(tài)(與x軸分別成0°、45°、90°和135°)各采集了5組磁異常數(shù)據(jù)，并將樣本數(shù)據(jù)以3∶2的比例隨機分成訓(xùn)練樣本和測試樣本，建立磁異常數(shù)據(jù)庫。實驗測量前，先測量背景場，對儀器進行校正。

圖6 實測磁異常場分布

以長度為1m，直徑為10cm的鐵質(zhì)圓柱體為例分析結(jié)果。將其放置在實驗臺的正中位置，與實驗臺x軸夾角為45°，測量其地磁總場后減去背景場，即可得到異常場。由圖6可見，實測數(shù)據(jù)包含很強的干擾和噪聲，主要來自大地、實驗儀器和其他地下大型鐵磁性物體，干擾了對目標(biāo)體方位、輪廓、邊緣等細節(jié)的分析。實驗地點處于中緯度地區(qū)，磁傾角約為55°，磁偏角約為-9°，斜磁化使磁異常信號畸變和偏移，提高了信號分析的難度。但是根據(jù)圖6的磁異常分布仍然可以較清晰地分辨出目標(biāo)的方位，也可大致辨別異常體的部分輪廓信息，且表現(xiàn)出與圖1理論仿真信號相似的形態(tài)特征。

圖6中屬性參量TMI、I1、I2和ME均在一定程度上包含目標(biāo)體的形狀和姿態(tài)信息。對比圖1的仿真結(jié)果可以看出，這4個參量的分布形態(tài)雖然發(fā)生了偏移和扭曲，但是等值線的形態(tài)仍與目標(biāo)體保持了一定的相似性，可反映出目標(biāo)的部分輪廓。而另外的5個張量矩陣獨立元素受到斜磁化影響較大，不能直觀地反映目標(biāo)形態(tài)，但是等值線峰值的連線卻可大致反映出目標(biāo)體的姿態(tài)。以上9個屬性參量分布圖表明，通過實測磁異常信號并不能直觀地識別目標(biāo)，但是實際信號與目標(biāo)之間仍保持一定的相關(guān)性。因此，不失一般性地，認為這9個屬性參量可以用于實際目標(biāo)分類。這就是本方法的基本思想，對于不能直觀反映目標(biāo)信息的實測磁場，利用支持向量機對信號的特征進行分類，以達到識別目標(biāo)的目的。

最后采用PSO -SVM進行分類，結(jié)果見圖7，識別正確率約為87.5%。

由實驗可知，本方法對小型鐵磁性物體的識別取得了一定的效果。由于地磁背景場比較強，異常體或目標(biāo)體的磁異常信號易被背景磁場淹沒，降低了正確率。如圖7所示，識別錯誤主要集中在標(biāo)簽1和標(biāo)簽2，有兩個標(biāo)簽1被識別為2，一個標(biāo)簽2被識別為1，這是因為這兩種識別對象的形狀相似度較高。因此關(guān)于細節(jié)識別有待進一步的研究。

圖7 實驗分類結(jié)果

5 結(jié)論

(1)利用矩陣最大特征值，構(gòu)造模式識別的特征向量，分析并提取張量矩陣的9個屬性參量聯(lián)合識別磁性目標(biāo)體，一定程度上彌補了磁測數(shù)據(jù)精度不足的缺點，也避免了大量磁測數(shù)據(jù)的計算。

(2)將機器學(xué)習(xí)的方法引入地下磁性目標(biāo)體的識別。利用量子粒子群改進的支持向量機識別地下小型目標(biāo)體的形狀，取得了較好的效果。

(3)本文僅研究磁性異常體形狀的識別，但此方法還可以識別磁目標(biāo)的姿態(tài)、深度等。由實驗結(jié)果可知，本文方法對磁異常信號的測量和特征分析仍有很大的研究空間。