压缩感知(Compressed Sensing, CS)描述的是可从高度不完备的线性测量中高精度重构未知目标。创造性的把L1范数最小化和随机矩阵理论有机结合, 可得到稀疏信号重建的最佳效果。利用压缩感知重构一个连续信号, 不再与香农-奈奎斯特采样定理所要求的“频带”有关, 而是与未知信号的“稀疏度”有关。CS由DONOHO, CANDÈS, TAO三位著名数学家提出[1-2], 并入选了2007年美国十大科技进展, 早期的主要几篇文献的引用已超过5×104次。
从数学上来理解, CS就是求解一个线性方程组Y＝AX。Y是已知的观察数据, X是要重构的目标。特别强调的是, 这里重构信号A是一个非方阵, 其矩阵的行数远少于列数。这相当于这个线性方程组的方程个数远少于未知数的个数。显然, 这样的方程组有无穷多个解, 是一欠定反问题。但是CS告诉我们, 如果先验知道未知目标X有一定的稀疏性或者具有可压缩性, 且A的列之间具有一定的非相关性(比如随机矩阵或勘探中的随机接收器布点), 那么可以用稀疏约束优化的算法X=min‖Y-AX‖2+a‖TX‖0精确求解X。其中, T表示稀疏变换, 零范数表示稀疏约束正则化。
满足这个精确重构所需要的测量数M可以表示为M=Slg(N/S), 其中N表示待重构信号的长度, S表示该信号的稀疏度。即便在有噪声的情况下, 即Y＝AX+Noise, 也可以进行高精度重构。我们补充上面提到的两个概念:稀疏和可压缩性。稀疏是指一个信号如果其大部分振幅值为零, 只有少数元素的振幅为非零, 则称该信号为稀疏信号; 可压缩性是指在某种变换域下的稀疏, 即该信号在时空域下不稀疏, 但可以找到一种稀疏变换, 如傅里叶变换、小波变换、曲波变换等, 使得该信号经过变换后的系数只有少数非零值。我们可以从如下几个角度来理解CS:从数学上讲, CS本质是降维, 从低维空间去研究高维空间; 从信号上讲, CS本质是采样, 从香农采样的频率相关到稀疏度相关; 从工程上讲, CS本质是成本, 从物理测量成本转移到数学计算成本, 用数学计算来弥补实际勘探的采集不足。
从20世纪40年代开始, 统治信号和信息领域的一直是香农-奈奎斯特采样定理, 即要想从离散信号中恢复连续信号, 其离散信号的采样率至少应该是此连续信号最高频率的两倍。可以说大部分相关信号采集的设备都是基于香农定理设计的, 如相机、雷达、核磁共振等医学成像、通讯、信息等众多工程领域的设备。CS理论出现后, 其采样率不再和信号的频带有关, 而是和信号的稀疏度有关, 从而打破了传统理论的束缚。因此很多基于传统理论设计的方法、软件和设备都可以得到升级换代。2007年, 莱斯大学提出的单像素相机就是一个很好的例子, 不再去追求千万像素的分辨率, 而是基于压缩感知理论用一个像素的时间序列成像就可重构出高分辨率图像。
压缩感知是一个理论框架, 不是单一的技术, 其三要素是随机测量、目标的稀疏表示和稀疏促进的优化重构算法。地震勘探中, 降低野外数据采集的成本而又能保证勘探精度是一个很重要的课题。反过来说, 从现有采集的地震数据中能否挖掘出更多的宝贵信息来提高勘探精度？由于CS理论的出现, 这个问题由香农采样转变到稀疏重构的思路上, 带来了很多变革。CS的成功, 在于改变了信号的采集方式。CS所涉及的技术, 在地震勘探的采集、处理、正演模拟、成像反演等方面都会带来变革, 但其核心问题之一还是地震数据和地质目标的稀疏表示。
1 地震勘探中压缩感知的最新研究进展学术期刊Geophysics已多次组织了与压缩感知地震勘探相关的专题, 如2010年关于采样和波场表示的专题、2015关于地震数据表示的专题。The Leading Edge也于2017年8月刊出了压缩感知和高效采集的专刊。压缩感知框架已经渗入了传统地震勘探的各个环节, 例如地震数据的采集、压缩、重构、模拟、偏移、反演等, 并且在很多新的勘探技术中也发挥了作用, 例如同时震源勘探、延时勘探等。不管应用于何种技术, 地震数据的稀疏表示和随机采样始终是压缩感知应用的重要理论基础。本章节首先介绍地震数据稀疏表示和随机采样的最新进展, 然后分别从传统勘探技术和新的勘探技术两个方面介绍压缩感知在地震勘探应用中的最新进展。
1.1 稀疏变换地震数据的稀疏表示按来源可分为3类。第一类是从图像处理技术中借用已有的稀疏变换, 如曲波变换[3-4]、剪切波变换[5-6]等等。曲波变换与剪切波变换均为多尺度多方向小波变换, 适用于表示图像中的线状特征, 因此也适用于表示地震数据的同相轴特征。第二类是利用字典学习方法从地震数据本身自适应构造稀疏变换, 如K-均值奇异值分解(K-SVD)方法[7]、数据驱动紧框架方法(data-driven tight frame, DDT)[8]等等。这类方法通常在块尺度下进行操作, 适用于处理具有复杂结构的地震数据。第三类是直接针对地震数据设计稀疏变换。物理小波变换[9]分别考虑了地震数据在时间上的小波状和空间上的抛物状特征, 能够对叠前地震数据进行稀疏表示。地震小波变换[10-11]考虑了地震同相轴斜率变化的光滑性, 实现了一种针对地震数据的提升小波。非对称C-子波变换[12]考虑了地震子波的非对称性, 能够在对地震数据稀疏表示的同时提取出相关信息, 如走时、衰减参数、非对称性参数等。
除了传统意义上的向量的稀疏性, 矩阵的低秩性(即奇异值的稀疏)也是一种广义上的稀疏性, 并且结合随机采样已经应用于地震数据重构等[13]。此外, 目前还出现了如拓展稀疏以及深度稀疏等新的稀疏性概念。拓展稀疏考虑了能量呈集中分布的情况; 深度稀疏来源于深度学习技术, 即将稀疏基进一步进行稀疏表示, 实现一个递归的稀疏表示过程。
1.2 随机采样采样过程的设计需要满足一定的随机性, 同时满足从高维到低维的采样。随机降维采样在压缩感知地震勘探如何体现呢？首先体现在观测系统的设计上。我们可以采用不满足奈奎斯特采样定律的检波器或者炮点分布来获取采样数据, 以此重构完整数据。SHAHIDI等[14]引入了一种抖动(jitter)随机采样方法, 即检波器随机分布在规则网格所限定的区域内部, 可控制相邻检波器的最大间距, 避免数据大片丢失。NAGHIZADEH等[15]总结了几种采样方法, 包括规则采样、随机采样以及随机加规则采样, 并研究了不同采样的频谱特征。NAGHIZADEH[16]还提出了一种双编织三维地震数据采样方法, 这种采样方法可以配合傅里叶方法进行数据重构。HERRMANN[17]论述了利用随机采样甚至可从更少的数据中获得更多的信息。国内很多学者也做出了突出贡献, 例如陈生昌等[18]提出了几种随机与均匀结合的高效采集方法, 并在数值试验中取得较理想的效果。
除了观测系统的下采样设计, 还可以利用随机同时震源实现下采样。具体来讲, 传统勘探中不同炮是按次序激发, 现在可以同时激发, 减少了采样时间, 但采集中不同炮点的数据会混叠在一起。可以在同时激发中添加随机性来解决这类问题, 比如炮点激发时间有随机延迟(炮编码[19]), 另外可以采用随机的激发波形(波场模拟过程中可以实现[20])。这样在共炮点域具有相关性的数据, 转换到共检波点域就是不相关的, 进而可以进行分离。已有石油公司提出了有效的压缩感知采集系统, 并已在实际勘探中进行了应用测试[21-22]。
1.3 模拟、成像、混合震源、延时勘探有限差分正演模拟的计算量取决于模型尺度以及观测系统的分布, 而不是取决于波场的复杂度。对于大尺度细密度的数据生成, 即使波场非常简单, 仍需要大量的计算资源。HERRMANN等[23]提出了基于压缩感知的同时全波形正演模拟, 通过少量的独立震源同时激发, 得到混叠波场, 然后利用压缩感知恢复算法得到分离后的波场。NEELAMANI等[20]考虑了波动方程格林函数的稀疏性, 在震源同时激发时采用随机带限脉冲波形, 利用压缩感知理论和曲波变换来恢复原始波场。
VERA RODRIGUEZ等[24]实现了基于压缩感知的微地震偏移成像技术。LI等[25]利用压缩感知和随机优化技术实现了地震数据的快速全波形反演。作者使用了两种采样策略:随机震源编码和随机震源位置。LI等[26]和李翔[27]利用稀疏约束的Gauss-Newton方法进行基于下采样数据的全波形反演, 减弱了下采样所产生的人为噪声, 改善了反演结果。王华忠等[28]提出了基于特征波的成像方法, 不是使用全波场而是使用压缩感知框架提取出特征波场, 保持走时特征, 然后进行成像。
为了提高勘探效率, 可以采用同时震源技术得到混叠的波场, 然后利用压缩感知恢复出没有混叠的波场。ABMA等[29]利用傅里叶变换实现了混叠波场的分离, 然后用于偏移成像, 得到了与直接使用无混叠波场相似或更好的结果。原因是使用波场混叠技术可以在相同时间内获取更多炮点位置的波场记录。LIN等[19]利用复杂的震源编码技术以及曲波变换将不同的混合震源数据分离开来。
OGHENEKOHWO等[30]提出了基于压缩感知的低成本延时地震勘探技术。不同时间的炮点布置在同一区域, 但具体位置每次都不同, 不同时间的勘探可以实现信息共享。另外, 检波器点是下采样布置的, 因此勘探成本较低。
2 讨论与结论CS如何在地震勘探开发中发挥作用？首先是针对数据采集, 在CS理论框架下, 可以用比传统方法更少的采集成本实现同等(甚至更高)精度的反演。这个“更少”到底是多少, 它与前面提到的反演目标的稀疏度有关。这个稀疏度(包括变换域的稀疏度), 对一般的地震数据而言, 是有一定的先验知识; 对实际数据的反演精度而言, 只能有一个统计上的结论, 即高概率重构。此外, CS也可以深入挖掘现有勘探采集数据的信息, 从而提高数据重构(道加密)和偏移成像的反演精度。
如何才能用好CS？无非是在实际应用中灵活把握CS的三要素:随机采集(包括炮点和检波器点两方面的随机)、目标的稀疏表达(明确具体问题的目标是高密数据还是地下介质速度或其它物性参数？确定目标后, 再找到或设计一个稀疏变换, 进行目标的稀疏表达)和稀疏约束优化重构的快速算法。一般来说, 重构更高维的目标, 需要用的采集数据(百分比)可更少。因为更高维的目标, 有更多的先验信息可以利用。CS在五维数据的插值重构比在二维或三维数据的重构更能发挥出优势, 比如它有更多的断层信息可利用, 所以如何在高维目标的稀疏表达中更好利用结构的几何信息是一个很关键的课题。
经过近十年的发展, 压缩感知理论已开始深入到地球物理勘探的方方面面, 为降低勘探成本、提高采集效率、提高成像分辨率提供了新的理论支撑。相信在不久的将来, 压缩感知能够真正应用于实际油气勘探, 加速同时震源、页岩油、深层勘探等非常规勘探发展。
压缩感知可以理解为一种浅层的学习, 结合当前的深度学习技术是未来发展的一个方向。
致谢: 感谢李幼铭、张捷、赵邦六、陈生昌、李成博等老师给予的巨大帮助, 感谢课题组于四伟博士在此稿件的撰写过程中做出的很多辅助工作。| [1] | DONOHO D. Compressed sensing[J]. IEEE Transactions on Information Theory, 2006, 52(4): 1289-1306 DOI:10.1109/TIT.2006.871582 |
| [2] | CANDÈS E, ROMBERG J, TAO T. Robust uncertainty principles:exact signal reconstruction from highly incomplete frequency information[J]. IEEE Transactions on Information Theory, 2006, 52(2): 489-509 DOI:10.1109/TIT.2005.862083 |
| [3] | MA J, PLONKA G. The curvelet transform[J]. IEEE Signal Processing Magazine, 2010, 27(2): 118-133 DOI:10.1109/MSP.2009.935453 |
| [4] | CANDÈS E, DONOHO D. New tight frames of curvelets and optimal representations of objects with piecewise C2 singularities[J]. Communications on Pure and Applied Mathematics, 2004, 57(2): 219-266 DOI:10.1002/cpa.v57:2 |
| [5] | EASLEY G, LABATE D, LIM W Q. Sparse directional image representations using the discrete shearlet transform[J]. Applied and Computational Harmonic Analysis, 2008, 25(1): 25-46 DOI:10.1016/j.acha.2007.09.003 |
| [6] |
冯飞, 王征, 刘成明, 等. 基于Shearlet变换稀疏约束地震数据重建[J].
石油物探, 2016, 55(5): 682-691 FENG F, WANG Z, LIU C M, et al. Seismic data reconstruction based on sparse constraint in the shearlet domain[J]. Geophysical Prospecting for Petroleum, 2016, 55(5): 682-691 |
| [7] | BECKOUCHE S, MA J. Simultaneous dictionary learning and denoising for seismic data[J]. Geophysics, 2014, 79(3): A27-A31 DOI:10.1190/geo2013-0382.1 |
| [8] | YU S, MA J, ZHANG X, et al. Interpolation and denoising of high-dimensional seismic data by learning a tight frame[J]. Geophysics, 2015, 80(5): V119-V132 DOI:10.1190/geo2014-0396.1 |
| [9] | ZHANG R, ULRYCH T J. Physical wavelet frame denoising[J]. Geophysics, 2003, 68(1): 225-231 DOI:10.1190/1.1543209 |
| [10] | FOMEL S, LIU Y. Seislet transform and seislet frame[J]. Geophysics, 2010, 75(3): V25-V38 DOI:10.1190/1.3380591 |
| [11] |
刘财, 李鹏, 刘洋, 等. 基于seislet变换的反假频迭代数据插值方法[J].
地球物理学报, 2013, 56(5): 1619-1627 LIU C, LI P, LIU Y, et al. Iterative data interpolation beyond aliasing using seislet transform[J]. Chinese Journal of Geophysics, 2013, 56(5): 1619-1627 DOI:10.6038/cjg20130519 |
| [12] | BOβMANN F, MA J. Asymmetric chirplet transform for sparse representation of seismic data[J]. Geophysics, 2015, 80(6): WD89-WD100 DOI:10.1190/geo2015-0063.1 |
| [13] | YANG Y, MA J, OSHER S. Seismic data reconstruction via matrix completion[J]. Inverse Problems and Imaging, 2013, 7(4): 1379-1392 DOI:10.3934/ipi |
| [14] | SHAHIDI R, TANG G, MA J, et al. Application of randomized sampling schemes to curvelet-based sparsity-promoting seismic data recovery[J]. Geophysical Prospecting, 2013, 61(5): 973-997 DOI:10.1111/gpr.2013.61.issue-5 |
| [15] | NAGHIZADEH M, SACCHI M D. On sampling functions and Fourier reconstruction methods[J]. Geophysics, 2010, 75(6): WB137-WB151 DOI:10.1190/1.3503577 |
| [16] | NAGHIZADEH M. Double-weave 3D seismic acquisition, part 1:sampling and sparse Fourier reconstruction[J]. Geophysics, 2015, 80(6): WB143-WB162 DOI:10.1190/geo2015-0161.1 |
| [17] | HERRMANN F. Randomized sampling and sparsity:getting more information from fewer samples[J]. Geophysics, 2010, 75(6): WB173-WB187 DOI:10.1190/1.3506147 |
| [18] |
陈生昌, 陈国新, 王汉闯. 稀疏性约束的地球物理数据高效采集方法初步研究[J].
石油物探, 2015, 54(1): 24-35 CHEN S C, CHEN G X, WANG H C. The preliminary study on high efficient acquisition of geophysical data with sparsity constraints[J]. Geophysical Prospecting for Petroleum, 2015, 54(1): 24-35 |
| [19] | LIN T T Y, HERRMANN F J. Designing simultaneous acquisitions with compressive sensing[J]. Expanded Abstracts of 71st EAGE Annual Conference, 2009: 4-9 |
| [20] | NEELAMANI R N, KROHN C E, KREBS J R, et al. Efficient seismic forward modeling using simultaneous random sources and sparsity[J]. Geophysics, 2010, 75(6): WB15-WB27 DOI:10.1190/1.3509470 |
| [21] | LI C, MOSHER C, MORLEY L, et al. Joint source deblending and reconstruction for seismic data[J]. Expanded Abstracts of 83rd Annual Internat SEG Mtg, 2013: 82-87 |
| [22] | MOSHER C, LI C, JANISZEWSKI F, et al. Operational deployment of compressive sensing systems for seismic data acquisition[J]. The Leading Edge, 2017, 36(8): 661-669 DOI:10.1190/tle36080661.1 |
| [23] | HERRMANN F J, ERLANGGA Y A, Lin T T. Compressive simultaneous full-waveform simulation[J]. Geophysics, 2009, 74(4): A35-A40 DOI:10.1190/1.3115122 |
| [24] | VERA RODRIGUEZ I, KAZEMI N. Compressive sensing imaging of microseismic events constrained by the sign-bit[J]. Geophysics, 2016, 81(1): KS1-KS10 DOI:10.1190/geo2015-0216.1 |
| [25] | LI X, ARAVKIN A Y, VAN LEEUWEN T, et al. Fast randomized full-waveform inversion with compressive sensing[J]. Geophysics, 2012, 77(3): A13-A17 DOI:10.1190/geo2011-0410.1 |
| [26] | LI X, ESSER E, HERRMANN F J. Modified Gauss-Newton full-waveform inversion explained—Why sparsity-promoting updates do matter[J]. Geophysics, 2016, 81(3): R99-R112 DOI:10.1190/geo2015-0431.1 |
| [27] |
李翔. 基于压缩感知技术的全波形反演[J].
石油物探, 2017, 56(1): 20-25 LI X. Full-waveform inversion from compressively recovered updates[J]. Geophysical Prospecting for Petroleum, 2017, 56(1): 20-25 |
| [28] |
王华忠, 冯波, 王雄文, 等. 特征波反演成像理论框架[J].
石油物探, 2017, 56(1): 38-49 WANG H Z, FENG B, WANG X W, et al. The theoretical framework of characteristic wave inversion imaging[J]. Geophysical Prospecting for Petroleum, 2017, 56(1): 38-49 |
| [29] | ABMA R, HOWE D, FOSTER M, et al. Independent simultaneous source acquisition and processing[J]. Geophysics, 2015, 80(6): WD37-WD44 DOI:10.1190/geo2015-0078.1 |
| [30] | OGHENEKOHWO F, WASON H, ESSER E, et al. Low-cost time-lapse seismic with distributed compressive sensing, part 1:exploiting common information among the vintages[J]. Geophysics, 2017, 82(3): P1-P13 DOI:10.1190/geo2016-0076.1 |