Химическая физика, 2019, T. 38, № 7, стр. 44-48

1. ВВЕДЕНИЕ

Основу многочисленных исследований наблюдавшихся перед сильными сейсмическими событиями возмущений параметров ионосферной плазмы в качестве возможных ионосферных предвестников землетрясений (ИПЗ) составил обширный массив данных наземных и спутниковых средств наблюдений. Основными источниками данных на начальных этапах были станции вертикального, наклонного и возвратно-наклонного зондирования ионосферы. В целях более детального изучения был осуществлен анализ данных ряда проектов, таких как OGO-6, AEC, ISIS-2, в том числе специальных миссий DEMETER, QuakeSat, CSES-1. Тем не менее, более чем за 40 лет исследований не было выработано консенсуса в научном сообществе по вопросу существования ИПЗ. Развитие применения спутниковых навигационных систем (GPS – Global Positioning System, ГЛОНАСС – Глобальная Навигационная Спутниковая Система) для картирования полного электронного содержания (ПЭС) ионосферы Земли, а также публикация NASA (National Aeronautics and Space Administration) глобальных карт ПЭС ионосферы в свободном доступе (ftp://cddis.nasa.gov.gnss/products/ionex/) сделало этот параметр основным (как минимум наиболее частым) объектом для исследований как обеспечивающий глобальное покрытие по пространству (разрешение 2° по широте и 5° по долготе) и практически непрерывные наблюдения по времени (шаг – 2 ч). Тем не менее, большинство посвященных данной тематике работ ограничивается исследованием частных случаев – так называемых “case studies”. Такой подход позволяет выявить некоторые особенности ПЭС ионосферы в периоды, предшествующие конкретным сильным сейсмических событиям, но требует отдельного обоснования возможности распространения установленных особенностей для общего случая. Количество публикаций, анализирующих значительное количество событий, достаточно мало. Это должно быть связано с большой трудоемкостью ручного анализа глобальных карт ПЭС ионосферы. При этом анализ пространственных характеристик предполагаемых ПЭС-предвестников в большинстве работ не автоматизирован.

2. МЕТОДЫ

Обзор существующих публикаций по проблеме ИПЗ не позволяет выявить единый или общепринятый метод, применяемый для расчета фоновой (невозмущенной) вариации в целях выявления возмущений ПЭС ионосферы литосферного происхождения. Тем не менее, можно выделить ряд наиболее часто используемых подходов. Приведем их.

1. Выбор ближайшего магнитоспокойного дня (или осредненных значений таких ближайших дней) [1] в качестве опорного.

2. Использование эмпирических моделей, например, различных версий справочной модели ионосферы IRI (International Reference Ionosphere), для расчета фоновой невозмущенной вариации [2].

3. Использование различных фильтров, включая “скользящие” средние [3], вейвлет-анализ, калмановскую фильтрацию, интерквартильный анализ [4], а также метод анализа главных компонент [5], искусственные нейронные сети [6], генетические алгоритмы [7] и некоторые другие. Наибольшее распространение исторически в качестве фоновой вариации получили различные средние, вероятно, в силу простоты их реализации и невысоких требований к вычислительным ресурсам. При их расчете разными коллективами авторов использовалось различное количество дней для определения спокойных значений. Обычно за ширину окна “скользящего” среднего принимают 7, 14, 15, 27, 30 или 31 день [8–11]. Обоснования выбора именно такого количества дней для расчета невозмущенной вариации в статьях обычно не приводится, зависимость получаемой вариации от количества дней, использованных при ее определении, не исследуется и не обсуждается. При этом обычно рассчитывается невзвешенное среднее, т.е. вклад каждого дня учитывается с одинаковым весовым коэффициентом (вне зависимости от возмущенности и его близости к рассматриваемому событию). Сам интервал (окно сглаживания) может быть размещен для отсчетов, предшествующих, следующих или же центрально относительно рассчитываемого момента.

3. ПРОБЛЕМЫ И ЗАДАЧИ

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

1. Отсутствие обоснования применяемой методики предобработки исходных данных. Например, в работе [12] используют линейную интерполяцию глобальных карт ПЭС для увеличения временного разрешения с 2 ч в исходных данных до 1 ч в анализируемой выборке. Обоснование того, какие преимущества это дало, не приводится. Отсутствие искажений в статистических свойствах полученной выборки по сравнению с исходной не проверено. Таким образом, массив отсчетов увеличен в два раза, что может привести к увеличению расхода машинного времени, а в худшем – привнести связанные с процедурой или накоплением численных ошибок особенности, не присутствующие в исходных данных.

2. Использование неявных предположений о свойствах анализируемых данных. На статистические свойства выборки накладываются некоторые ограничения (наиболее частое из которых – гауссов закон распределения величины), которые явным образом не приводятся, выполнение этих ограничений не проверяется и в лучшем случае заменяется некоторыми общими рассуждениями (о центральной предельной теореме). Таким образом, применимость метода и полученные с его помощью результаты оказываются под вопросом.

3. Фоновая вариация рассчитывается без обоснования выбранного метода вычисления опорной вариации. В расчетах разных авторов (или разных работах одного коллектива авторов) используется разное количество дней для определения спокойного невозмущенного состояния, но при этом не исследуется вопрос зависимости получаемых результатов от данного параметра. Пример такой зависимости приведен на рис. 1.

Рис. 1.

Карты отклонений ПЭС ионосферы от невозмущенного состояния для 9 января 2010 г., 5 LT. Величина возмущения рассчитана в процентах (верхний ряд) и с нормировкой абсолютных единиц на выборочное среднеквадратическое отклонение, т.е. в количестве СКО (нижний ряд). Фоновая вариация определялась как “бегущее” среднее по предшествующим моментам времени. Правый столбец соответствует ширине окна “скользящего” среднего, равной 7 дням, левый – 15 дням. Ось абсцисс – геомагнитная долгота, ось ординат – геомагнитная широта. Черная кривая – линия положения терминатора, звезда – положение эпицентра землетрясения, ромб – соответствующая магнитосопряженная точка.

4. Несогласованное определение критериев аномальности возмущения. Эта проблема напрямую связана с упомянутой выше проблемой определения фоновой вариации. Авторы обычно не исследуют зависимость морфологических свойств выявленных ПЭС-предвестников от, например, количества дней, учтенных в расчете фоновой вариации. Эта проблема обозначена в работах [13] и [14] и требует дальнейших исследований. При этом сами отклонения определяются как в абсолютных (TECU – Total Electron Content Units, 1 TECU = 10¹⁶ электрон/м²), так и относительных (проценты, количество среднеквадратичных отклонений) единицах. При этом не исследуется вопрос зависимости выявленных особенностей от используемого подхода (например, характерные времена проявления ИПЗ для “абсолютных” и “относительных” вариантов различаются). Для обнаружения в рассчитанных возмущениях аномальных значений требуется задать количественный критерий “аномальности” данных. Обычно в качестве последнего задается некоторое фиксированное значение, например, 50%, или определяется по аналогии с доверительным интервалом. В последнем случае аномальным считаются значения, выходящие за границы интервала M(X) ± k ⋅ СКО, где M(X) –выборочное среднее (или математическое ожидание) совокупности, используемой в вычислениях; СКО – выборочное среднеквадратическое отклонение; k – определяемый экспертным путем множитель, обычно принимающий значения 1, 1.5, 2, 2.5 или 3. При этом в известных нам работах не исследуется тот факт, что при одном и том же значении k, но разных размерах окна “бегущего” среднего, используемого для расчета фоновой вариации (например, 15 или 30 дней), будут получены разные ПЭС-предвестники – разной амплитуды и линейных размеров. Из этого следует, что значение k должно зависеть от количества дней, учитываемых при расчете невозмущенной опорной вариации. Вопрос калибровки используемых методов определения порога “аномальности” данных многими авторами не рассматривается. Более того, использование исследователями несогласованных критериев для расчета и детектирования аномалий делает эти результаты несопоставимыми, по крайней мере для прямого количественного сравнения, что также создает сложности для выявления общих закономерностей на базе анализа конкретных отдельных сильных сейсмических событий.

5. Многие исследования ограничиваются рассмотрением только одномерных данных. Даже те работы, в которых анализируются данные глобальных карт ПЭС, чаще всего рассматривают их как совокупность одномерных временных вариаций (независимый временной ряд в каждом узле сетки). Такое рассмотрение не позволяет исследовать существенные особенности ПЭС-аномалий, связанные с их линейными размерами и положением в пространстве. Те же работы, которые тем или иным образом учитывают совместные вариации в нескольких точках (см., например, [15, 16]), тоже имеют ряд ограничений. В работе [15] требуют (см. стр. 698), чтобы приемники GPS-сигналов располагались на одной и той же (или очень близкой) геомагнитной широте, а также и долгота не различалась слишком сильно. Давиденко [16] фактически выполняет анализ данных в узле (см. рис. 8 в работе [16], нижнюю панель, где по оси абсцис – количество дней до события, по оси ординат – LT – местное время) или паре узлов, т.е. осуществляет кросс-корреляционный анализ между парой соседних приемников.

6. Условия гелиогеофизической активности не учитываются при анализе. В большинстве работ “учет” геомагнитных возмущений осуществляется путем исключения возмущенных дней из рассмотрения (в том числе из расчета фоновой вариации). Из немногих работ, в которых используется иной подход к учету возмущений гелиогеофизической активности, можно отметить работы [16, 17]. В работе [17] с помощью вейвлет-преобразований (в русскоязычной литературе вместо термина “вейвлет” иногда используется термин “всплеск”) исключил из вариаций ПЭС ионосферы, предшествующих М9.0 землетрясению, составляющую, определяемую изменениями солнечной активности. Тем не менее в работе [17] не обосновывается выбранный базис вейвлет-разложения, его применимость для разных гелиогеомагнитных условий, а также не проведено исследование зависимости получаемых результатов от выбранного базиса, например от широко распространенных – Хоара, Гаусса, Добеши, “Мексиканская шляпа” и др. В работе [17] приведен пример использования этого метода только для одного события; не сообщалось об успешном применении данного метода для обработки большой выборки событий. Аналогичные замечания справедливы и для работы [16].

7. Ни одна из известных нам работ не использует для анализа поставляемые НАСА вместе с глобальными картами ПЭС ионосферы соответствующие им карты выборочных среднеквадратичных отклонений (RMS). Исследователи, которые восстанавливают ПЭС самостоятельно (а не использовали данные NASA), не приводят соответствующих оценок СКО для своих данных. Таким образом, не удается оценить, как полученные значения соотносятся с величиной СКО (т.е. с оценкой естественной изменчивости исходных данных для того же периода).

4. ОБСУЖДЕНИЕ

Систематическому количественному анализу накопленных и представленных в научных публикациях материалов по наблюдавшимся перед сильными землетрясениями возмущениям ПЭС ионосферы препятствует ряд факторов. К наиболее существенным из них следует отнести недостаточным образом обоснованные преобразования исходных данных и отсутствие учета зависимости получаемых результатов от применяемого метода определения опорной и возмущенной вариаций. В случае работ, анализирующих поставляемые NASA глобальные карты ПЭС ионосферы Земли (в ionex-файлах ПЭС-данные маркированы как TEC), отдельно стоит отметить отсутствие учета исследователями поставляемых вместе с этими ПЭС-картами соответствующих карт выборочного СКО (в ionex-файлах маркированы как RMS, т.е. “root-mean-square”). Вместе эти факторы приводят к невозможности количественного сопоставления полученных амплитуд “аномальных” вариаций ПЭС с уровнем естественной изменчивости ионосферы.

Определение аномалий из соотношений X ‒ kL < X < X + kL (Х – некоторая оценка фоновой вариации, например математическое ожидание, калиброванная эмпирическая кривая, прогноз и т.п.) выглядит разумным, но требует дополнительных исследований. В приведенном соотношении L является некоторой заданной характеристикой изменчивости вариации, k – некоторый калибровочный множитель для учета размера выборки и желаемого доверительного интервала, возможно, определяемый экспертным путем. Не затрагивая проблему выбора опорной вариации Х, можно предложить несколько критериев для определения L-параметра на основе данных о выборочном среднеквадратичном отклонении исходных данных для заданного момента времени:

В оценках (1)–(5) индекс “i” соответствует позиции (номеру) в соответствующем временном ряде, а при группировке анализируемых значений по местному времени LT индекс “i” соответствует номеру дня в году. Суммирование (или поиск максимума) осуществляется по всему окну “бегущего” среднего, ширина которого задается как N. Отметим, что в случае анализа поставляемых NASA глобальных карт ПЭС ионосферы Земли предполагается прямое использование RMS-данных (соответствующих глобальных карт выборочного СКО) из соответствующих ionex-файлов.

Оценки L₁ и L₂ являются пессимистическими и, скорее всего, окажутся малопригодными для практического применения. Параметр L₁ представляет собой прямую аналогию правилу сложения абсолютных погрешностей для значений, заданных с известной точностью, и при этом является оценкой наихудшего случая, когда все погрешности накапливаются и не происходит их частичной взаимной компенсации. Параметр L₂ по факту представляет собой эвклидово расстояние в пространстве ошибок (считая каждую из ошибок независимой). Оценка L₂ является менее пессимистичной, но возможность ее использования на практике также подлежит отдельной проверке. Оценка L_∞ может быть неустойчивой из-за чувствительности к наличию выбросов, и, таким образом, перед ее применением может потребоваться фильтрация анализируемых данных для их сглаживания. Оценки ${{L}_{{{\text{*}}E}}}$ и ${{L}_{{*\zeta }}}$ являются по форме математическим ожиданием и дисперсией. Множитель M в уравнении (5) принимается за N или N – 1 (N – размер выборки) по аналогии со смещенной и несмещенной дисперсией. Мы ожидаем, что оценки ${{L}_{{*E}}}$ и ${{L}_{{*\zeta }}}$ более робастны, но это предположение требует дальнейших исследований.

5. ЗАКЛЮЧЕНИЕ

Несмотря на более чем 40-летнюю историю исследований по проблеме ионосферных предвестников землетрясений (ИПЗ) к настоящему времени не существует общепринятой методики расчета опорной вариации и возмущений ПЭС в целях поиска и выявления ИПЗ. В данной работе предпринята попытка систематически рассмотреть методы расчета фоновой вариации ПЭС, возмущений и определения критериев их аномальности применительно к задачам поиска ИПЗ.

Показано, что распространение в качестве опорной вариации получили следующие данные: 1) ближайшего магнитоспокойного дня или их комбинации, 2) справочной модели ионосферы IRI различных версий или 3) применение различных методов фильтрации. К последним можно отнести вейвлет-анализ, калмановскую фильтрацию, интерквартильный анализ, метод главных компонент, искусственные нейронные сети, генетические алгоритмы, а также всевозможные “скользящие” средние. При расчете последних обычно ширину окна “скользящего” среднего принимают в диапазоне 7–31 день. Такое разнообразие методов приводит к зависимости получаемого результата (т.е. амплитуды и пространственных размеров ИПЗ) от способа вычисления возмущения и опорной вариации и к невозможности прямого количественного сопоставления и обобщения результатов отдельных конкретных ситуаций.

Отдельной проблемой при анализе поставляемых NASA глобальных карт ПЭС ионосферы является отсутствие учета исследователями имеющихся карт выборочных СКО (RMS) для тех же периодов. Для решения последней проблемы нами предложены лебеговы метрики L₁, L₂, L_∞, ${{L}_{{*E}}}$, ${{L}_{{*\zeta }}}$ для оценки “аномальности” вариаций ПЭС, основанные на использовании доступных данных выборочного СКО. На основе анализа свойств предложенных оценок мы ожидаем, что оценки ${{L}_{{*E}}}$ и ${{L}_{{*\zeta }}}$ более робастны.

Авторы благодарят рецензентов за полезные дискуссии, NASA за предоставленные в свободном доступе данные глобальных карт ПЭС ионосферы Земли (ftp://cddis.gsfc.nasa.gov/gnss/products/ ionex/).