Многоструктурная модель геофизического сигнала

УДК 519.6

МНОГОСТРУКТУРНАЯ МОДЕЛЬ ГЕОФИЗИЧЕСКОГО СИГНАЛА

© 2007 г О.В. Мандрикова

A new method of processing of signals having complicated structure is offered in the article. The described in the article technique suggests combination of autoregression method - integrated moving average and fast wavelet transformation. The way of construction of the obtained multi-structural model of signal is considered on the example of processing of registration data of the critical frequency foF2.

Предлагается новый метод обработки геофизических сигналов, имеющих сложную структуру. Сложная структура сигналов регистрации различных геофизических параметров затрудняет процесс построения моделей для таких временных рядов с целью исследования динамических процессов на всех уровнях атмосферы и в различных областях магнитосферы Земли. Для поиска устойчивых изменений применяют сглаживание, что автоматически приводит к потере значимых коротко-периодных вариаций. Для оценки периодических изменений компонент временного ряда получили большое признание ряды Фурье, но наличие в геофизических сигналах резких пикообразных особенностей влечет большие погрешности в расчетах, а также невозможность идентификации точного местоположения особенности. Научный интерес к проблемам исследования взаимодействия Солнца, атмосферы и литосферы определяется тем, что результаты таких исследований крайне важны для решения целого ряда фундаментальных научных задач физики атмосферы, ионосферы, магнитосферы, распространения радиоволн. Одновременно с регулярными изменениями параметров ионосферы, соответствующих их суточному и сезонному ходу, в сейсмоактивных регионах Земли может наблюдаться также и аномальное поведение распределения характеристик ионосферы, наблюдаемых накануне землетрясения и определяемых механизмами лито-сферно-ионосферных взаимодействий. Однако регистрация надежного прогностического признака сейсмического события в вариациях ионосферных параметров сталкивается с серьезными трудностями, которые заключаются в том, что они проявляется на фоне возмущений, определяющихся активностью Солнца. В статье предлагается метод обработки сигналов регистрации ионосферных параметров, в основе которого лежит конструкция вейвлет-разложения. Эта математическая конструкция позволяет исследовать сигнал c достаточной степенью локализации и является одним из наиболее эффективных методов детального изучения скрытых закономерностей нестационарных временных рядов.

Метод исследования изменений в параметрах ионосферы рассматривается в статье на примере обработки сигналов критической частоты /сР2. Вариации критической частоты имеют сложную структуру. На фоне регулярных изменений, обусловленных суточным и сезонным ходом, могут возникать резкие одиночные «пики» длительностью от несколько десятков минут до нескольких часов. В свою очередь, их вели-

чина и длительность зависят от многих факторов. В процессе решения поставленной задачи был разработан новый подход к построению алгоритмов обработки геофизических сигналов, который предполагает совмещение модели авторегрессии - проинтегрированного скользящего среднего (АРПСС), традиционно используемой для изучения закономерностей нестационарных временных рядов, и быстрого вейвлет-преобразования (конструкция кратномасштабного анализа (КМА)). Данная задача становится реализуемой благодаря возможности вейвлет-разложения произвольного сигнала на составляющие. Эта конструкция позволяет работать с дискретными данными, дает информацию о структуре сигнала, позволяет отделить высокочастотные локальные нестационарности и провести дополнительный анализ компонент, а также получить информацию о масштабах процесса. Образующаяся после применения конструкции КМА сглаженная составляющая ряда имеет более простую структуру и может быть аппроксимирована моделью АРПСС. Описание временных рядов на основе параметрических методов позволяет получить точные оценки статистических характеристик процесса и имеет высокое спектральное разрешение. Получаемая при этом модель процесса достаточно легко интерпретируется физически, что играет важную роль в задачах анализа экспериментально полученных данных. Но модели данного класса имеют ряд ограничений, и это сужает область их применения (в первую очередь делается допущение о том, что временной ряд может быть описан линейным стохастическим дифференциальным уравнением, вторым важным моментом является стационарность временного ряда). В том случае, когда известна физическая модель процесса, модификация параметрической модели, как правило, не является сложной задачей и может быть решена с применением широкого класса известных нелинейных моделей. Но в большинстве задач обнаружения зависимость между различными переменными неизвестна. Предлагаемый подход к построению модели позволяет не только справиться с описанными выше проблемами, но и значительно улучшает качество обработки подобных сигналов.

Многоструктурная модель сигнала

1. Разделение сигнала на составляющие в пространстве Ь2(Я).

Пусть /&12{К) - данные регистрации ионосферного параметра; /0 - некоторая проекция / на

¥0 =с1оя г Ш2 [1, 2]; У0 в этом случае

Ь

является пространством выборки; /0 - данные регистрации / в V). Нулевой уровень детальности соответствует интервалу взятия отсчетов. Уровень детальности в вейвлет-преобразовании принято ассоциировать с частотным дискретом в классическом дискретном Фурье-преобразовании. Ограничиваясь пятым уровнем разложения, мы получаем

К, = W_x

-W_ з

-w_5

Тогда /0 имеет единственное разложение:

/о (0 = g-i(0 + g-2(О + g-з (0 + «-4(0 + g-5 (0 + f-5 (0 ,

fi(x) = ^cJkul(2Jx-k)eVi

= 4

keZ,

gJ(x) = 'ZdJky/(2J х-к)(

к

J J k>

W;

= ^

где ц- - скейлинг-функция; ц/ - соответствующий ей вейвлет.

Тогда /0(х) единственным образом определяется последовательностями с] и dJ . Последовательности коэффициентов 3J являются детализирующими.

Коэффициенты 5-го уровня разложения с] соответствуют сглаженной составляющей сигнала. Если в разложении используется ортогональный вейвлет, то прямая сумма становится ортогональной и разложение / ( х) не только единственно, но и компоненты разложения взаимно ортогональны.

Определенная свобода выбора базисных функций позволяет подобрать вейвлеты, наиболее подходящие для обработки каждого конкретного сигнала.

2. Идентификация модели АРПСС для сглаженной компоненты сигнала.

Для сглаженной компоненты сигнала, полученной после применения конструкции кратномасштабного разложения и представленной в виде коэффициентов

с], идентифицируется параметрическая модель из класса моделей АРПСС. Процесс идентификации состоит из следующих этапов [3, 4]:

1. определение порядка разности 3 модели, обеспечивающее стационарность процесса;

2. идентификация типа модели из класса моделей АРПСС;

3. оценка параметров модели;

4. диагностическая проверка полученной модели.

Последовательность коэффициентов с\ после выполнения данных этапов аппроксимируется моделью: фФУёск-> =в{В)ак\ ф(В) = 1-ф1В-....-фрВр;

-oqBq-,

R"cJ - cJ

п ск ~ к-п

няющую кратномасштабное представление компонент временного ряда и модель АРПСС. В итоге мы получаем нелинейную параметрическую модель нестационарного временного ряда, не выдвигая при этом каких-либо предположений по поводу вида ее нелинейности и учитывая нестационарности лишь определенного типа. В этом заключается основное преимущество данного метода.

Модель имеет вид

z =

j,k<=z

£ (zt,ЧЬо;

ot=Vdzt, 4Kjk{t) = 2JI24>{2Jt-kb0)

-e(t), - базисный

вейвлет, - соответствующий ему двойственный; фу,...,фр- коэффициенты авторегрессии.

Моделирование сигнала критической частоты

При обработке использовались файлы с часовыми данными критической частоты /ср2. Всего для экспериментов было взято 19 файлов, содержащих результаты измерений за период 1968 - 1986 гг. Обработка производилась в среде МайаЬ, представляющей собой универсальную интегрированную систему компьютерной математики. Учитывая сезонный характер поведения временного ряда, процедура его разложения с последующей идентификацией модели АРПСС проводилась отдельно для каждого времени года. Так как предоставленные данные имеют разрывы, то с целью уменьшения погрешности получаемых результатов в качестве критерия выбора сезона выберем сезоны с наименьшим количеством пропусков. Анализ данных показал, что таким сезоном является зима. Также для устранения краевого эффекта в начале и в конце сигнала временные ряды дополним значениями, захватывающими конец осени и начало весны (рис. 1, зима 1971 г., 2664 отсчета).

6{В) = 1 -вхВ-

Таким образом, процесс построения модели включает два основных этапа. Первый этап заключается в применении конструкции кратномасштабного анализа. Второй этап базируется на идентификации и оценке модели АРПСС. Получаемая при этом модель сигнала содержит многоуровневую структуру, объеди-

Рис. 1. Результаты измерения данных критической частоты 1^2, 15 ноября - 15 марта. 1971 г., 2664 отсчета

На рис. 2 показаны результаты разложения: детализирующие компоненты с 1 по 4-й уровни и аппроксимирующая компонента 5-го уровня. Анализируя полученные результаты, видно, что детализирующие компоненты 3-го, 4-го уровней разложения содержат ярко выраженный суточный ход.

k

Рис. 2. Результаты дискретного вейвлет-разложения сигнала

критической частоты ^2, 15 ноября-15 марта. 1971г.

На основе методики построения многоструктурной модели далее необходимо идентифицировать модель АРПСС для полученных после преобразования сглаженных компонент сигнала с]. В соответствии с разработанной методикой не указывается, на каком именно уровне разложения необходимо выделить сглаженную компоненту сигнала, подлежащую идентификации АР-модели. Эта задача должна быть решена в процессе проведения экспериментов и полностью зависит от особенностей конкретного временного ряда и целей моделирования. Важным здесь является выбор информативных компонент исследуемого процесса. В каждом отдельном случае информативными компонентами могут быть как детализирующие компоненты определенных уровней разложения сигнала, так и сглаженные компоненты. В процессе вейвлет-разложения, с одной стороны, мы получаем сглаженную компоненту все более простой структуры, что соответственно упрощает как процесс идентификации АР-модели, так и понижает ее порядок, с другой - часть полезной информации, содержащейся в сигнале, переходит в компоненты разложения. Поэтому выбор уровня разложения является важной задачей, и от него во многом зависят результаты исследования.

Поскольку основной целью построения модели сигнала критической частоты является выявление возможных аномальных эффектов, связанных с повышенной сейсмической активностью, в первую очередь АР-модель была идентифицируема для периодов лет, не содержащих сейсмических событий рассматриваемого класса. Идентификация модели была проведена для сглаженных компонент 3, 4 и 5-го уровней разложения. В процессе было выявлено, что наилучшим уровнем

разложения сигнала критической частоты является 4-й. На 3-м уровне сглаженная компонента имеет еще достаточно сложную структуру, что затрудняет процесс идентификации модели и показывает плохие оценки полученных моделей при их диагностике (существенная автокорреляция остатков, первые значения автокорреляций остатков частично выходят за доверительные границы). На 5-м уровне разложения для каждого отдельного года была получена модель либо второго, либо третьего, либо четвертого порядка. Модели имели хорошие характеристики, но общий анализ ионосферного процесса в данном случае провести затруднительно. В отличие от описанного компоненты 4-го уровня разложения показали при идентификации удивительно хорошие результаты. Для всех сейсмически спокойных лет были идентифицированы АР-модели пятого порядка, имеющие хорошие характеристики (остатки являются белым шумом; допустимые значения первых коэффициентов автокорреляции). Более того, параметры полученных моделей достаточно близки между собой (наблюдаемая разница не превышает по значению 0,2). В табл. 1 приведены результаты моделирования сигналов нескольких лет. Особо следует отметить момент яркого проявления в результатах моделирования солнечного одиннадцатилетнего цикла (1972, 1984 гг.). Модели этих лет имеют наибольшее сходство друг с другом: константы, характеризующие средний уровень процесса, и параметры модели. На рис. 3-5 показаны соответственно модели компонент сигналов за 1971, 1984, 1985 гг.

Таблица 1

Параметры моделей сигналов

Модель компоненты сигнала

1971 г. 1972 г. 1984 г. 1985 г.

Константа = Константа= Константа= Константа=

=0,03183 =0,04083 =0,04617 =0,01554

Коэффициенты

АР(1) -0,6477 -0,7416 -0,7442 -0,8846

АР(2) -0,7928 -0,5505 -0,6598 -0,7797

АР(3) -0,2275 -0,09206 -0,2324 -0,2702

АР(4) -0,539 -0,3662 -0,4047 -0,3374

АР(5) -0,3564 -0,5355 -0,4494 -0,3738

Анализ результатов моделирования

Для проведения анализа полученных результатов моделирования с целью выявления возможных аномальных процессов и явлений в сигнале критической частоты накануне наиболее сильных сейсмических явлений были использованы данные каталога землетрясений п-ова Камчатка. Поскольку идентификация проводилась для зимних периодов указанных лет, не содержащих сейсмических событий, землетрясения анализируемого класса наблюдаются только в весенние месяцы, они приведены в табл. 2. На рисунках, содержащих модели компонент, моменты произошедших землетрясений отмечены стрелками.

Анализ этих компонент показывает, что во всех случаях накануне сейсмических событий наблюдается резкое увеличение амплитуды колебаний. Выявленные аномальные области показаны на рисунках пунктирными линиями.

Таблица 2

Данные сейсмического каталога п-ова Камчатка

Год Дата Энергия землетрясений

1971 14 марта 12,9

1984 5 марта 13

1985 3 апреля 13,2

Т-Г  '1-1-1-1-1-г——г

Рис. 3. Модель компоненты 4-го уровня вейвлет-разложения сигнала критической частоты, 15 ноября-15 марта, 1971 г.

Рис. 4. Модель компоненты 4-го уровня вейвлет-разложения сигнала критической частоты, 15 ноября-15 марта, 1984 г.

Таким образом, на основе результатов обработки геофизических данных были наглядно показаны преимущества разработанной многоструктурной модели данных, которая на основе совмещения известной модели временных рядов - АРПСС и сравнительно новой математической конструкции - вейвлет-преобразования позволяет моделировать сложные

геофизические процессы. Модельное представление сигнала критической частоты показало, что временные периоды без повышенной сейсмической активности, с предварительным их разделением на сезонные и разложением в компоненты четвертого масштабного уровня на основе конструкции дискретного вейвлет-преобразования, имеют линейную структуру и поддаются моделированию.

Рис. 5. Модель компоненты 4-го уровня вейвлет-разложения сигнала критической частоты, 15 ноября-15 марта, 1985 г.

Также анализ полученных результатов обработки временных рядов различных районов и их сопоставление с данными сейсмического каталога позволил выявить разномасштабные аномальные эффекты в ионосферном процессе, которые предшествуют наиболее сильным сейсмическим событиям (энергетического класса с k>12,5 в радиусе R~200 км от П.-Камчатского).

Литература

1. Daubechies I. Ten Lectures on Wavelets: Пер. с англ. Ижевск, 2001.

2. Chui K.C. An Introduction to Wavelets: Пер. с англ. М., 2001.

3. Бокс Дж., Дженкинс Г. Анализ временных рядов прогноз и управление: Пер. с англ. М., 1974.

4. Никифоров И.В. Последовательное обнаружение изменения свойств временных рядов. М., 1983.

Институт космофизических исследований и распространения радиоволн ДВО РАН_27 февршш 2007 г.