基于GA和PSO智能算法的大地电磁一维反演分析

陈 杰，杨 磊

(1.中水珠江规划勘测设计有限公司，广东 广州 510610；2.中铁第四勘察设计院集团有限公司 地质路基设计研究院，湖北 武汉 430063；

1 引 言

大地电磁测深方法是利用大地中频率范围分布宽广的天然电磁波来探测不同地层深度的电性结构的地球物理勘探方法。大地电磁测深方法具有多项优点：不受高阻屏蔽，对低阻异常分辨率高；不用人工供电，勘探成本低且工作方便；勘探深度范围大。近几十年来，这一系列优势使得大地电磁方法在矿产勘探及普查、地壳岩石圈电性结构研究、海洋地球物理勘探、地热勘探和隧道隐伏岩溶结构勘察等都扮演着至关重要的角色[1-4]。

在实际铁路勘察项目中，复杂山区往往需要建设多条深埋长隧道。高频大地电磁法具有探测深度大、对低阻异常体敏感和对复杂地形适应能力强等优点，因而逐渐被用于铁路隧道勘察设计任务[5-8]。为了更直观地获取地下深部响应特征，大地电磁法数据反演解释工作尤为重要。近年来，国内外学者将许多非线性反演方法应用于地球物理反演问题中，比如粒子群算法、遗传算法、人工神经网络法等，均取得了一定的效果[9]。粒子群优化算法(Particle Swarm Optimization,PSO)是在描述鸟儿寻找食物的Boid模型基础上题设出来的。其核心思想是通过整个群体中每一个个体之间的相互协助以达到整个群体的信息共享，从而在杂乱无序的求解空间中迅速找到最佳求解路线。该算法主要的优势是原理易懂、容易实现，并且不需要调节太多的参数。目前，粒子群优化算法也开始应用于地球物理反演的相关研究。例如，韩瑞通等[10]利用改进的粒子群算法对不同类型的大地电磁一维层状地电模型进行了反演测试；崔益安等[11]利用粒子群优化算法实现了双频激电数据的联合反演；程久龙等[12]利用粒子群优化算法实现了全空间条件下的矿井瞬变电磁法反演；戴前伟等[13]在探地雷达全波形一维反演中采用改进型粒子群算法，发现该反演算法对多参数反演具有良好的抗噪性和有效性。此外，遗传算法(Genetic Algorithm,GA)是模拟大自然中自然选择和遗传学中生物进化理论的一类仿生算法，该算法通过跟寻自然进化和遗传规律来搜索求解域内的最优解；在计算机的模拟运算过程中，其求解过程转变为生物进化中染色体基因的交叉、互换和变异等过程，在处理较为复杂的组合优化问题中具有较好的全局搜索能力，其在信号处理和地球物理反演问题的研究中备受关注。李帝铨等[14]利用遗传算法实现了CSAMT(Controllable Source Audio Frequency Magnetotelluric Method)的最小构造反演；童孝忠[15]利用遗传算法实现了大地电磁测深的正则化反演。李强强[16]针对大地电磁的反演特点,设计了大地电磁的遗传反演程序。李玲[17]采用改进的遗传算法进行一维可控源音频大地电磁法反演研究，其结果反映出该方法具有一定的抗噪能力和实用性。综上所述，PSO算法和GA算法已开始应用于地球物理反演问题研究中，然而目前两类方法的反演效果在不同测试条件下却缺少一定的比较分析。

因此，本文的主要构思是基于大地电磁法的基本理论，通过编写PSO和GA两类非线性反演算法，对大地电磁勘探典型地电模型进行两种算法的反演测试对比分析，为未来两类智能算法在地球物理反演中的应用提供有价值的参考。

2 大地电磁理论和优化算法

2.1 大地电磁法正演理论

大地电磁测深方法是一种利用宇宙中的雷暴、太阳风等垂直入射到大地介质的天然电磁场作为场源信号的电磁勘探方法，其勘探前提也是需要满足卡尼亚经典理论中的先验假设[18]。

在均匀半空间中，当深度为z处时，电场强度可写成为

(1)

此外，电磁波在地下传播过程中，趋肤深度δ基本反映了地层不同深度的电性结构和各频段有效穿透深度之间的关系：

(2)

式中，ρ为电阻率(Ω·m)。趋肤深度反映出了地下岩性不同电阻率ρ情况下的有效穿透深度。一般而言，测点的高频信息反映出浅部的电性特征，低频信息则反映出地层深部电性特征。

视电阻率概念是从均匀介质中电阻率和波阻抗关系引申出来的。假设地电剖面是水平分层均匀的，假设大地构造是由N层水平层状介质组成，其中每一层的电阻率为ρ1，ρ2，ρm…,ρN-1,ρN(单位为Ω·m),厚度为h1,h2,hm,…,hN-1(单位为m)。对于一维地电模型，计算视电阻率ρa和相位角φ为：

式中,Z1为第一层地面的波阻抗,kg/(m2·s)；E和H为空间总场中相交的电场值(V/m)和磁场值(A/m);Re[Z]和Im[Z]为Z的实部和虚部。

地面阻抗递推公式为：

(4)

阻抗Z在每一层地电面具有连续性。式(4)中Zm+1表示的就是m+1层顶面的波阻抗;Z0m表示第m层的特征阻抗;km为第m层复波数;hm为第m层的厚度。对于最下面的第N层，其波阻抗等于介质中的特征阻抗ZN=-iωμ/km。因此，通过最底层特征阻抗的逆向递推求解，进而可以求出整个地层模型的波阻抗信息。

2.2 粒子群优化算法原理及其改进

粒子群算法是模拟鸟儿寻找食物的过程，在该算法中是把鸟群简化成粒子群，把寻找食物的过程简化为寻找最优解。每一个粒子通过跟踪自身路径的最优值和群体路径中的最优值来不断更新自己，第i个粒子的速度更新公式为：

Vi=ωVi+c1ξ(pbestxi-xi)+c2η(gbestxi-xi)

(5)

位置更新公式为：

xi=xi+rVi

(6)

公式中,Vi为第i个粒子的速度;xi为第i个粒子的位置;ω为惯性权重;c1,c2为学习因子，分别代表了自身认知和社会认知;ξ,η为0到1之间的随机数;pbestxi为自身所经过的最优位置;gbestxi为整个粒子群所经过的最优位置;r为约束因子。目前对粒子群优化算法的改进主要集中在对速度构造方式的改进上，本文选择调试后，采用公式(5)进行速度更新。

2.3 遗传算法原理及其改进

2.3.1 模型编码

遗传算法的初始状态是从可能为潜在解集的一个种群开始，其种群则是由许多不同种类和数量基因编码的个体组成，进而组成了遗传物质的载体染色体。其核心是通过将问题模型参数简化为相应的二进制编码，实现了将模型参数转换成基因型的编码工作。对应于大地电磁法反演，参与反演的参数就是每一层的厚度和该层的电阻率值，最后一层层厚度为无穷大，不参与反演。

2.3.2 生成初始模型集群并计算适值

初始模型集群应均匀随机地分布在求解区间内，种群规模应综合考虑计算成本。适值应与目标函数成反比关系，即目标函数值越小，适值越大，模拟 “适者生存”的原理。本文选取的目标函数为：

k=‖log2ρobs-log2ρmod‖2

(7)

式中，k表示为目标函数的拟合误差；ρobs和ρmod分别表示观测的视电阻率值(Ω·m)和模型反演的视电阻率值(Ω·m)。

2.3.3 交叉

交叉是指两个染色体的某一相同位置被切断，染色体个体将重新交叉组合、变换位置等，从而形成两个新的染色体。交叉的方式有许多种，评定交叉方式优劣的标准是能否增加种群的多样性，新的优秀个体是否得到增长。

2.3.4 变异

变异是指染色体个体的某一部位产生了突变，变异为一个新的个体，呈现出一组新的二进制编码，产生新的优秀个体。该过程同时丰富了种群的多样性。

2.3.5 精英选择

精英选择是改进原有的自然选择方式，计算其个体适应度；对性能优秀的个体进行保留，不让其参与到交叉、变异环节，直接保留到新的下一代。精英选择方法可以保证新老优秀个体在一系列变换过程中得到保留，保证了遗传算法在最优化问题的收敛性。

3 模型反演测试

为了探究这两种算法对地球物理非线性反演结果的影响，本文比较测试了简单二层地质模型、复杂三层和四层地质模型，选取了合适的地层参数，其反演测试效果基本相同。在模型测试展示环节，选择具有典型代表性的三层H型地电模型进行细致分析，其中绘出的三层地层模型如图2所示。

图2 大地电磁H型地电模型Fig.2 Magnetotelluric H-type geoelectric model

选取三层H型地层模型为ρ1=50 Ω·m,ρ2=10 Ω·m,ρ3=100 Ω·m，h1=100 m,h2=200 m。分别用PSO，GA两种算法进行理论正演模型数据的反演对比，然后给正演模型数据添加10 %的高斯噪声模拟实际野外数据。根据观测数据估测地电断面各层视电阻率，并将计算得到的趋肤深度作为各层的层厚度参考值，通过设定反演参数的大概范围来约束反演参数区间大小，以此加快迭代反演收敛进程。此次数值模拟计算的计算机配置为：CPU，i5-10210U 四核1.6 GHz,内存8 GB,其工作性能能较好地完成此次计算任务。

从结果图(图3～图6)和反演计算耗时统计(表1)来看，在理论模型分析时，PSO和GA算法均能收敛到较小误差值，几乎完美地拟合了理论模型曲线，也能寻找到最优的模型解。在反演迭代过程中，PSO算法在整个寻优过程是跟随当前最优解的过程，其收敛迭代速度较快，反演结果也较为精确；GA算法是整个种群整体参与寻优过程，其收敛迭代速度较慢，相对于PSO算法而言迭代次数较少。然后对正演模型加入10 %高斯噪声模拟实际大地电磁测深数据，发现PSO算法的收敛迭代速度明显要快于GA算法，将反演结果与理论模型进行比较分析，PSO算法反演结果在一定程度上能较好地拟合正演模型数据。而GA算法的收敛迭代速度较慢，但其是种群整体参与寻优过程，具有优良的全局搜索能力[19,20]，抗干扰寻优效果优于PSO算法。两类反演方法进行耗时统计比较，GA算法采用二进制编码，而PSO算法采用实数编码，PSO反演算法耗时明显少于GA算法，GA反演算法的耗时受到实际计算机内存大小的制约。

图3 PSO优化算法反演结果与理论模型对比Fig.3 Comparison of PSO optimization algorithm inversion results with theoretical model

图4 GA优化算法反演结果与理论模型对比Fig.4 Comparison of GA optimization algorithm inversion results with theoretical model

图5 PSO优化算法反演结果与加噪10 %模型对比Fig.5 Comparison of PSO optimization algorithm inversion results with 10 % noise addition model

图6 GA优化算法反演结果与加噪10 %模型对比Fig.6 Comparison of GA optimization algorithm inversion results with 10 % noise addition model

表1 反演计算耗时统计比较

4 结 论

PSO和GA智能算法都属于完全非线性反演算法，本文利用这两类算法实现了大地电磁一维典型三层地电模型的大地电磁法反演，并进行了抗干扰能力测试，结果表明：

1)PSO算法与GA算法的寻优能力都很强，在不加入误差的时候，两者均能较为准确地找到真实解。加入10 %的高斯误差测试比较，GA算法是整个种群整体参与寻优过程，PSO算法的整个搜索更新过程紧跟当前最优解的粒子，在复杂模型影响下，GA算法的寻优抗干扰能力优于PSO算法。

2)与PSO算法相比，GA算法收敛速度慢，具有优良的全局搜索能力，只要设置参数得当,一般不会陷入局部极小值；PSO算法在整个寻优过程是跟随当前最优解的过程，因而计算收敛速度更快，但更容易陷入局部极小值。

3)GA算法采用二进制编码，中间迭代过程复杂，程序运行时占用的内存较多，导致迭代计算耗时严重。PSO算法采用实数编码，中间迭代过程简单，程序运行时占用的内存明显少于GA算法。

综上所述，两类智能算法各有优劣，今后应该结合实际的观测数据、地质条件、计算能力进行综合的分析、选择。