CC BY
149
УДК 553.9; 51-74
ЧАСТОТНО-ВРЕМЕННОЙ АНАЛИЗ НЕСТАЦИОНАРНЫХ ГЕОФИЗИЧЕСКИХ ПРОЦЕССОВ НА ОСНОВЕ ВЕЙВЛЕТОВ И ЭМПИРИЧЕСКИХ МОД
© А.Е. Филатова, А.Н. Павлов, А.Е. Храмов, А.В. Иванов,
В.В. Грубов, И.А. Яшков, Е.Н. Егоров, М.В. Храмова
Ключевые слова: автоматическая диагностика; вейвлетный анализ; эмпирические моды; цифровой сигнал; сейсморазведка; звуковые и поверхностные волны.
В работе представлены результаты частотно-временного анализа сигналов с локализованными помехами на примере сложных многокомпонентных и зашумленных процессов, регистрируемых при проведении сейсморазведочных работ. Для проведения частотно-временного анализа используются как непрерывный вейвлетный анализ, так и метод эмпирических мод. Показано, что метод разложения сигнала на эмпирические моды является перспективным инструментом исследования структуры экспериментальных данных, который расширяет возможности изучения динамики систем с меняющимися во времени характеристиками.
Существенный прогресс в развитии методов цифровой фильтрации и анализа структуры нестационарных процессов в последние десятилетия был связан с появлением новых инструментов исследования, к числу которых относится концепция вейвлетов [1-5]. Первоначально вейвлет-анализ был предложен в качестве альтернативы классическому спектральному анализу, основанному на преобразовании Фурье [1, 6]. В настоящее время теория вейвлетов является, возможно, единственным новым математическим направлением, которое очень быстро нашло практическое применение почти во всех естественных науках и многих областях техники.
Большую роль играют вейвлеты и при анализе структуры многокомпонентных сигналов с меняющимися во времени характеристиками. Классические вероятностные и спектральные методы анализа [7] являются инструментами исследования стационарных случайных процессов; их применение для обработки нестационарных данных приводит к различным проблемам интерпретации полученных результатов. В частности, обнаружение двух пиков в спектре мощности с некратными частотами может соответствовать принципиально разным случаям: одновременному присутствию двух независимых ритмов колебаний в динамике изучаемой системы или переключению частоты, при котором в каждый момент времени существует только один ритмический процесс [8]. Известно несколько подходов, позволяющих проводить исследование процессов, не ограничиваясь требованием стационарности, включая концепцию аналитического сигнала, анализ флуктуаций относительно тренда, распределение Виг-нера-Вилля и т. д., однако, в соответствии с выводами работ [2, 7-8], вейвлет-анализ обладает наиболее широкими возможностями и целым рядом преимуществ по сравнению с упомянутыми подходами [5].
В последние годы был разработан еще один вариант частотно-временного анализа нестационарных
процессов, получивший название метода эмпирических мод [9]. Данный подход уже хорошо зарекомендовал себя при решении широкого круга задач [10-11]. Более того, в ряде работ приводятся примеры, свидетельствующие о том, что метод эмпирических мод может превосходить вейвлет-анализ с точки зрения частотновременного разрешения [9, 11]. Главное преимущество данного метода по сравнению с вейвлетами заключается в том, что он не требует выбора базиса и настройки параметров преобразования, поэтому его проще использовать на практике при проведении вычислений.
Как вейвлет-анализ, так и метод эмпирических мод представляют собой перспективные инструменты цифровой обработки экспериментальных данных, позволяющие извлекать информацию о меняющихся во времени характеристиках процесса в условиях нестационарной динамики и при наличии помех [5, 12]. Практическое применение этих подходов в разных областях естествознания позволило исследователям существенно расширить возможности адекватной расшифровки информации, содержащейся в экспериментальных данных. Примером могут служить камеральные геофизические работы [13-16]. Наиболее активно развивающимися областями геофизики в настоящее время являются отрасли, связанные с поиском и разведкой залежей различных природных ископаемых и, в первую очередь, нефтяных и газовых месторождений. Максимальный объем ежегодных полевых поисковых и разведочных работ с последующими камеральными исследованиями зарегистрированных данных приходится на сейсмическую разведку методом отраженных волн. При этом обеспечивается большая глубинность исследований земной коры и высокая детальность информации. Так как проведение полевых исследований является дорогостоящим методом получения информации о литологическом строении земной коры, ключевой задачей является повышение качества обработки экспериментальных материалов для извлечения максимально
возможной информации о структуре записанных сигналов. В этой связи привлечение методов анализа, использующих вейвлет-преобразование и разложение сигнала на эмпирические моды, представляется одним из наиболее перспективных вариантов расширения инструментария для адекватной расшифровки такой информации из сложных многокомпонентных и зашумленных сигналов, регистрируемых при проведении сейсморазведочных работ. Наряду с получением дополнительных сведений о литологическом строении земной коры, анализ данных сейсморазведки является также и хорошим тестом на эффективность применяемых методов цифровой обработки сигналов, позволяющим совершенствовать используемый инструментарий и, тем самым, определять дальнейший прогресс в развитии методов анализа нестационарных данных.
Для частотно-временного анализа многокомпонентных сигналов сейсморазведки было использовано непрерывное вейвлетное преобразование (НВП) [1-5], которое, применительно к рассматриваемой задаче, представляет собой свертку исследуемого сигнала, х(ґ), и набора базисных функций ф^ т:
ад /Л
Ж(5, т) = | х(0ф5,т (ґ)Л ; ф^ (ґ) = -1 Фо 1^^), (1)
— ад
где і - временной масштаб, определяющий растяжение или сжатие материнской функции; т - временной сдвиг вейвлетного преобразования; уо('л) - материнский вейвлет. В данной работе был использован материнский Морле-вейвлет:
ц0 (^) = л—1/4е-/'ИоЛе_л2/2 , (2)
поскольку он, как было показано ранее [5, 17], является оптимальным базисом для частотно-временного представления нелинейных нестационарных сигналов различной природы, <в0 - параметр вейвлета Морле (центральная частота).
В силу локализации вейвлета у0 во временной области, коэффициенты преобразования зависят от двух параметров - масштаба і (или частоты / = 1/і) и момента времени ґ. Эти особенности НВП приводят к ряду отличий от классического спектрального анализа, основанного на преобразовании Фурье, которые способны устранить ряд принципиальных недостатков классического подхода [5, 14, 17].
Именно по этой причине НВП воспринимается экспериментаторами как естественное развитие методов спектрального анализа процессов с меняющимися во времени характеристиками, позволяющее получать результаты, интуитивно понятные специалистам в области радиофизики, например, описание многокомпонентных нестационарных процессов в терминах мгновенных амплитуд, мгновенных частот и фаз ритмических составляющих. Для НВП хорошо проработана теория [1-6, 8], позволяющая использовать строгие математические результаты при рассмотрении вопросов спектрального разрешения метода и т. д. Недостатки данного подхода не столь очевидны, и их выявление возможно лишь при проведении различных тестирований, которые позволяют выявить следующие ограни-
чения вейвлет-анализа: 1) существование краевых эффектов [8]; 2) наличие интерференций, если спектральное разрешение недостаточно для разделения спектральных «пиков» в частотной области, причем, этот эффект зависит от степени нестационарности [6]; 3) возникновение эффектов ложной амплитудной (или частотной) модуляции [18], которое тесно связано с наличием интерференций между боковыми спектральными составляющими модулируемого процесса.
Альтернативным НВП методом частотно-временного анализа нестационарных процессов является метод эмпирических мод (ЭМ), который лишен вышеупомянутых недостатков.
Остановимся подробнее на процедуре разложения сигнала по ЭМ. На отрезке сигнала х(Г) между двумя последовательными экстремумами (например, между двумя минимумами на сигнале в моменты времени Г и Г+) можно формально ввести локальную высокочастотную (ВЧ) составляющую {^(Г), Г- < Г < Г+}. Аналогично можно формально ввести низкочастотную (НЧ) составляющую т(Г) следующим образом: т(Г) = х(Г) - <^(Г). Таким образом, сигнал х(Г) на отрезке Г- < Г < Г+ представляется в виде суммы ВЧ <^(Г) и НЧ т(Г) компонентов. Повторяя процедуру п раз для низкочастотной составляющей, можно разложить сигнал на набор мод.
Метод ЭМ предполагает следующую процедуру выделения НЧ и ВЧ составляющих (эмпирических мод) сигнала:
1) нахождение всех экстремумов сигнала х(Г);
2) интерполяция сигнала между минимумами (максимумами) и построение огибающей е^Г) (ещах(Г));
3) расчет средней НЧ составляющей т(Г) = (е^Г) + + вщах(Г))/2;
4) выделение ВЧ составляющей <^(Г) = х(Г) - т(Г);
5) повторение п. 1-4 для НЧ части сигнала т(Г). В п. 1 вместе х(Г) необходимо использовать сигнал т(Г).
Очевидно, что первая эмпирическая мода d^(Г) будет самая высокочастотная, а с ростом номера моды частота будет уменьшаться, поэтому высшие эмпирические моды Зп(Г) (п > 1) окажутся низкочастотными. Среднее значение каждой эмпирической моды нулевое.
В соответствии с выводами статьи [10], метод эмпирических мод позволяет анализировать модулированные по амплитуде (или частоте) сигналы даже в том случае, если частоты модулируемого и модулирующего процессов близки (или почти совпадают). Очевидно, что данный случай представляет практический интерес, т. к. применение НВП будет приводить к сильной интерференции мгновенных частот соответствующих процессов [19]. Как недостаток метода ЭМ можно отметить то обстоятельство, что физическая интерпретация эмпирических мод, на которые раскладывается анализируемый сигнал, часто бывает затруднена из-за наличия большого числа идентифицированных мод с близкими частотами. Само разложение на ЭМ интуитивно менее понятно по сравнению с частотно -временным анализом в рамках НВП. Поэтому мы считаем, что оба рассматриваемых подхода, как НВП, так и метод ЭМ, необходимо рассматривать как потенциально мощные и при этом взаимодополняющие инструменты исследования структуры сигналов с меняющимися во времени характеристиками.
В данной работе наряду с проведением расчетов с использованием собственных программ для выполне-
ния НВП, авторы применяли программу выделения ЭМ, разработанную в научной группе П. Фландрина (http://perso.ens-lyon.fr/patrick.flandrin/emd.html). Перед применением данного метода к экспериментальным данным вначале проводилось его тестирование с использованием стандартных вариантов сигналов с меняющимися характеристиками - переключения частоты колебаний и процессов с монотонным изменением частоты (линейных чирпов) [19]. Результаты, полученные при тщательном подборе параметра <в0 базисной функции Морле и при проведении расчетов методом эмпирических мод, в проводимых нами исследованиях оказались очень похожими и сопоставимыми с точки зрения ошибки вычисления.
Проведя тестирование и убедившись в том, что метод эмпирических мод и НВП позволяют получить сходные результаты с близкой величиной погрешности, дальнейшие расчеты нами проводились для первого из этих подходов, учитывая его преимущества с точки зрения автоматизации процедуры вычисления, т. е. уменьшения роли субъективных факторов, способных повлиять на результаты вычислений. С помощью данного подхода были проанализированы трассы сейсмограммы общего пункта возбуждения, подвергнутые предварительной процедуре фильтрации от волн-помех. Анализ трассы, изображенной на рис. 1, позволил идентифицировать 7 эмпирических мод, частотновременной анализ каждой из которых может быть проведен с помощью классического преобразования Гильберта. Представление получаемых при этом результатов на плоскости время-частота аналогично рассмотрению «хребтов» энергетического спектра непрерывного вейвлет-преобразования (см. подробнее работы [5, 19]). По аналогии с НВП, для более информативного представления результатов представляется целесообразным изображать соответствующие зависимости оттенками цвета, которые будут отражать информацию о мгновенных амплитудах эмпирических мод, а не только об их мгновенных частотах.
Рис. 2 иллюстрирует результаты анализа частотновременной динамики трассы сейсмограммы общего пункта возбуждения (после проведения предварительной фильтрации, подробно описанной в работе [20]) на основе метода эмпирических мод. Данный подход позволяет идентифицировать ряд характерных ритмических процессов в структуре полезного сигнала сейсмо-
хЮ6
2.0 ----------------1----------------1-------■--------
1.0
х 0.0
-1.0
-2.0 ----------------------'---------------------1---------------------
0.0 1.0 2.0 3.0
Г, С
Рис. 1. Вид типичной временной реализации трассы сейсмограммы (трасса 4) общего пункта возбуждения
60
О1-----------1------------1-----------1-----------
О 0.5 1 1.5 2
І, С
Рис. 2. Результаты частотно-временного анализа трассы сейсмограммы после проведенной вейвлет-фильтрации помех (см. работу [20]) на основе метода ЭМ
граммы, идентифицируя особенности литологического строения земной коры в конкретном пункте приема. В табл. 1 показаны характерные частоты ритмических процессов, наблюдаемых в сейсмограмме, выделенные как методом ЭМ (частоты мод), так и с помощью НВП для нескольких трасс анализируемых материалов сейсмограмм. Видно хорошее согласие результатов, полученных с использованием обоих методов, однако, метод ЭМ является менее затратным с точки зрения машинного времени и более простым с точки зрения реализации и модификации алгоритмов.
Следует также отметить, что сопоставление структуры сигналов, регистрируемых при проведении сейсморазведочных работ в разных пунктах, позволяет диагностировать и выявлять отличия в строении земной коры. Учитывая сложность регистрируемых сигналов, привлечение высокоточных инструментов исследования нестационарных зашумленных сигналов, таких как НВП или метод ЭМ, позволяет надеяться на то, что с их помощью удастся существенно повысить качество расшифровки информации из экспериментальных данных, регистрируемых при проведении сейсморазведочных работ.
Итак, в данной работе обсуждается проблема частотно-временного анализа процессов с меняющимися во времени характеристиками на примере данных геофизических сейсморазведочных исследований. Отметим, что именно изучение сложных процессов в геофизике, биологии и медицине является мощным стимулом к развитию и совершенствованию методов цифровой обработки экспериментальных данных. Фактически, большинство новых инструментов анализа временных рядов, активно вовлеченных в научные исследования за последние десятилетия, появилось в результате исследования данных натурных экспериментов в этих областях. Например, первые работы, в которых применялась концепция вейвлетов [1, 6], появились в ходе изучения геофизических процессов. Метод анализа флуктуаций относительно тренда [7], ставший одним из новых методов корреляционного анализа нестационарных процессов, впервые был предложен в качестве инструмента диагностики сердечно-сосудистых патологий. Метод мультифрактального анализа, использующий вейвлет-преобразование [18], рассматриваемый в настоящее время в качестве, возможно, наиболее мощного инструмента статистического анализа неоднородных процессов, приобрел значительную популярность только
Таблица 1
Частоты и моменты времени наиболее интенсивного проявления характерных ритмических процессов, выделенные методом ЭМ и НВП, для различных трасс сейсмограммы
Характерная частота (частота ЭМ/частота метода НВП), Гц
Момент времени, когда данный ритм максимально интенсивен (ЭМ/НВП), с
№ ритма (номер ЭМ) 1 2 3 4 5 6 7
Трасса 4 4.62/4.53 1,31/1,32 5,18/5,18 1,26/1,25 11,23/11,27 1,41/1,43 16,78/16,77 1,02/1,04 23,12/23,18 0,97/0,95 35,2/35,3 0,85/0,82 51,6/51, 5 0,85/0,82
Трасса 8 4,69/4,70 1,38/1,39 5,21/5,22 1,24/1,26 11,24/11,27 1,45/1,46 16,78/16,78 1,04/1,04 23,12/23,16 0,97/0,97 35,3/35,3 0,85/0,82 51,6/51, 5 0,86/0,84
Трасса 16 4,78/4,76 1,42/1,43 5,32/5,30 1,26/1,28 11,26/11,28 1,49/1,49 16,78/16,79 1,06/1,08 23,12/23,17 0,97/0,98 35,4/35,4 0,87/0,86 51,7/51, 5 0,91/0,90
Трасса 32 4,82/4,81 1,52/1,54 5,41/5,39 1,26/1,29 11,27/11,30 1,51/1,52 16,79/16,79 1,12/1,10 23,18/23,18 0,98/0,99 35,5/35,5 0,86/0,86 51,8/51, 6 0,95/0,94
Трасса 48 4,95/4,92 1,60/1,61 5,45/5,46 1,31/1,29 11,30/11,34 1,55/1,56 16,78/16,78 1,14/1,14 23,18/23,17 0,99/1,05 35,5/35,4 0,85/0,84 51,8/51,7 0,98/0,98
после появления работ по анализу сердечного ритма, а до этих работ данный метод почти не применялся в задачах цифровой обработки экспериментальных данных, представляя интерес, в основном, в теоретических исследованиях.
Геофизические данные сейсморазведки представляют собой один из тех примеров, когда задача извлечения максимально возможной информации о структуре земной коры является стимулом к развитию и совершенствованию методов цифровой обработки сигналов. Эти методы могут базироваться на новых математических алгоритмах (включая аппарат НВП или синтетические алгоритмы, которые используют вейвлет-анализ как одну из ступеней анализа численных данных) или алгоритмах, уже апробированных в других областях естественных наук (например, в радиофизике). И здесь одной из наиболее актуальных является задача исследования структуры сигналов сейсморазведки. Так как полевые работы являются наиболее дорогостоящими этапами сейсмогеологических исследований, необходимо решать задачу прогнозирования литологического состава изучаемого участка недр путем максимально детальной и качественной обработки и интерпретации данных экспериментов. Рассмотренный в данной работе метод эмпирических мод является новым инструментом исследования, который ранее не применялся для решения подобного рода задач. Учитывая значительный потенциал данного метода при анализе многокомпонентных процессов в условиях нестационарности и наличия помех, а также принимая во внимание ряд преимуществ данного подхода, например, по сравнению с вейвлет-анализом (применение адаптивных базисов и возможность реализации сравнительно простой процедуры автоматизации расчетов), данный инструмент исследования представляет несомненный интерес как перспективный инструмент получения новой информации о литографической структуре земной коры, позволяющий эффективнее анализировать структуру данных сейсморазведки.
ЛИТЕРАТУРА
1. DaubechiesI. Ten lectures on wavelets. Philadelphia: S.I.A.M., 1992.
2. Meyer Y. Wavelets and operators. Cambridge: Cambridge University Press, 1992.
3. Grossman A., Morlet J. Decomposition of Hardy Functions Into Square Integrable Wavelets of Constant Shape // SIAM Journal Mathematical Analysis. 1984. V. 15. № 4. P. 273.
4. Meyer Y. Wavelets: Algorithms and applications. Philadelphia: S.I.A.M., 1993.
5. Короновский А.А., Храмов А.Е. Непрерывный вейвлетный анализ и его приложения. М.: Физматлит, 2003.
6. Morlet J., Arens G., Fourgeau I., Giard D. Wave propagation and sampling theory. Part 1. Complex signal and scattering in multilayer media // Geophysics. 1982. V. 47. P. 203-221.
7. Бендат Дж., Пирсол А. Прикладной анализ случайных данных. М.: Мир, 1989.
8. Mallat S. G. A wavelet tour of signal processing. N. Y.: Academic Press, 1998.
9. Huang N.E., Shen Z., Long S. The empirical mode decomposition and Hilbert spectrum for nonlinear and non-stationary time series analysis // Proc. R. Soc. London. Ser. A. 1998. V. 454. P. 903-995.
10. Wu Z., Huang N.E. A study of the characteristics of white noise using the empirical mode decomposition method // Proc. R. Soc. London. Ser. A. 2004. V. 460. P. 1597-1611.
11. Huang N.E., Shen Z., Long S.R. A new view of nonlinear water waves: the Hilbert spectrum // Annu. Rev. Fluid Mech. 1999. V. 31. P. 417.
12. Addison P.S. The illustrated wavelet transform handbook: applications in science, engineering, medicine and finance. Philadelphia: IOP Publishing, 2002.
13. Wavelets in geophysics / eds. E. Foufoula-Georgiou, P. Kumar. N. Y.: Academic Press, 1994.
14. Филатова А.Е., Артемьев А.Е., Короновский А.А., Павлов А.Н., Храмов А.Е. Успехи и перспективы применения вейвлетных преобразований для анализа нестационарных нелинейных данных в современной геофизике // Известия вузов. Прикладная нелинейная динамика. 2010. Т. 18. № 3. С. 3-23.
15. Филатова А.Е., Овчинников А.А., Короновский А.А., Храмов А.Е. Применение вейвлетного преобразования для диагностики волн-помех звукового и поверхностного типов по цифровым данным наземной сейсморазведки // Вестник Тамбовского университета. Серия Естественные и технические науки. 2010. Т. 15. № 2. С. 561565.
16. Павлов А.Н., Филатова А.Е., Храмов А.Е. Частотно-временной анализ нестационарных процессов: концепции вейвлетов и эмпирических мод // Известия вузов. Прикладная нелинейная динамика. 2011. T. 19. № 2. С. 141-157.
17. Sitnikova E., Hramov A.E., Koronovskii A.A., Luijtelaar E.L. Sleep spindles and spike-wave discharges in EEG: their generic features, similarities and distinctions disclosed with Fourier transform and continuous wavelet analysis // J. of Neuroscience Methods. 2009. V. 180. P. 304.
18. Павлов А.Н., Анищенко В.С. Мультифрактальный анализ сложных сигналов // Успехи физических наук. 2007. T. 177. № 8. С. 859-876.
19. Pavlov A.N., Makarov V.A., Mosekilde E., Sosnovtseva O. V. // Briefings in Bioinformatics. 2006. V. 7. P. 375.
20. Филатова А.Е., Павлов А.Н., Короновский А.А., Храмов А.Е. Диагностика и фильтрация различных волновых компонент цифровых данных наземной сейсморазведки на основе вейвлетного анализа //
Вестник Тамбовского университета. Серия Естественные и технические науки. 2011. Т. 16. № 2. С. 468-475.
БЛАГОДАРНОСТИ: Работа была поддержана Министерством образования и науки РФ в рамках реализации Федеральной целевой программы «Научные и научно-педагогические кадры инновационной России» на 2009-2013 годы.
Поступила в редакцию 11 июля 2012 г.
Filatova A.E., Pavlov A.N., Khramov A.E., Ivanov A.V., Grubov V.V., Yashkov I.A., Egorov E.N., Khramova M.V. FREQUENCY AND TIME ANALYSIS OF NON-STATIO-
NARY GEO-PHYSICAL PROCESSES BASED ON WAVELETS AND EMPIRICAL MODES
We have presented the results of time-frequency analysis of the non-stationary signals on the example of the complex experimental geophysical multi-component processes recorded during seismic sounding studies. We have used both continuous wavelet analyses and empirical mode methods. We have shown that empirical modes concept is the perspective tools to study the structure of non-stationary multi-component processes and broadens the possibilities of study of dynamics of the complex systems with time-dependent characteristics.
Key words: automatic diagnostics; wavelet analysis; empirical modes; digital signals; prospecting seismology; sonic and surface waves.
