基于改進粒子群深度神經網絡的頻率域航空電磁反演

廖曉龍，張志厚，范祥泰，路潤琪，姚禹，曹云勇，徐正宣,

(1.西南交通大學地球科學與環境工程學院，四川成都，611756；2.中鐵二院成都地勘巖土工程有限責任公司，四川成都，610000)

1954年，頻率域航空電磁法在加拿大應用于礦產勘探，并發揮了重要作用。經過60多年發展，該方法已廣泛應用于基礎地質調查、油氣資源勘探及環境工程等領域[1]。頻率域航空電磁的常用成像方法包括質心深度成像法[2]、微分電阻率成像法[3]等，此類方法可快速獲得地下典型參數，但成像分辨率較低。為獲得更加準確、定量化的地下電性信息，諸多學者對反演進行了研究，包括一維層狀介質反演[4-6]、擬二維橫向約束反演[7-9]和整體反演[10-11]等，這些研究極大推動了航空電磁反演的實用化，但此類方法依賴初始模型的選取，反演解容易陷入局部極小，對反演的準確性和穩定性有一定影響。近年來，完全非線性反演方法在地球物理數據處理解釋中受到高度重視，此類方法對初始模型依賴小，無需計算靈敏度矩陣，有著良好的全局尋優性能。朱凱光等[12-13]研究了時間域航空電磁的神經網絡反演；劉力輝[14]等采用深度神經網絡進行地震巖相反演；奚先等[15]實現了散射波場的深度學習反演成像；師學明等[16]研究了大地電磁阻尼粒子群優化反演；柳建新等[17]將遺傳算法引入大地電磁測深二維反演；胡祖志等[18-19]實現了大地電磁模擬退火反演及人工魚群反演；殷長春等[20]采用變維數貝葉斯方法對頻率域航空電磁資料進行反演。以上研究均取得了很好的成效，可見完全非線性算法在地球物理資料處理中發揮著重要作用，但這些算法的求解精度和收斂速度仍有提升空間。每種反演方法都有其優點及局限性，將不同反演方法的優勢結合起來，通常能達到更好的反演效果。高明亮等[21]將免疫遺傳算法與BP神經網絡結合，實現了二維高密度電法數據反演；陳東敬等[22]使用帶模擬退火的擬BP 神經網絡反演重力資料；戴前偉等[23]結合粒子群優化算法與BP神經網絡，進行了電阻率成像非線性反演；王鶴等[24]結合遺傳算法與BP 神經網絡，實現了大地電磁反演，取得了極好的效果。以上研究表明，混合優化算法能兼顧各方法優點，提高反演的精度與速度，保證反演的穩定性。基于此，本文將深度神經網絡出色的局部尋優性能與改進粒子群優化算法良好的全局尋優能力相結合，提出了基于改進粒子群深度神經網絡的頻率域航空電磁反演方法。首先，闡述頻率域航空電磁一維正演理論；其次，介紹改進粒子群優化算法優化深度神經網絡的過程；第三，利用該方法對理論模型進行反演測試；最后，對實測航空電磁數據進行反演，并將反演結果與Occam 反演結果進行對比，以檢驗此方法的有效性。

1 頻率域航空電磁1D正演理論

目前，頻率域航空電磁系統主要采用水平共面裝置[20]，水平層狀介質的回線源激發和接收模型如圖1所示。根據NABIGHIAN[25]，垂直磁偶源接收線圈的垂直磁場分量表達式為

其中：

式中：Hz為總場；m為發射磁矩；rTE為TE極化的反射數；h為飛機的飛行高度；J0為零階貝塞爾函數；r為接收線圈與發射線圈的距離；Yn為本征面波阻抗；Y?n為本征阻抗；μn與σn分別為第n層的磁導率和電導率，通常認為大地磁導率等于自由空間的磁導率，即μn=μ0；ω為發射線圈發射角頻率；kn為第n層的波速。

反演處理中所使用的是總場去除一次場并歸一化后的數據[20]：

式中：HZ0為一次場；d為歸一化垂直磁場分量。層狀介質電磁響應公式含有零階貝塞爾函數的廣義積分，為保證計算速度與精度，本文采用快速漢克爾變換數值積分進行計算[26]。

2 改進粒子群深度神經網絡反演

2.1 深度神經網絡

深度神經網絡(deep neural network，DNN)，是近年興起的一種學習能力極強的機器學習方法，在語音識別、圖像識別及自然語言處理等領域應用廣泛。它起源于傳統的人工神經網絡，具備更深的網絡深度，能夠融合各種非線性特征，具有強大的非線性學習和擬合能力[27]。

圖1 頻率域航空電磁層狀介質模型示意圖Fig.1 Schematic diagram of airborne EM layered medium in frequency domain

DNN 由多層自適應的非線性單元組成，每個神經層之間完全連接，通過逐層學習不斷減小模型預測誤差。圖2所示為DNN 示意圖，由圖2可見：輸入層節點數M由輸入向量維數決定；隱含層節點數目與輸入輸出節點個數聯系緊密；輸出層節點數N由輸出向量維數決定。

DNN的正向傳播通過眾多的神經元計算實現，工作過程可用式(6)和(7)來表示：

式中：xi為輸入層第i個節點的輸入；wji為輸入層第i個節點與隱含層第j個節點的連接權值；θj為隱含層第j個節點的閾值；yj為隱含層第j個節點的輸出；f為激活函數。常用有以下幾種。

DNN 的反向傳播過程即為網絡的反向傳播及自我修正過程，代價函數常有2種。

圖2 深度神經網絡示意圖Fig.2 Schematic diagram of DNN

式中：di是訓練數據的正確輸出；M是節點的個數。式(11)為均方誤差函數，式(12)為交叉熵函數。由于DNN 具有較多的隱含層，訓練難度較大，為使代價函數達到最小值而避免陷入局部極值，常用隨機梯度下降和Adam等優化方法[28]。

本文采用6 層全鏈接DNN 模型，代價函數選擇均方誤差函數，隱層激活函數使用ReLU 函數(能有效改善“梯度彌散”現象)，輸出層采用sigmoid激活函數。

2.2 粒子群優化算法

粒子群優化算法(particle swarm optimization，PSO)是KENNEDY 等[29-30]提出的一種智能優化算法，其基本思想源于對鳥類覓食的研究。假設在1個D 維搜索空間中，有m個粒子組成的集群X={X1,X2,…,Xm}，第i個粒子的位置為Xi=(Xi1,Xi2,…,XiD)，第i個粒子的速度為Vi=(Vi1,Vi2,…,ViD)。在粒子群優化算法迭代進化過程中，粒子會根據個體最優位置和群體最優位置不斷調整自身移動方向，再到達下一位置。個體最優位置稱為個體極值(Gbest)，群體最優位置稱為全局極值(Pbest)，每個粒子移動的迭代公式如下[31]：

式中：t為進化次數；w是慣性權重；c1和c2分別為學習因子，一般取常數為2；r1和r2為[0，1]之間的隨機數。慣性權重描述了上一代粒子速度對下一代粒子速度的影響，改變其大小可調節算法的全局尋優與局部尋優的能力。w較大有利于提高算法的全局尋優能力；w較小可增強算法的局部尋優性能。學習因子決定著算法的“自我認知”和“社會認知”能力，影響著種群靠近最優區域的能力。

2.3 改進粒子群優化算法優化深度神經網絡

邵洪濤等[32-33]研究了粒子群優化算法的改進措施，取得了較好的效果。本文提出一種振蕩衰減改進措施，即將慣性權值設置為

式中：ws和we分別為慣性權值的初值和終值；α為衰減因子；慣性權值w(t)隨著迭代次數t振蕩衰減。

圖3 慣性權重迭代曲線Fig.3 Inertia weight iteration curve

慣性權重迭代曲線如圖3所示，由圖3可見：當迭代剛開始時，粒子的慣性權重取值較大，粒子將以較大的飛行速度尋找到全局最優解的范圍，保證粒子群優化算法的全局搜索能力；在搜索的后期慣性權重取值較小，粒子的飛行速度較小，進而可以對最優解領域范圍進行精確搜索，保證算法的局部搜索能力。同時，在搜索過程中，慣性權重振蕩調整，可以不斷改變粒子的速度與位置，增強粒子之間的信息交互，避免慣性權重減小導致算法陷入局部最優。

在求解復雜的非線性問題時，深度神經網絡在設置初始權值時依賴于對樣本空間反復試驗，易產生網絡不穩、收斂速度慢以及陷入局部最優等問題。而改進的粒子群優化算法具有良好的全局尋優能力，將其與深度神經網絡結合，理論上可以尋找到求解問題的全局最優解。改進粒子群優化算法優化深度神經網絡基本流程如圖4所示。本文利用改進粒子群優化算法對深度神經網絡優化步驟如下。

圖4 改進粒子群優化深度神經網絡流程圖Fig.4 Flow chart of DNN optimized by IPSO

1)依據待解決問題，初始化深度神經網絡，設定網絡的輸入層、隱含層和輸出層神經元個數，并確定傳遞函數、學習函數和訓練函數。

2)根據深度神經網絡模型確定粒子的維數，使粒子每個維度分量對應神經網絡的連接權值與閾值，并隨機初始化粒子群。

3)計算種群中每個粒子的適應度fi。將每一輪最小適應度記為個體極值Gbest，將種群中所有粒子中最小適應度記為全局極值Pbest。

4)更新個體極值Gbest和全局極值Pbest。在每一輪迭代過程中，若粒子適應度fi小于Gbest，則Gbest=fi；若Pbest小于Gbest，則Pbest=Gbest。

5)根據式(13)和(14)更新粒子的速度與位移，采用式(15)慣性權重進行改進。按照粒子位置區間[Xmin，Xmax]和速度區間[Vmin，Vmax]規范超出限定范圍的粒子。若滿足迭代的終止條件，則進行下一步，否則返回步驟3)。

6)當達到迭代次數或目標精度時，選擇最佳粒子作為深度神經網絡的初始連接權值和閾值，然后進行網絡訓練，直到達到最大迭代次數或滿足精度要求，算法結束。

3 模型試驗

本文分別采用改進的粒子群深度神經網絡(IPSO-DNN)、標準粒子群深度神經網絡(PSODNN)和深度神經網絡(DNN)對層狀地電模型進行反演測試。設置參數如下：DNN 為6 層全鏈接方式(包含輸入、輸出及隱含層)，其中隱含層4 層，節點數隨訓練數據集改變；粒子群規模20，進化次數200，學習因子c1=c2=2.0，PSO-DNN 算法慣性權值取0.8，IPSO-DNN算法慣性權值見式(15)。

3.1 2層地電斷面模型

將正演得到的d值的實分量等距取20 個頻點作為網絡的輸入，則輸入層節點數為20；將2層地電斷面模型的電阻率及厚度作為網絡的輸出，則輸出層節點(ρ1，ρ2，h1)個數為3。為減少計算量，將電阻率與厚度值固定在一定范圍，ρ1和ρ2取值范圍為100～1 000 Ω?m，h1的取值范圍為100～1 000 m，每個參數等距選擇10個值，即ρ1與ρ2為100，200，…，1 000 Ω?m，h1為100，200，…，1 000 m。則模型總量為1 000組，將其中800組模型作為訓練樣本，200組模型作為測試樣本，并對數據進行歸一化處理。

利用IPSO-DNN，PSO-DNN 和DNN 算法分別對800組已知模型進行學習訓練，圖5所示為訓練過程中均方誤差(MSE)收斂曲線。由圖5可見：DNN，PSO-DNN 和IPSO-DNN 算法起始反演結果均方誤差分別為0.102，0.090 和0.079，說明改進粒子群優化算法降低了網絡的初始誤差；在訓練過程中，IPSO-DNN算法在前20輪收斂速度極快，明顯優于其余2 種算法，當訓練至100 輪以后，PSO-DNN和DNN算法收斂速度極慢，計算已陷入局部極小，而IPSO-DNN 算法避開了此問題，100輪以后反演結果均方誤差已小于0.010；當訓練至200輪時，DNN算法的反演結果均方誤差為0.038，PSO-DNN 算法的反演結果均方誤差為0.025，IPSO 算法的反演結果均方誤差為0.008，說明IPSO-DNN 算法在收斂速度及精度上均高于DNN和PSO-DNN算法。

對200組測試數據分別利用IPSO-DNN，PSODNN 與DNN 算法進行反演，并從中隨機抽取10組數據的反演結果如表1所示，反演結果相對誤差如表2所示。從表1和2可見：隨機抽取的10組數據中IPSO-DNN算法反演相對誤差明顯低于其余2種反演方法的相對誤差。為了更直觀地展示反演結果，抽取一組G 型與D 型地電斷面模型繪制電阻率-深度圖，如圖6所示。從圖6可見：相比于PSO-DNN和DNN算法，IPSO-DNN算法反演出的電阻率和深度更接近真實模型。

圖5 不同方法訓練2層模型反演結果的均方誤差收斂曲線Fig.5 MSE convergence curve of inversion for training two-layer model with different method

3.2 3層理論模型

采用3層水平模型進行研究，則網絡輸出層節點(ρ1，ρ2，ρ3，h1，h2)個數為5。令ρ1，ρ2和ρ3的取值范圍100～1 000 Ω·m，h1和h2的取值范圍為100～1 000 m，每個參數等間隔取10 個值，則模型總量為10×104。將其中80%(即8×104組)數據模型作為訓練樣本，其余20%(即2×104組)數據模型作為測試樣本。利用IPSO-DNN，PSO-DNN 和DNN 算法分別對8×104組已知模型進行學習訓練，圖7所示為訓練過程中均方誤差(MSE)收斂曲線。由圖7可見：DNN，PSO-DNN 和IPSO-DNN 算法起始均方誤差分別為0.093，0.060和0.040，IPSODNN算法降低了網絡的起始誤差；在訓練過程中，IPSO-DNN 算法前10 輪均方誤差就已下降到較低值，相比之下，DNN和PSO-DNN算法收斂速度較慢；當訓練至200 輪時，IPSO 算法的均方誤差為0.016，PSO-DNN 算法的均方誤差為0.024，DNN算法的均方誤差為0.029，說明在3 層地電斷面模型中，IPSO-DNN算法的收斂速度及精度也高于其余2種算法。

表1 2層模型10組數據不同方法反演結果對比Table 1 Comparison of inversion results of 10 sets of data for the two-layer model with different methods

表2 2層模型10組數據不同方法反演結果相對誤差對比Table 2 Comparison of relative errors of 10 set of data for the two-layer model with different methods %

圖6 2層理論模型不同方法反演結果對比Fig.6 Comparison of inversion results of the two-layer model with different methods

對2×104組測試數據分別利用IPSO-DNN，PSO-DNN 和DNN 算法進行反演，并從中抽取10組數據的IPSO-DNN反演結果如表3所示，相對誤差如表4所示。由表3和4 可知：IPSO-DNN 算法的反演結果與理論模型的平均相對誤差僅為4.53%，而PSO-DNN 和DNN 的平均相對誤差均超過了10%。從表3中分別抽取A，H，K 和Q 型地電斷面模型繪制電阻率-深度圖，如圖8所示，IPSO-DNN 算法的反演結果能更準確地接近真實模型。

3.3 不同程度干擾下的反演結果對比

為測試該算法在加入干擾后的適用性，本文在正演得到的測試集數據中加入了均值和方差均為零的正態分布隨機噪聲，含噪等級分別為3%，5%和10%。分別利用IPSO-DNN，PSO-DNN 和DNN 算法對3 種不同程度噪聲水平的數據進行反演，并從中隨機抽取了15 組反演結果，計算出與理論數據平均相對誤差，如表5所示，由表5可見：在2層及3層模型的含噪數據反演結果中，不論噪聲程度如何，IPSO-DNN算法反演結果相對誤差均低于PSO-DNN 和DNN 算法；此外，噪聲等級3%對各個方法反演效果的影響并不大。但隨著噪聲等級增加至5%和10%，PSO-DNN 和DNN 反演結果已出現較大誤差，IPSO-DNN反演誤差雖然也有增大，但誤差增加幅度比PSO-DNN 和DNN反演算法的小，表明IPSO-DNN算法對于含噪數據具有更高的反演精度。

表4 3層模型10組數據不同方法反演結果相對誤差對比Table 4 Comparison of relative errors of inversion results of 10 set of data for the three-layer model with different methods %

4 實測數據反演

利用IPSO-DNN方法對某地區近地表的實測航空電磁數據進行反演，驗證該方法的有效性。實測數據為HCP裝置測量的垂直磁場分量，共觀測5個頻點，分別為386，1 538，6 257，25 790 和100 264 Hz。根據該地區已知的地質資料及巖礦物性測定結果，確定了適合本地區的IPSO-DNN算法參數取值范圍，并進行學習訓練；再利用訓練好的參數對實測數據進行反演，與Occam 反演結果進行對比。

圖8 3層理論模型不同方法反演結果對比Fig.8 Comparison of inversion results of the three-layer model with different methods

表5 不同程度干擾下各方法反演相對誤差對比Table 5 Comparison of relative error inversion of various methods under different degrees of interference%

圖9所示為航空電磁實測數據Occam 反演結果。由圖9可知剖面地層的連續性較好，地下電性結果大體分為3 層：第1 層電阻率為8～18 Ω·m，厚度約為10 m，經實際勘查確定為砂泥層；第2層電阻率為3～8 Ω·m，厚度約為10 m，確定為黏土層；第3 層電阻率為8～18 Ω·m，確定為含砂質地層。

圖9 航空電磁實測數據Occam反演結果(據殷長春等[20])Fig.9 Inversion results of survey data with Occam's method

圖10 航空實測數據IPSO-DNN與Occam反演結果對比圖Fig.10 Comparison of inversion results of IPSO-DNN and Occam data from survey data

圖10所示為4 個測點的IPSO-DNN 與Occam反演結果對比圖，圖10中的虛線為測點A，B，C和D(距離分別為25.7，157.7，288.2 和466.1 m)進行IPSO-DNN 反演的結果，實線為Occam 反演結果。由圖10可見：2種方法反演結果變化趨勢基本吻合，IPSO-DNN 算法的有效性得以驗證。圖10(d)中測點D的第1和3層反演電阻率與Occam反演結果有一定差異，說明此方法的預測結果與實際工區結果仍存在誤差，需要考慮更多實際影響因素。

5 結論

1)將深度神經網絡與改進粒子群優化算法結合，提出了改進粒子群神經網絡反演方法，將該方法運用于頻率域航空電磁反演，理論和實際結果均證明該方法的可行性。

2)改進粒子群深度神經網絡算法結合了改進粒子群優化算法的良好的全局尋優性能與深度神經網絡出色的局部尋優能力，有效避免了神經網絡訓練學習訓練陷入局部極小等問題。相比于粒子群神經網絡和單一神經網絡，該方法具有更高精度與收斂速度；此外，改進粒子群深度神經網絡算法對含干擾數據的反演效果也優于粒子群深度神經網絡算法和單一的深度神經網絡算法。

3)對于實測數據，受計算量的影響，需要了解測區相關地質資料以確定網絡參數控制范圍；為了減小實際資料反演誤差，應考慮更多實際影響因素，建立符合實際地質情況的數據集。