Исследование низкочастотных вариаций во вращении Земли

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

Получены следующие основные результаты, выносимые на защиту: организация и создание базы данных астрометрических наблюдений ПВЗ для исследования низкочастотных вариаций вращения Землиисследование квазишестилетних вариаций в неравномерности вращения Земли на протяжении более 100 летв рамках этого исследования обнаружено изменение характера этих вариаций в конце XX векамодель возбуждения свободной нутации Земли (ЧДП), позволяющая на основе проведенного численного эксперимента объяснить вариации амплитуды и фазы ЧДП с помощью нелинейного механизма передачи углового мохмента от внешних оболочек Земли к внутреннимисследование геофизического возбуждения низкочастотных колебаний среднего полюса Земли (волн Марковича) вариациями барического поля Северо-Атлантического колебанияисследование возможностей сингулярного спектрального разложения для целей прогноза параметров вращения Земли.

Основные результаты данного исследования докладывались на следующих конференциях и семинарах:

1. Международный коллоквиум в рамках комиссии № 5 MAC «Международное сотрудничество в распространении астрономических данных», Россия, СПб, Пулково, 2−9 июля 1996 г.

2. Всероссийская конференция «Внутреннее ядро Земли — 2000», Россия, Москва, ОИФЗ, 27−29 ноября 2000 г.

3. Всероссийская конференция с международным участием «Астрометрия, геодинамика и небесная механика на пороге XXI века», ИПА, СПб, Россия, 19−23 июня 2000 г.

4. Всероссийская астрономическая конференция (ВАК-2001), Россия, СПб, СПбГУ, 6−12 августа 2001 г.

5. Международная конференция Journees 2003, «Astrometry, Geodynamics and Solar System Dynamics: from milliarcseconds to microarcseconds», St. Petersburg, Russia, 22−25 September 2003.

6. Международная конференция «Astronomy in Ukraina: Past, Present and Future», 15−17 июль 2004, Голосеево, Киев, Украина.

7. Всероссийская астрономическая конференция (ВАК-2004), Россия, Москва, МГУ-ГАИШ, 3−10 июня 2004 г.

8. Сагитовские чтения, Москва, МГУ-ГАИШ, 30−31 января 2005 г.

9. Международный симпозиум European Geosciences Union G9 Symposium «Geodetic and Geodynamic programmes of the Central European Initiative», Vienna, Austria, 25−30 April 2005.

Ю.Восьмой съезд AO и международный симпозиум «Астрономия -2005: состояние и перспективы развития», Москва, ГАИШ-МГУ, 1−6 июня 2005 г.

11. Всероссийская конференция «Внутреннее ядро Земли — 2005», Россия, Москва, ИФЗ, 16−17 ноября 2005 г.

12. На заседаниях астрометрических семинаров ГАО РАН.

Результаты работы опубликованы в следующих статьях:

1. Gorshkov V.L., 1995, A statistical analysis of the Pulkovo time service database, A&A Transactions, v.8, pp. 307−309.

2. Горшков B.JI., Губанов B.C., Танюхина B.B., 1996, Приливные деформации Земли по астрономическим наблюдениям за период 1962;1996 гг., В кн.: «Современные проблемы и методы астрометрии и геодинамики», Труды конф., ИПА РАН, СПб., с. 349 — 352.

3. Горшков B. JL, Щербакова Н. В, 1998, О Каталоге служб времени 2, Изв. Г АО, № 213, с. 25−35.

4. Горшков В. Л., Миллер Н. О., Воротков М. В., Баушев А. Н., 1998,.

О персистентности параметров ориентации Земли, Изв. ГАО, № 213, с.269−272.

5. Горшков В., Миллер И. О., Персиянинова Н. Р., Прудникова Е. Я., 2000, Исследование геодинамических рядов методом главных компонент, Изв. ГАО РАН. № 214, с. 173−180.

6. М. В. Воротков, В. Л. Горшков, 2000, Эволюционные портреты геофизических рядов: Приливной процесс и сейсмические ряды, Изв. ГАО РАН, № 214, с.429−439.

7. М. В. Воротков, В. Л. Горшков, Н. О. Миллер, ЕЛ. Прудникова, 2002, Исследование основных составляющих в движении полюса Земли, Изв. ГАО РАН, № 216, с.406−414.

8. Горшков В. Л., Щербакова Н. В., 2002, Изменение долготы Пулкова и дол-гопериодические вариации скорости вращения земли, Изв. ГАО РАН. № 216. с. 430−437.

9. Горшков В. Л., Воротков М. В., 2002, Динамика движения полюса и долго-периодические вариации скорости вращения Земли, Изв. ГАО РАН. № 216. с. 415−425.

10.Gorshkov V., Shcherbakova N., Miller N., Prudnikova E., 2003, Tidal variations from local astrometric EOP sets, in the Proceedings of the Journees «As-trometry, Geodynamics and Solar System dynamics: from milliarcseconds to microarcseconds», A. Finkelshtein and N. Capitane (ed), IAA and Obserava-toire de Paris, p.236−238.

11.Горшков В. Л., Воротков M.B., 2004, Динамика движения полюса и сейсмический процесс, Изв. ГАО РАН. № 217. с. 379−387.

12.Горшков В. Л., 2004, О методах прогнозирования в геодинамике, Изв. ГАО РАН. № 217. с. 365−378.

13.Chapanov, Ya., Vondrak J., Gorshkov V., Ron C., 2005, Six-year cycles of the Earth rotation and gravity, in Proc. EGU G9 Symp. «Geodetic and Geody-namic programmes of the СЕ1», Vienna, Austria, April 2005, pp. 221 — 230.

14.Gorshkov V., Miller N., Naumov V., Prudnikova E., Shcherbakova N., 2005, Pulkovo coordinates from astrooptical observations, Kinematic and Physics Celestial Bodies, Suppl. Ser., № 5, 311−315.

Помимо этих публикаций участие в международном проекте по созданию базы данных астрометрических определений ПВЗ отражено в работах руководителя проекта д-ра Я. Вондрака (Vondrak etal., 1995; 1998).

Вклад автора в работы, выполненные в соавторстве равный.

Введение

В.1. База данных.

В данной работе исследованы взаимодействия некоторых геофизических процессов и их проявления в динамике вектора вращения Земли. Исследования сложных взаимосвязей в геофизической среде невозможны без создания и расширения баз данных (БД) геофизических и геодинамических параметров. Извлечение геофизической информации из наблюдений служб определения параметров вращения Земли (ПВЗ) классическими астрометрическими средствами весьма затруднительно. Частично это обусловлено тем, что в рамках отдельной обсерватории собрать и единообразно переработать весь поток исходных наблюдений невозможно, а в координационные центры поступает итоговая, не всегда единообразно редуцированная информация. Переработка такой неоднородной (в широком смысле) информации и, как правило, сглаживающая модель редукций приводят к тому, что в шумовые остатки переходит многое из того, что в будущем могло бы стать предметом научного исследования.

Создававшаяся в лаборатории вращения Земли (JTB3) Пулковской обсерватории с 80-х годов прошлого века БД исходных наблюдений ПВЗ на классических инструментах оказалась востребованной после запуска КА HIPPARCOS. Рабочей группой MAC «Исследование вращения Земли в системе каталога HIPPARCOS» под руководством д-ра Я. Вондрака (Астрономический Институт АН Чешской Республики) был организован сбор исходных наблюдательных данных ПВЗ классическими средствами определения. При этом предполагалось уточнение долгопериодических вариаций в ПВЗ для следующих геодинамических исследований (Capitane, 1991):

— исследования возбуждения чандлеровской составляющей движения полюса (ЧДП) на основе более точных данных об амплитудных и фазовых вариациях этого колебания,.

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

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

— исследование амплитуды главного члена нутации и свободной нутации ядра для разрешения проблемы его сплюснутости и в целом для построения моделей Земли.

Некоторые из этих проблем исследованы в настоящей работе.

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

В результате этой работы создан однородный ряд ПВЗ по наиболее точным и продолжительным наблюдениям за период с 1899 по 1992 годы, основанный на 3.3 миллионах наблюдений (Vondrak et al., 1995). При этом, помимо решения чисто методических задач по совместному уравниванию разнородных астрометрических наблюдений, была оценена взаимная ориентация систем HIPPARCOS и VLBI (Vondrak et al., 1996). На основании Пулковской части этой международной БД ПВЗ.

— был создан сводный каталог координат (прямых восхождений) и собственных движений КСВ-2 (Горшков, Щербакова, 1998), в котором был реализован новый метод уравнивания (Горшков, 1983), устраняющий регулярную часть суточных ошибок наблюдений,.

— предложена рефракционная модель так называемого эффекта ветра в наблюдениях ПВЗ, определены и учтены в БД вносимые им систематические погрешности (Горшков, 1988),.

— получены новые оценки чисел Лява к и, А (Горшков и др., 1996; Gorshkov et al., 2003),.

— уточнены координаты Пулкова и оценена скорость их изменения в течение XX века (Горшков, Щербакова, 2002; Gorshkov et al., 2005).

В связи с вышеупомянутой задачей исследования взаимосвязей в геофизической среде БД в значительной степени расширена интегрированными в нее геофизическими данными, доступными через интернет или опубликованными в открытой печати. Ассимилированные данные представлены в виде удобном для развиваемого пакета прикладных программ по анализу БД и частично доступны на сайте (http://www.gao.spb.ru/russian/1 g/database.html) лаборатории геодинамики ГАО РАН.

В.2. Обзор исследований низкочастотных вариаций во вращении Земли.

Возросшие точности и плотности рядов ПВЗ активизировали исследования высокочастотных вариаций вектора вращения Земли, в частности, исследования свободной нутации ядра. Однако в области его низкочастотных вариаций остались вопросы, дискуссии по которым продолжаются до сих пор. В основном они касаются проблем возбуждения и амплитудно-частотной модуляции чандлеровского движения полюса (ЧДП) и декадных вариаций продолжительности суток (ПС). Как заметил Р. Гросс (Gross, 2000), «хотя ЧДП исследуется более столетия, механизм его возбуждения остается ускользающим (elusive)».

Предлагается широкий спектр возможных механизмов возбуждения ЧДП вплоть до рассмотрения различных механизмов, объясняющих это колебание как вынужденное (Авсюк, 1996; Курбасова, Рыхлова, 2001, 2002; Акуленко и др., 2002). Если в 1980;х годах влияние атмосферы и океана на возбуждение ЧДП оценивалось в 25−30% (Wahr, 1983), то из последних работ следует, что совместного влияния все более точно учитываемых геофизических факторов (атмосферный, океанический и гидрологический угловые моменты) достаточно для возбуждения ЧДП. При этом учитываются:

— возбуждение океаническими течениями Эль-Ниньо (Сидоренков, 1997) и Северо-Атлантическим колебанием (Chao, Zhou, 1999),.

— атмосферное возбуждение с учетом вклада ветра (Brzezinski, Petrov, 1995; Furuya et al., 1996),.

— возбуждение океаническими течениями и донными нагрузками (Ponte et at, 1998; Ponte, Stammer, 1999),.

— совместное влияние океанических и атмосферных процессов (Celaya et al., 1999; Gross, 2000; Brzezinski, Nastula, 2002; Liao et al., 2003).

Рассматривается также возможность возбуждения ЧДП гироскопическим моментом сил со стороны внутреннего ядра (Dumberry, Bloxham, 2002; Спиридонов, Акименко, 2003).

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

Для описания возмущенного вращения Земли и, в частности, для решения проблем возбуждения ЧДП используются дифференциальные линеаризованные уравнения Лиувилля, выражающие сохранение момента импульса во вращающейся вместе с Землей системе координат. Применяемый в данной работе обобщенный вид уравнений относится к Земле, испытывающей упругие деформации под переменными атмосферно-океаническими нагрузками и вращательные деформации. Детальное описание вывода этих уравнений можно найти в монографиях (Манкк, Макдональд, 1964), (Lambeck, 1980), (Мориц, Мюллер, 1992). Для координат мгновенного полюса вращения, не сильно удаляющегося от полюса фигуры, это уравнение в компактной комплексной форме записывается следующим образом (Gross, 2000):

РхУ + (i/Ocv,)dpxy/dt = (1 в).

Здесь р — мгновенные координаты полюса в комплексной плоскости р — х + iy Осу/ ~ &>cw (l+i/2Qc") — комплексная частота свободной нутации, со^ = 0.84 цик.

1 Система координат (начало — в центре масс Земли, координата л- — в направлении гринвичского меридиана, у — в направлении 90° в. д) вращается вместе с Землёй. При задании у в направлении 90° з.д. координаты полюса запишутся следующим образом р~х — iy. В данных IERS положение мгновенного полюса вращения (в соглашении IERS (2003) Celestial Intermediate Pole — CIP, что соответствует учёту нутаций до 2 суток) относительно полюса фигуры (принятого в ITRS референц-полюса Земли — IRP) определяется именно таким образом.

Система координат для возбуждающих функций стандартна и, следовательно, Хху=Хх + iXy. лов в год {сру в дальнейшем) — наблюдаемая частота чандлеровского движения полюса, Qcw — добротность системы на резонансной частоте со^ (для упрощения в дальнейшем обозначим Q = Qcw). Данные прошлых астрометрических наблюдений приводили к оценкам Q меньше 30, что обусловлено шириной полосы спектральной линии ЧДП, определяемой точностью наблюдений. В настоящее время значения фактора добротности оцениваются в пределах от 0=121 (Brzez-inski, Petrov, 1995) до (2=179 (Wilson, Vicente, 1990). В работе (Спиридонов, Акименко, 2003), основанной на прямом вычислении координат полюса по данным моментов импульса атмосферы и океана (ААМ и ОАМ) при варьировании со^ и в одной из оптимальных моделей получено Q ~ 30.

Для описания возбуждения неравномерности вращения Земли (осевая компонента вектора вращения Земли) уравнение Лиувилля записывается следующим образом:

ALOD/LODo = -Хт (2в), где LOD0 = 86 400 сек, a ALOD — разность между наблюдаемой продолжительностью суток и LODo. Для упрощения в дальнейшем будем применять обозначение LOD = ALOD.

Правые части уравнений (1в) и (2в) содержат компоненты возбуждающей функции момента импульса Хху— Хх + iXy И Хт, Эти функции были введены в работе (Barnes et al., 1983) и названы возбуждающими функциями углового момента. Функции Хху удобны для вычисления экваториальных эффективных угловых моментов, а ХТ — для вычисления осевого углового момента по атмосферным (ААМftp://ftp.aer.corn/pub/collaborations/sba/') и океаническим (ОАМftp://euler.jpl.nasa.gov/sbo//) данным. Каждый из компонентов Хху и Хт можно представить в виде суммы функций, определяемых переменной частью давления атмосферных и океанических масс (ХР), и функций, определяемых движе.

JV нием этих масс (X). Имея геофизические данные в виде рядов давления, скоростей атмосферных ветров и океанических течений, можно оценить их вклад в динамику вектора вращения Земли от перераспределения масс (У) и относиW тельных движений (X), иначе говоря, за счёт изменения тензора инерции ЛI (t) и относительного момента количества движения h (t). В случае исследования сезонных и более низкочастотных составляющих вектора вращения Земли членами с производными от Al (t) и H (t) можно пренебречь ввиду их малости (1:365 для сезонных вариаций) (Сидоренков, 2002) и тогда Х^ и Хт для жесткой не-деформируемой Земли запишутся следующим образом:

Zxy=(&I + Qh)/(C-A).

Z г = (AI33 + Qh3)/C (Зв).

Здесь С = 8.0365 и, А — 8.0101 — осевой и экваториальный моменты инерции Земли (С-А = 0.26 398) в единицах 1037 кг-м2 (Groten, Е., 1999), Q = 7.292 115* 10~5 радиан/сек — средняя угловая скорость вращения Земли.

Для работы с конкретными атмосферными и океаническими данными Al (t) и H (t) удобнее выразить в сферических координатах. Тогда компоненты функ.

Р VV ций возбуждения Хил с учётом упругих деформаций Земли запишутся следующим образом (Barnes et al., 1983): р -1.098/г4 9х. 2 .

X х = 77—тт— J J Ps sm 9 cos (poosXdXd.

C-A)g J J p= 1.098 Л 4 я.

УС у.

— I .иу"я f г. 2. ,, , ,.

—-— ря sin q> cos (p sin Л a A, d (p.

C — A) g I J p 0.7 5 57?4 ^f^j 2 .

Хт =—-p, — I л cos cp аласр

CS I I (4b).

X&trade——1−5913Я f r f^MSjn^C0S^C0S д vcos^sinХ)дЛд (рАр

1.5913Л3.

Xу =——г—J j J (u sin cp cos qy sin Я + vcos (pco s Xdcpdp.

R 3×7 = 777T-J J «COS2?» dXdcpdp.

Здесь для географических координат использованы стандартные обозначения ф и R = 6371 км — средний радиус Земли, g = 9.81 м/сек2 — среднее ускорение силы тяжести, ps (q>, X, t) — атмосферное давление на поверхности Земли, u ((p, X, p, t) и v ((p, X, p, t) — восточная и северная составляющие скорости ветра на высоте, где атмосферное давление составляет р (10 < р < ps мбар в рядах ААМ). В работе Н. С. Сидоренкова (2002) приведены более точные формулы, учитывающие сжатие Земли и изменение ускорения силы тяжести g ((p, h), где h — высота слоя над поверхностью Земли.

Численные коэффициенты в (4в) (трансформирующая функция в терминах работы (Манк, Макдональд, 1964)) отражают реакцию упругой Земли на нагрузку и вращательную деформацию. Возникающие при этих деформациях дополнительные потенциалы вызывают добавочные возбуждения через изменения тензора инерции и зависят от соответствующих чисел Лява. Помимо этого, жидкое ядро Земли, имея очень узкую околосуточную резонансную частоту, практически не участвует в вариациях вектора вращения Земли для периодов от нескольких дней до нескольких лет. Поэтому в этих коэффициентах также учтена редукция за счет перехода от моментов инерции всей Земли (С и А) к моментам инерции мантии и коры (Ст и Ат) (Eubanks, 1993). Заметим, что эти коэффициенты необходимо изменить, если возбуждающий процесс не связан с атмосферными или океаническими нагрузками (например, землетрясения). Кроме того, в разных источниках эти численные коэффициенты слегка различаются в зависимости от выбранных значений чисел Лява.

Вычислив по различным временным рядам геофизических данных возбуждающие функции с помощью уравнений (4в), можно, подставив их в (1 В, 2в), после численного интегрирования получить оценки «геофизических» координат полюса р (х, у), которые затем можно сравнить с астрономическими наблюдениями параметров вращения Земли (ПВЗ). Именно по такой схеме проведено моделирование и оценка низкочастотных составляющих в движении полюса в данном исследовании.

Большинство исследований возбуждения ЧДП геофизическими процессами основаны на перевычислении (reanalysis) метеорологических данных с 1948 года Национальным Центром предсказаний окружающей среды (NCEP/NCAR) США (http://www.ncep.noaa.gov/reanalysis/). Более точные данные, основанные на метеорологии в эпоху ИСЗ, приведены с 1979 года в reanalysis2. Данные по океаническим угловым моментам имеются с 1980 года (Ponte et al., 1998; Gross et al., 2003). Однако самый интригующий с точки зрения тестирования любых моделей момент в зарегистрированной истории движения полюса приходится на конец двадцатых годов прошлого столетия. Тогда произошло почти полное затухание ЧДП с последующей потерей фазы этого колебания без каких-либо заметных возмущений годичного компонента, а, стало быть, и без серьезных метеорологических аномалий. Предпринимаются некоторые «реставрационные» метеорологические усилия в этом направлении по данным о поверхностной температуре (Rosen, Salstein, 2000), однако надо отчетливо понимать, что полученные ряды будут носить в значительной степени модельный характер, так как данные по всему профилю атмосферы в то время не регистрировались. Более того, даже имеющиеся метеорологические данные в немалой степени носят модельный характер.

Что касается годичного компонента в движении полюса, то его возбуждение совместным действием атмосферного и океанического угловых моментов (ААМ и ООМ) не вызывает сомнений. В последнее время появились работы, посвящённые влиянию сезонного компонента водного баланса атмосферы и океана на возбуждение вариации вектора вращения Земли (Jochmann, 1999; Wunsch, 2000). Годичный компонент движения полюса, будучи не строго периодичным и не повторяющимся по амплитуде год от года, может в свою очередь служить возбуждающим фактором для свободной нутации Земли. На это обратил внимание ещё Джеффрис (I960), указав на «возможность возбуждения свободной нутации нелинейным откликом на него внутренних слоев Земли».

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

Помимо резонансного ЧДП в движении полюса обнаруживаются более низкочастотные, нерегулярные и слабые колебания (волны Марковича). В данной работе проведено совместное исследование колебаний полюса, СевероАтлантического колебания (North-Atlantic Oscillation — NAO) и уровня Балтийского моря. Показана возбуждающая роль мощного атмосферно-океанического явления NAO в образовании волн Марковича, а также синхронизм колебаний уровня Балтийского моря и Ладожского озера с индексом NAO.

В области неравномерности вращения Земли самые мощные процессы, имеющие квазирегулярный характер с характерным временем порядка десятков лет (так называемые декадные вариации), также не имеют однозначного геофизического объяснения. Дискутируются, в основном электромагнитное, обусловленное топографией границы ядро-мантия, гравитационное (Kuang, Bloxham, 1997; Rubincam, 2003) и связанное с ним конвективное (Pais, Hulot, 2000; Гох-берг, 1995) взаимодействия мантии и ядра для их объяснения. Однако помимо декадных и довольно мощных приливных вариаций во вращении Земли присутствуют квазирегулярные вариации с характерным временем порядка 5−8 лет, которые практически не исследованы. В данной работе рассматривается возможность возбуждения этих вариаций взаимодействием внутренних оболочек Земли, нелинейно реагирующих на фазовые соотношения сезонных и чандле-ровских колебаний полюса.

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

В.З. Обзор применяемых методов исследования.

В.3.1. Сингулярный спектральный анализ.

Помимо традиционных методов статистических исследований (различные фильтрации, спектральные методы), одним из основных методов исследования в данной работе был метод сингулярного спектрального анализа (SSA) в программной реализации СПбГУ (Данилов, Жиглявский, 1997). Так как его использование пока не столь широко, как активно применяемые различные варианты Фурье анализа и вэйвлет анализа, приведем его краткое описание.

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

Задавшись числом М < N/2, значениями исходного ряда длины N последовательно заполняют строки матрицы X. При этом первая строка содержит первые М элементов ряда, вторая — со второго элемента по М+1 и так далее, пока ряд не исчерпается. После центрировки по столбцам и соответствующей нормировки вычисляется корреляционная матрица R = ХХТ, сингулярное разложение которой R = PLPT даёт диагональную матрицу собственных чисел L и ортогональную матрицу собственных векторов Р матрицы R. В программной реализации этого метода (http://vvvv.gistatgroup.coin/gus/). названного авторами «Гусеница» (М — длина гусеницы), главные компоненты исходной матрицы Y-XP могут быть исследованы, визуализированы и упорядочены по возрастанию их вклада в исходный ряд. Это позволяет интерактивно производить непосредственный поиск гармонических компонентов, фильтрацию или сглаживание ряда, выбирая соответствующие значимые компоненты Y. Ввиду ортогональности матрицы г можно восстановить матрицу X — УР используя при атом, выбранные главные компоненты Y.

L-ооственные вектора г корреляционной матрицы выступают и роли передаточных функций соответствующих фильтров. Ширина полосы пропускания зависит от формы передаточной функции фильтра и определяется как видом собственного вектора, так и длиной интервала усреднения М. Чем больше М, тем уже полоса фильтра. Выбор нескольких главных компонентов Y эквивалентен параллельному соединению соответствующих фильтров, что позволяет управлять формой спектральной характеристики. При выборе М, которое значительно меньше характерной ритмичности исследуемого ряда (в пределе при А/=2), фактически происходит его сглаживание. Периодические, но не обязательно гармонические составляющие исследуемого ряда образуют графически хорошо различимую пару соседних компонентов Y (типа фигур Лиссажу).

Данный метод имеет определённые аналогии с вэйвлет анализом и динамическим Фурье анализом (если ряд состоит из набора строго гармонических компонентов, то, фактически, осуществляется разложение в ряд Фурье). Наиболее важные преимущества как самого метода, так и его реализации для данного исследования состоят в том, что:

— базовые функции метода порождаются самим исследуемым рядом, так как являются собственными векторами корреляционной матрицы R. (В Фурье анализе разложение происходит на гармонические компоненты, а в вэйвлет анализе — на выбранные локально симметричные базовые функции, наиболее подходящие к исследуемому процессу.);

— имеется возможность восстановления ряда по информативным, необязательно гармоническим компонентам, исследование которых интерактивно доступно;

— возможна оценка не только мгновенных частоты и амплитуды периодических компонентов анализируемого ряда, но и их фазы;

— по выбранным информативным компонентам возможен прогноз ряда.

Исследование временных рядов с помощью SSA, тем не менее, не является «автоматической» процедурой. Имеется ряд обстоятельств, которые обязательно должны быть приняты во внимание.

Помимо стандартных ограничений (равномерная сетка, стабильность дисперсии), необходимо соизмерять продолжительности всего ряда N и характерный период исследуемой составляющей Т. В случае достаточной мощности этой составляющей соизмеримость не имеет значения, но чем ближе уровень исследуемого процесса к уровню шумов, тем всё более важным становится це-лочисленность отношения N/T. Из этого следует, что длину «гусеницы» М следует выбирать по возможности близкой к Т, при этом приходится жертвовать (укорачивать) длиной исходной реализации. Практика работы с рядами разной природы с помощью SSA показывает, что для каждого ряда существует свое оптимальное отношение соизмеримости М, N и исследуемого процесса Т. Теоретические особенности и возможности метода приведены в ряде статей авторов на том же сайте (http://www.gistatgroup.com/gus/books.htrril).

В реальных физических процессах зачастую происходит мультипликативное взаимодействие составляющих. Большинство методов позволяют выделить аддитивные составляющие, однако в случае их не строгой регулярности появляются фиктивные компоненты. Все компоненты, выделяемые методом «Гусеница», аддитивны, но, тем не менее, использование метода в мультипликативной ситуации оказывается продуктивным.

Нами в работе (Воротков и др., 2002) проведено исследование метода на модельных рядах с нестабильной амплитудой и фазой. На рис. 1 В из этой работы представлен ряд длиной 1560 точек с проявлением нестабильности в двух местах (160, 860 точки) продолжительностью 110 точек. Модельный ряд состоит из гармонической составляющей с периодом 11 точек и шумовой составляющей (<т в пределах 20 — 80% от амплитуды гармоники) с нормальным распределением. Гармоническая составляющая была подвергнута амплитудной модуляции таким образом, чтобы амплитуда на указанных участках нестабильности плавно уменьшалась вдвое и затем восстанавливалась. На этих же участках таким же образом изменялась на половину периода фаза. Использованная модель создавалась «похожей» на ряд ЧДП. В одной из работ на вышеупомянутом сайте разработчиков исследовалась модель с экспоненциальной модуляцией амплитуды гармоники с последующей автоматизацией процедуры идентификации и группировки восстанавливаемых компонент. Там же отмечено, что визуальный способ идентификации остается самым гибким. В данном исследовании всегда использовался именно такой способ выбора восстанавливаемых компонент.

Рис. 1 В. SSА-декомпозиция и восстановление зашумленного гармонического ряда с амплитудно-фазовой модуляцией. Вверху — модельный ряд, а) — главный компонент восстановленной гармонической составляющей, б) — компонент, отражающая амплитудно-фазовую модуляцию, в) — сумма (а+б).

На рис. 1 В второй график (а) отображает компонент, сопоставимый с исходной гармоникой модельного ряда. Третий график (б), отражающий амплитудно-фазовую «интервенцию», в сумме со вторым даёт гармоническую составляющую исходного ряда (в). Эти компоненты не являются тождественными исходным неаддитивным составляющим. Тем не менее, можно выявить ха актер взаимодействия этих составляющих в различные моменты времени. График (а) демонстрирует наличие интервалов нестационарности ряда, а размах колебаний графика (б) может быть интерпретирован, как интенсивность модулирующего процесса. Видно, что метод позволяет анализировать ряды, порожденные неаддитивным взаимодействием процессов. При этом оказалось, что результат практически не зависит от уровня шума.

Метод продолжает работать при одновременном исследовании нескольких дискретных временных рядов одинаковой длины (многомерный SSA). В этом случае в многомерную выборку преобразуется многомерный же временной ряд. Главные компоненты, полученные этим методом, являются общими для всей системы рядов, в то время как собственные вектора состоят из частей, соответствующих отдельным рядам. Опыт показывает, что при совместном исследовании нескольких временных рядов сначала необходимо убедиться с помощью одномерного варианта SSA, что в каждом из этих рядов присутствует искомая компонента. Как отмечают разработчики продукта теорию одномерного метода SSA можно считать в основном законченной, в то время как достаточно полной теории для многомерного SSA (MSSA — Multi-channel SSA) не существует. В частности для MSSA не разработаны методы прогноза и нахождения моментов разладки многомерных рядов.

В.3.2. Анализ персистентности (R/S — статистика) гой мето п о аммно нами еализованный и использ емый в анно" работе, связан с исследованием динамики поведения показателя Харста (Я), полученного по скользящему фрагменту исследуемых временных рядов.

В 50-х годах Харстом (Hurst, 1957) был предложен метод нормированного размаха для исследования временных рядов. Суть метода состоит в исследовании накопленных за определенный интервал времени отклонений процесса x (t) от среднего значения. Пусть для дискретного ряда Хпсреднее процесса за этот интервал, тогда накопленное отклонение процесса в момент tj равно: j i, n ~ (xi ~ Xn). i=i.

Размах этой величины для п значений ряда равен:

R (n) = max X, «- min X, «.

1< j<).

Харст исследовал нормированный размах г = R/S, где S — дисперсия ряда. Поэтому метод иногда так и называют R/S — статистикой.

За прошедшее с тех пор время этим методом были исследованы десятки эмпирических рядов, характеризующих самый широкий спектр атмосферных, гидрологических, геологических и астрономических (числа Вольфа) рядов. Процессы такого рода относятся к обобщенным броуновским движениям (фрактальным, дробным — в отечественной литературе нет устоявшегося аналога английскому варианту fractional Brownian motion). Они имеют определенную «память», статистически выражающуюся в независимой от времени корреляции последующих приращений процесса от предшествующих.

C (0 = 22/M-l .

Параметр Н (0 < Н < 1) определяет степень персистентности. Такие процессы обладают свойством статистического самоподобия (масштабной инвариантностью) (Кроновер, 2000). Это значит, что преобразование масштаба времени в р раз (например, равномерное прореживание), а масштаба приращений в рн раз, не меняет вид распределения процесса, то есть.

X (t+pAt)-X (t) = pH (X (t+&)-X (t)). Отсюда следует, что приращение процесса пропорционально AtH.

Из того, что приращение процесса пропорционально AtH, следует также пропорциональность r (At) ~ At н, на основании чего после логарифмирования и производится оценка параметра Я. Этот параметр является устойчивым (ро-бастным) по отношению к варьированию типа распределения приращений процесса X (t). На этом основании по оценке Я, полученной для некоторого временного ряда, делают вывод о наличии персистентности (Я > 0.5), её отсутствии (Я" 0.5) или о наличии антиперсистентности (Я< 0.5) в данном временном ряде. Значению Я = 0.5 соответствует процесс с независимыми приращениями (обычное броуновское движение). Если же Я^ 0.5, то приращения процесса стохастически зависимы, причем значениям Я > 0.5 отвечает положительная корреляция приращений, а значениям Я < 0.5 — отрицательная.

Персистентность свидетельствует о наличии в системе процессов, поддерживающих (в целом) наметившиеся тенденции к изменению её состояния, а антиперсистентность, наоборот, о процессах, препятствующих изменению состояния системы.

Для многих природных процессов этот показатель в среднем равен Я= 0.73, то есть является некоторой универсальной постоянной. При исследовании длительных эмпирических рядов значения Я > 0.5 свидетельствуют о наличии в них долговременных самоподдерживающих тенденций (персистентности), то есть текущее состояние процесса в значительной степени зависит от его предыдущих состояний. Можно говорить об эффектах памяти в таких рядах или об их фрактальности. Мандельбротом (1982) действительно была показана линейная зависимость между фрактальной размерностью ряда D и показателем Харста {D-2H).

При использовании метода Харста необходимо помнить о некоторых условиях, которым должны удовлетворять исследуемые ряды:

— необходима большая их продолжительность, т.к. существует отмеченная Мандельбротом возможность асимптотической сходимости Я—> 0.5 при п —> оо;

— необходимо отсутствие значимых периодичностей и трендов в них, что приводит к завышенным значениям Н.

Последнее обстоятельство не столь существенно, т.к. в работе (Peters, 1994) продемонстрирована возможность приближенно оценивать эти неявные периодичности в ряде по виду кумулятивной кривой и, более того, определять периоды, в которых фрактальная размерность ряда может оказаться различной. Однако для выделения любых не стохастических (детерминированных) составляющих существуют более эффективные методы, например, сингулярный спектральный анализ, описанный выше. Поэтому для надежного суждения об эффектах долговременной памяти во временных рядах совершенно необходимо предварительно отделить детерминированные компоненты от стохастических.

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