周辉
（青岛海洋大学海洋地球科学学院，青岛 266003）
何樵登
（长春地质学院地球物理系，长春 130026）
摘要 各向异性介质参数反演通常为非线性优化问题。非线性反演方法可以分为两大类：随机搜索方法，如Monte Carlo法、模拟退火和遗传算法及基于非线性最小平方理论的非线性解析反演方法。遗传算法能寻找到全局最优解，但它为一种较费时的方法。非线性解析反演方法能给出一个与初始模型有关的局部最优解。然而，这种方法具有较快的收敛速度。遗传算法与非线性解析反演方法相结合的反演方法利用这两种反演方法的优点而克服其缺点。因此，结合反演方法既能快速收敛，又能寻找到全局最优解。如何合理地将遗传算法和非线性解析反演方法结合是十分重要的。本文提出一种结合方案，即在连续若干次遗传算法迭代后作一次非线性解析反演。理论算例表明结合反演方法具有上述特点。
关键词 遗传算法 非线性解析反演 非线性结合反演 各向异性介质
1 引言
遗传算法为随机搜索类方法之一，它以概率论为理论基础，用于求解多极值复杂优化问题[9]。遗传算法不要求已知模型空间中后验概率密度的形状并能广泛搜索模型空间。遗传算法模拟自然选择和遗传规律，并遵循适者生存的原则。
遗传算法由Holland在1975年提出[4]。Berg首先将遗传算法应用于地球物理优化问题[1]。Stoffa等系统地研究了种群大小、交叉概率、选择概率和变异概率对多参数优化问题收敛性和收敛速度的影响[11]。Sen等讨论了在选择概率中引入温度参数的作用并提出一些退火方案[10]。周辉等则研究了目标函数与收敛速度和解的精度的关系[16]。
基于最小平方优化理论的非线性反演方法是两大类反演方法之一。当给定的初始模型位于目标函数全局最优解所在的峰谷附近时，这种下降类方法能给出正确解而与初始模型位置无关。下降类算法研究得较深入，应用较广。
Tarantola提出一种基于广义最小二乘法的多维多偏移距声波地震波形解释的一般性非线性地震波形反演方法[12]。随后，Tarantola将该理论推广于各向同性介质的弹性波反演[13]。Gauthier等用理论数据验证了Tarantola提出的方法的正确性[2]。稍后，Tarantola研究非线性解析法反射波弹性反演的策略，指出以纵横波的波阻抗和密度作为反演参数，才尽可能使反演参数之间相互独立[14]。Pan用τ—P变换研究层状声学介质中平面波地震记录非线性解析反演的理论和可行性[6]。为了更多地利用地震数据中的信息，包括VSP资料中反射和转换信息，Mora作了一些工作[5]。当仅用反射数据时反演主要解决引起反射的P波和S波的波阻抗突变。当利用转换数据时，则能分辨大尺度的P波和S波速度变化。Sambridge等改进了修改模型的方法[8]。在子空间中，可同时得到P波、S波波阻抗和密度。周辉等将非线性梯度反演方法推广于多维、多道、多分量任意弹性各向异性介质参数的反演[17]。
非线性解析反演方法和遗传算法结合的反演方法利用非线性解析反演和遗传算法的优点，克服它们的缺点。因此，结合反演方法不仅能搜索到全局最优解，而且能较快地收敛。Porsani等在遗传算法和广义线性反演方法相结合方面作了一些研究[7]。
本文讨论各向异性介质的非线性解析反演方法和遗传算法与非线性解析反演方法相结合的结合反演方法[17]。对于遗传算法读者可参考遗传算法的相关文献[3,9～11]。
2 各向异性介质参数非线性解析反演方法
2.1 共轭梯度法
反演的目的是利用地面或井中测得的位移场ui（xr,t）求取地下介质密度分布ρ（x）和弹性参数分布Cijkl（x）。ρ（x）、Cijkl（x）称为模型参数。x为研究介质中或边界上任一点，x=（x1,x2,x3）,xr为接收点。反演的目标是使目标函数
岩石圈构造和深部作用
取极小值。其中Cd、Cm分别为数据（波场）和模型参数的协方差算子。m0为先验模型参数，m为反演过程中求得的模型参数。由于模型参数有多个，故用向量表示。ucal为给定m的波动方程正演记录，uobs为观测波场，上角标t表示转置。地震记录u和模型参数m之间的函数关系为
岩石圈构造和深部作用
g为非线性算子，（2）式为波动方程的算子形式。记第n次迭代时的模型参数为mn，则有
岩石圈构造和深部作用
及共轭梯度法的迭代公式[15]
岩石圈构造和深部作用
其中Gn为g对mn的Frechet导数，ηn为一常数，可由多种方法计算[5,8]。
梯度
式（4）为梯度反演方法的基本公式。当该公式中的每一量都已知时，迭代就可进行。在这些变量中，最关键的是梯度向量。
2.2 目标函数
在最小二乘理论中，权函数是协方差算子逆的核。假设数据集中的误差是不相关的，它仅取决于时间或源和接收器的位置，那么有[14]
岩石圈构造和深部作用
其中σ为数据的均方差。
2.3 各向异性介质中的弹性波动方程
令fi（x,t;xs）是第s次激发的内体力密度，Ti（x,t;xs）是地球表面S的应力矢量分量，ni（x）是表面的单位法向分量。那么与第s次激发相应的位移由以下微分方程组给出[15]
岩石圈构造和深部作用
2.4 梯度向量
式（4）中梯度向量的分量为[17]
岩石圈构造和深部作用
其中，T为地震记录的长度，
岩石圈构造和深部作用
其中，t∈[T,0]，
3 结合反演方法
3.1 遗传算法和非线性解析反演方法的优缺点
遗传算法是利用概率论来求解多极值复杂优化问题的一种随机搜索方法，由一组随机选取的模型开始，不需要更多的先验信息，广泛而有效地对模型空间的最优部分采样。尽管遗传算法是基于自然选择、遗传规律，搜索模型空间的最优部分而求得最优解，但它是一种计算量很大的方法。由于地震模型空间大，用全局最优化方法估计各向异性介质参数的地震波形反演十分费时。
目标函数的梯度信息是非线性解析反演方法修改模型参数的依据，它能给出一个接近初始模型的一个局部最优解。如果初始模型选择得合适，即当初始模型处在全局最优解所在的目标函数低谷时，非线性解析反演方法能收敛于全局最优解。然而，恰好给出一个接近全局最优解的初始模型的概率是非常小的，尤其对没有模型参数的任何先验信息的情况。但应强调的是，非线性解析反演方法具有较快的收敛速度。
发挥非线性解析反演方法快速收敛和遗传算法能搜索到全局最优解的优点，而克服前者仅能寻找到局部最优解和后者运算量大的缺点是很有意义的。非线性解析反演方法和遗传算法相结合的反演方法可达到上述目的。在结合反演方法中，遗传算法的作用是提供接近全局最优解的模型，非线性解析反演的作用是尽快求出全局最优解。因此，结合反演方法具有搜索到全局最优解的能力和比遗传算法收敛速度快的特点。
3.2 结合方案
遗传算法在优化过程中连续不断地搜索整个模型空间。在每次迭代结束后，得到一个本代的最优模型。根据遗传算法的数学原理[3]，最优模型的数量在下一代中得以增加，同时经交叉和变异作用又有新的模型产生。在下一代种群中，最优模型可能与前一代的相同，也有可能劣于前一代的最优模型。所有这些最优模型可能在目标函数的同一低谷处，也有可能在其它低谷处。遗传算法寻找最优模型要经过多次迭代才能确定一个极值。遗传算法的随机性导致遗传算法是一种费时的方法。然而正是遗传算法的这种随机性保证了它能搜索到全局最优解。
如果将每次遗传算法迭代的最优解作为非线性解析反演的初始模型，非线性解析反演可以找出与初始模型毗邻的局部最优解。由于非线性解析反演是一种确定性的方法，它按目标函数的梯度方向修改模型，所以非线性解析反演方法只需几次迭代即可收敛。非线性解析反演求得的解是否为全局最优解，非线性解析反演方法本身是无法得以保证的。只有当遗传算法提供接近全局最优解的初始模型时，非线性解析方法反演才能收敛到全局最优解。
结合反演方法中遗传算法和非线性解析反演方法的匹配方式是十分重要的。非线性解析反演方法得到接近遗传算法提供的初始模型的局部最优解后，在以后若干代中因遗传算法的随机性而使其最优解与该局部最优解相同。如果每次遗传算法迭代后作非线性解析反演，那么结合反演的结果在几代内都是相同的。显然其中的一些非线性解析反演是没有必要的。因此，结合方式应为在连续多次遗传算法迭代后作一次非线性解析反演，然后将非线性解析反演的结果作为下一代种群中的一个母本模型。图1为结合反演的框图。
图1 结合反演框图
4 算例
为了验证结合反演方法的优越性，对一维多层横向各向同性介质参数的反演理论实例作了分析。
图2是目标函数值与迭代次数的关系图。在该结合反演算例中每次遗传算法迭代后就作一次非线性解析反演迭代。结合反演的误差在开始几次迭代中下降很快，尤其在前3次。结合反演方法在第10次迭代达到的较小误差，遗传算法在第42次迭代才达到。结合反演的误差比遗传算法的跳跃得严重。这是因为非线性解析反演得到的模型在遗传算法中作为母代参加繁衍。这个模型因遗传算法的随机性常常被新的模型替代。这两个模型可能位于目标函数两个不同的低谷中，因此非线性解析反演的结果不同。
尽管结合反演的目标函数有些振荡，但也存在连续几次迭代目标函数几乎不变的现象。这意味着这几次迭代的最优模型是很接近的。在这种情况下非线性解析反演不能提供较大的改进。所以，此时的非线性解析反演是没有必要的，否则只能增加计算量。
图2 结合反演（实线）和遗传算法（虚线）的误差与迭代次数的关系
结合反演中每次遗传算法迭代后作一次非线性解析反演迭代
图3是另一个例子。在该结合反演例子中，每五次遗传算法迭代作一次非线性解析反演。在这里遗传算法占主要地位。此时结合反演的误差函数明显比遗传算法的小。结合反演的误差在第5次迭代末突然下降，并在第10次迭代时的小误差，遗传算法在42代才达到。遗传算法始终没有到达结合反演的最小误差。结合反演的误差在后期迭代过程中平稳下降，这是遗传算法占主导地位的原因。
从该例可知，若遗传算法与非线性解析反演方法比较合理地结合，结合反演方法比遗传算法具有快得多的收敛速度。
5 结论
非线性结合反演方法扬遗传算法和非线性解析反演方法之长，抑其之短，它是一种具有较快收敛速度的全局反演方法。
在结合反演中遗传算法和非线性解析反演方法的结合方式是重要的。从算例可得出，五次遗传算法迭代后作一次非线性解析反演的结合反演的效果明显优于每次遗传算法迭代后都作非线性解析反演的结合反演的效果。但是在结合反演中连续作多少次遗传算法迭代及连续迭代次数在整个迭代过程中的可变性还有待于进一步研究。
图3 结合反演（实线）和遗传算法（虚线）的误差与迭代次数的关系
结合反演中每五次遗传算法迭代后作一次非线性解析反演迭代
在结合反演中遗传算法的作用是提供接近全局最优解的初始模型。结合反演的运算速度主要取决于遗传算法的运算速度。均匀设计理论可以应用于遗传算法以加快随机搜索的速度。
与遗传算法相同，其它随机搜索方法也可用来与非线性解析反演方法形成结合反演方法。
参考文献
[1]E.Berg.Simple convergent genetic algorithm for inversion of multiparameter data.SEG60 Expanded Abstracts,1990，Ⅱ，1126～1128.
[2]O.Gauthier,J.Virieux and A.Tarantola.Two-dimensional nonlinear inversion of seismic waveforms:Numerical results.Geophysics,1986,51,1387～1403.
[3]D.E.Goldberg.Genetic Algorithms in Search,Optimiztion,and Machine Learning.Addison-Wesley,Reading,MA,1989.
[4]J.H.Holland.Adaptation in Natural and Artifical Systems.The University of Michigan Press,Ann Arbor,1975.
[5]P.Mora.2D elastic inversion of multi-offset seismic data.Geophysics,1988,52,2031～2050.
[6]G.S.Pan,R.A.Phinney,and R.I.Odom.Full-waveform inversion of plane-wave seismograms in stratified acoustic media:Theory and feasibility.Geophysics,1988,53,21～31.
[7]M.J.Porsani,P.L.Stoffa,M.K.Sen,et al..A combined Genetic and linear inversion algorithm for seismic wave-form inversion.SEG63 Expanded Abstracts,1993,692～695.
[8]M.S.Sambridge,A.Tatantola and Kennet.An alternative strategy for nonlinear inversion of seismic waveforms.Geophysical Prospecting,1991,39,723～736.
[9]M.Sambridge,and G.Drijkoningen.Genetic algorithms in seismic waveform inversion.Geophys.J.Int.,1992,109,323～342.
[10]M.K.Sen,P.L.Stoffa.Rapid sampling of model space using genetic algorithms:examples from seismic waveform inversion.Geophys.J.Int.,1992,109,323～342.
[11]P.L.Stoffa,M.K.Sen.Nonlinear multiparametre optimization using genetic algorithms:Inversion of plane-wave seismograms.Geophysics,1991,56,1794～1810.
[12]A.Tarantola.Inversion of seismic reflection data in the acoustic approximation.Geophysics,1984（a）,49,1259～1266.
[13]A.Tarantola.The seismic reflection inverse problem.In:F.Santosa,Y.-H.Pao,W.W.System,and C.Holland Eds.Inverse problems of acoustic and elastic waves.Soc.Industr.Appl.Math.,1984（b）,104～181.
[14]A.Tarantola.A strategy for nonlinear elastic inversion of seismic reflection data.Geophysics,1986,51,1893～1903.
[15]A.Tarantola.Inverse problem theory:Methods for data fitting and model parameter estimation.Elsevier Science Publ.Co.Inc.,1987.
[16]周辉，何樵登.遗传算法在各向异性介质参数反演中的应用.长春地质学院学报，1995,25，增刊1,62～67.
[17]周辉.各向异性介质波动方程正演及其非线性反演方法研究.长春地质学院博士论文，1995.