基于粒子群優(yōu)化算法和F分布的微震源三步定位方法

龐聰，李查瑋*，馬武剛，程誠(chéng)，江勇，廖成旺

(1.中國(guó)地震局地震研究所，湖北 武漢 430071；2.地震預(yù)警湖北省重點(diǎn)實(shí)驗(yàn)室，湖北 武漢 430071；3.運(yùn)城學(xué)院數(shù)學(xué)與信息技術(shù)學(xué)院，山西 運(yùn)城 044000)

0 引 言

微震監(jiān)測(cè)技術(shù)在礦山巖體破裂、頁(yè)巖氣勘探、地震預(yù)警預(yù)報(bào)等防災(zāi)減災(zāi)領(lǐng)域[1-3]發(fā)揮著重要作用，而微震源定位精度的提高或算法穩(wěn)定性的改善是微震監(jiān)測(cè)理論研究的重要方向。

微震源定位算法研究已經(jīng)有許多顯著成果，例如傳統(tǒng)微震定位算法Geiger法、線性定位法、單純形法等，后人又在此基礎(chǔ)上結(jié)合智能算法、最優(yōu)化思想等創(chuàng)造出一系列有效的微震定位方法。例如，李邵紅等[4]為了消除微震定位中已知參數(shù)帶來(lái)的誤差與錯(cuò)誤，假定全部檢波器參數(shù)準(zhǔn)確，建立微震定位多目標(biāo)模型，應(yīng)用多目標(biāo)遺傳算法(NSGA-Ⅱ)進(jìn)行二次反演，剔除第一次反演獲得的定位異常樣本值，第二次反演可以獲得較為準(zhǔn)確的定位結(jié)果；朱權(quán)潔等[5]在四四組合法定位結(jié)果基礎(chǔ)上，利用聚類算法逐步過(guò)濾掉異常值，并將聚類中心曲線拐點(diǎn)與異常值離散度作為聚類持續(xù)與否的準(zhǔn)則；丁恩杰等[6]創(chuàng)新地應(yīng)用TDOA定位原理確立微震定位數(shù)學(xué)模型，該方法減少了速度因素帶來(lái)的誤差影響，但是研究?jī)?nèi)容僅局限于若干監(jiān)測(cè)點(diǎn)布置在單維直線上的情況，代表性不足，且微震源定位目標(biāo)函數(shù)存在主觀簡(jiǎn)化的現(xiàn)象，沒(méi)有給出相應(yīng)簡(jiǎn)化的理由。然后，上述微震定位方法雖然在異常值剔除、微震源定位模型、算法判定準(zhǔn)則等方面有一定的改進(jìn)和探索，但是也存在算法復(fù)雜度過(guò)高、運(yùn)算速度過(guò)慢、方法通用性差、穩(wěn)健性不好等問(wèn)題。

傳統(tǒng)微震定位算法的研究主要集中于精度的提高，忽視了安全生產(chǎn)管理、地質(zhì)活動(dòng)監(jiān)測(cè)中時(shí)間的重要性，以及某些智能算法尋優(yōu)結(jié)果隨機(jī)性大、穩(wěn)健性差的現(xiàn)狀。對(duì)于微震源定位結(jié)果的評(píng)價(jià)，主要依據(jù)誤差精度、算法執(zhí)行效率、算法穩(wěn)定性等方面。陳炳瑞等[7]發(fā)展了一種基于相鄰點(diǎn)到時(shí)差目標(biāo)方程和粒子群算法的微震源分層定位方法，利用礦山發(fā)震現(xiàn)狀反復(fù)校正定位結(jié)果，為了得到準(zhǔn)確、可靠的定位結(jié)果，計(jì)算過(guò)程過(guò)于復(fù)雜；李楠等[8]利用L1范數(shù)統(tǒng)計(jì)原理建立時(shí)間殘差公式，將監(jiān)測(cè)數(shù)據(jù)的中位數(shù)作為L(zhǎng)1范數(shù)統(tǒng)計(jì)的最佳估計(jì)，解決傳統(tǒng)微震源定位方法抗干擾性差、誤差較大等問(wèn)題。

針對(duì)微震震源定位研究存在的上述不足，本文嘗試發(fā)展一種結(jié)果穩(wěn)定、解算較快、精度更高的微震定位方法。在監(jiān)測(cè)點(diǎn)到時(shí)數(shù)據(jù)準(zhǔn)確的前提下，建立微震到時(shí)差目標(biāo)函數(shù)，為了獲得準(zhǔn)確初值和算法迭代停止閾值，設(shè)計(jì)由線性定位法、Geiger法、經(jīng)過(guò)優(yōu)化后的粒子群算法(Particle Swarm Optimization，PSO)等組成的三步反演方法，并將遺傳算法、非線性最小二乘算法、模擬退火算法、標(biāo)準(zhǔn)粒子群算法等作為定位方法參考對(duì)象，并討論P(yáng)SO隨機(jī)粒子分布對(duì)定位效果的影響。

1 關(guān)鍵技術(shù)

1.1 Geiger法與線性法

Geiger法具體步驟：先計(jì)算每個(gè)監(jiān)測(cè)點(diǎn)的走時(shí)、震源距離，再計(jì)算全體監(jiān)測(cè)點(diǎn)的偏微分矩陣，及時(shí)間差值與震源參數(shù)殘值，最終得到修正后的微震數(shù)據(jù)。

1.2 粒子群算法

粒子群算法(PSO)：它是一種常用的進(jìn)步型群體智能優(yōu)化算法，從隨機(jī)初值Pr(xr,yr,zr,tr)出發(fā)，通過(guò)更新粒子的飛行方向和飛行速度來(lái)尋找最優(yōu)解，然后通過(guò)適應(yīng)度Fitness(優(yōu)化目標(biāo)函數(shù))來(lái)評(píng)價(jià)解的質(zhì)量，從而獲得最優(yōu)近似解Ppso(x',y',z',t')。

粒子群算法核心是粒子的位置更新公式[7]:

Vid=w·Vid+c1r1(pid-xid)+c2r2(pgd-xid)

(1)

xid=xid+Vid

(2)

式中i=1,2,…,m；d=1,2,…,D；Vid為粒子在第d維的當(dāng)前速度；w為速度更新慣性權(quán)值；學(xué)習(xí)因子c1和c2為(0,2)區(qū)間內(nèi)的常數(shù)；r1和r2為[0,1]范圍內(nèi)的隨機(jī)數(shù)；pid為當(dāng)前單個(gè)粒子目前為止搜索到的最優(yōu)位置；pgd為當(dāng)前m個(gè)粒子所在粒子群搜索到的最優(yōu)位置[7]。

1.3 隨機(jī)粒子對(duì)粒子群的影響

在實(shí)際礦山開(kāi)采中，開(kāi)采環(huán)境較復(fù)雜，各種人類、機(jī)械的非正?；顒?dòng)或微震檢波器布置的局限性，都會(huì)導(dǎo)致監(jiān)測(cè)系統(tǒng)的隨機(jī)誤差無(wú)法遵循正態(tài)分布。粒子群算法的隨機(jī)粒子包含隨機(jī)位置和隨機(jī)速度，初值的隨機(jī)性與分布情況既影響到算法局部搜索與全局搜索的能力，也直接影響算法執(zhí)行效率。

傳統(tǒng)粒子群算法的隨機(jī)粒子分布情況遵循正態(tài)分布、均勻分布，初值的分布情況較為規(guī)律和簡(jiǎn)單，一定程度上增大了定位算法尋優(yōu)的難度。

(1)正態(tài)分布(N(μ,σ2))：該分布在自然界中較為常見(jiàn)，且具備對(duì)稱性、集中性、曲線面積一定、標(biāo)準(zhǔn)的數(shù)學(xué)描述方法等明顯特征，因而在工程計(jì)算中得到廣泛應(yīng)用。粒子群算法隨機(jī)解常采用標(biāo)準(zhǔn)正態(tài)分布(N(0,1))，即每一個(gè)隨機(jī)粒子都服從數(shù)學(xué)期望為0、方差為1的正態(tài)分布。

(2)均勻分布(U(a,b))：包含兩個(gè)上下限參數(shù)a和b，在幾何上可看作一個(gè)矩形，在概率統(tǒng)計(jì)學(xué)上描述為相同單位長(zhǎng)度上的分布概率相等，其中U(0,1)為標(biāo)準(zhǔn)均勻分布。

(3)F分布(F(m,n))：F分布是由兩個(gè)獨(dú)立的卡方分布組合后的分布，其自由度是m和n，實(shí)際應(yīng)用中多采用F(1,1)進(jìn)行數(shù)學(xué)運(yùn)算。

2 基于PSO和Geiger的微震源定位流程

(1)首先，獲取微震監(jiān)測(cè)到時(shí)數(shù)據(jù)與微震監(jiān)測(cè)系統(tǒng)監(jiān)測(cè)點(diǎn)坐標(biāo)，作為算法的原始輸入數(shù)據(jù)。

(2)根據(jù)線性定位法原理，建立線性矩陣關(guān)系式，隨機(jī)選取4組監(jiān)測(cè)點(diǎn)數(shù)據(jù)，計(jì)算得到初始反演結(jié)果，記為P2。線性法的矩陣關(guān)系式具體為:

P?=β

(3)

(4)

P=(x0,y0,z0,t0)

(5)

(6)

上式中，α為系數(shù)列向量，P表示微震求解參數(shù)，β為結(jié)果矩陣。

(7)

(5)異常值判斷。重復(fù)循環(huán)執(zhí)行上述方法n次(可取20或100)，獲得震源定位近似解集。異常值判定準(zhǔn)則為:

(8)

上式中，異常因子λ(λ>0)表示該法最優(yōu)震源解與Geiger解的各向最大差異閾值，可取為[5，20]；Dist(·)表示兩點(diǎn)之間的空間距離。異常值判定準(zhǔn)則在幾何意義上可描述為一個(gè)以為球心、除去球面的圓球狀幾何圖形，超出此范圍的點(diǎn)被視為異常值，予以去除，異常值數(shù)目記為ρ。

(6)對(duì)n-ρ個(gè)有效震源近似解按照質(zhì)心法求解最終微震源結(jié)果。

3 工程驗(yàn)證

3.1 微震數(shù)據(jù)描述

為了驗(yàn)證本文提出的定位方法的有效性和進(jìn)步性，應(yīng)用呂進(jìn)國(guó)等[9]公開(kāi)的某開(kāi)采礦井工程實(shí)驗(yàn)數(shù)據(jù)進(jìn)行驗(yàn)證。該礦一共安裝了30個(gè)微震監(jiān)測(cè)單元，爆破位置為[8 732.70，6 570.60，511.30]及起爆時(shí)刻為10時(shí)27分，P波傳播速度均值為5 700 m/s，爆破成功后觀測(cè)到8個(gè)P波到時(shí)數(shù)據(jù)，分別記為T(mén)1～T8，如圖1所示，實(shí)驗(yàn)仿真階段采用MATLAB 2019A數(shù)據(jù)處理工具以及相應(yīng)的標(biāo)準(zhǔn)函數(shù)工具箱。

圖1 微震監(jiān)測(cè)點(diǎn)位置

算法中PSO的具體參數(shù)為：粒子數(shù)目N為50，學(xué)習(xí)因子為1.5，學(xué)習(xí)因子為2.5，慣性權(quán)重為0.5，最大迭代次數(shù)D為100。

3.2 粒子群隨機(jī)粒子對(duì)定位的影響

為尋找適合本算法和微震定位模型的隨機(jī)分布，在STD循環(huán)判定準(zhǔn)則、到時(shí)差模型、循環(huán)重復(fù)100次等條件下，測(cè)試正態(tài)分布、F分布、均勻分布等不同隨機(jī)對(duì)象對(duì)算法定位結(jié)果的影響，分布參數(shù)同第一章節(jié)所述。

圖2 隨機(jī)粒子服從不同分布下的定位結(jié)果

表1 隨機(jī)粒子服從不同分布下的定位結(jié)果統(tǒng)計(jì)

從圖2和表1可看出，F(xiàn)分布的定位精度及尋優(yōu)速度皆優(yōu)于其他兩種分布，均計(jì)算耗時(shí)達(dá)到0.972 3 s，定位誤差均值為22 m，而基于均勻分布和正態(tài)分布的微震源定位效果圖都出現(xiàn)了較多的突跳現(xiàn)象，定位性能極不穩(wěn)定，不利于在實(shí)際礦山微震活動(dòng)監(jiān)測(cè)或頁(yè)巖氣勘探中應(yīng)用。

3.3 算法對(duì)比研究

為了進(jìn)一步驗(yàn)證三步反演法的優(yōu)越性，分別取模擬退火算法、粒子群優(yōu)化算法、遺傳算法、非線性最小二乘算法等方法對(duì)微震活動(dòng)監(jiān)測(cè)數(shù)據(jù)進(jìn)行震源反演，不同經(jīng)典算法與本文方法的性能對(duì)比結(jié)果如表2和圖3所示。

由表2和圖3可知，非線性最小二乘法雖然定位性能穩(wěn)定，但是震源求解精度不高，標(biāo)準(zhǔn)粒子群優(yōu)化算法受制于部分超參數(shù)設(shè)置困難，定位結(jié)果誤差極大，體現(xiàn)了PSO易陷入局部最優(yōu)的顯著特點(diǎn)，而本文方法吸取了Geiger法、群體智能優(yōu)化算法等方法的優(yōu)勢(shì)與不足，加入了異常值提出步驟，反演結(jié)果誤差較小，定位性能顯著穩(wěn)健。

圖3 微震源定位方法對(duì)比

表2 微震源定位方法對(duì)比結(jié)果

4 結(jié) 論

(1)本文設(shè)計(jì)了三步微震源反演方法，通過(guò)結(jié)合線性定位法、Geiger法、PSO等算法的優(yōu)劣勢(shì)，將線性定位結(jié)果作為Geiger法的初值，較好地解決了Geiger法受初值影響較大的情況；將Geiger法定位結(jié)果作為下一步的循環(huán)判定準(zhǔn)則，解決了粒子群算法在微震定位求解時(shí)最優(yōu)解隨機(jī)性大、誤差較高、穩(wěn)健性差的問(wèn)題。該算法結(jié)構(gòu)清晰、穩(wěn)定性好，誤差相對(duì)PSO、GA、SA等更小。

(2)對(duì)粒子群算法進(jìn)行改進(jìn)，通過(guò)隨機(jī)初解分布尋優(yōu)，較全面的對(duì)比了正態(tài)分布、均勻分布、F分布等隨機(jī)分布對(duì)算法定位精度與效率的影響，找出合適的隨機(jī)分布匹配PSO算法的隨機(jī)粒子。合理選擇隨機(jī)分布，對(duì)獲得較為準(zhǔn)確、高效的定位算法結(jié)果有重要意義。

鑒于研究時(shí)間有限，并未對(duì)正態(tài)分布、均勻分布、F分布等分布的分布參數(shù)與微震定位的關(guān)聯(lián)做更為詳細(xì)的剖析；可以從各類隨機(jī)分布規(guī)律出發(fā)，獨(dú)立設(shè)計(jì)一種隨機(jī)分布選擇函數(shù)或混合分布函數(shù)，適應(yīng)不同地質(zhì)條件、不同數(shù)據(jù)質(zhì)量下微震定位算法，以達(dá)到算法質(zhì)量的自動(dòng)優(yōu)化。