Физика Земли, 2022, № 1, стр. 100-117

ВВЕДЕНИЕ

Одной из важнейших характеристик сейсмического процесса является частота возникновения землетрясений разной силы. Зависимость частоты землетрясений от магнитуды называют законом повторяемости землетрясений. В сейсмологии наиболее широко известен закон повторяемости Гутенберга–Рихтера [Gutenberg, Richter, 1944]:

где: $N$ – число землетрясений с магнитудой $m \geqslant M$; $a$ – параметр сейсмической активности; $b$ – наклон графика повторяемости землетрясений.

Закон Гутенберга–Рихтера применим для описания сейсмичности различных тектонических режимов. Параметры $a$ и $b$ имеют принципиальное значение для анализа сейсмической опасности. Кроме того, временное снижение величины $b$ неоднократно наблюдалось перед возникновением сильных землетрясений [Gibowicz, 1973; Завьялов, 1984; 2006; Nuannin et al., 2005]. Временные снижения и нарастания сейсмической активности, известные как сейсмические затишья и форшоковые активизации, также рассматриваются в качестве возможных предвестников сильных землетрясений [Wyss et al., 1999; Zöller et al., 2002].

Определение значений a и b существенно осложняется тем, что обычно график повторяемости землетрясений соответствует теоретической зависимости (1) лишь на ограниченном интервале значений магнитуды (см. рис. 1). Выраженный излом кумулятивного графика повторяемости в области малых магнитуд обычно связывают с неполной регистрацией слабых землетрясений [Писаренко, 1989; Smirnov, 1998; Смирнов, Габсатарова, 2000; Wiemer, Wyss, 2000].

Магнитуда представительной регистрации ${{M}_{c}}$ – это минимальная магнитуда, для которой землетрясения в заданной области регистрируются в полном объеме. Значение ${{M}_{c}}$ соответствует началу линейной части графика повторяемости, поэтому знание ${{M}_{c}}$ становится критически важным для получения адекватных оценок параметров Гутенберга–Рихтера.

Из соотношения (1) следует, что магнитуда имеет экспоненциальное распределение с функцией плотности:

где $\beta = b\ln \left( {10} \right)$.

Существующие методы оценки параметра $b$ включают метод наименьших квадратов [Guttorp, 1987; Nuannin et al., 2005], метод максимального правдоподобия [Aki, 1965; Гусев, 1974; Bender, 1983] и метод, основанный на минимизации статистики критерия согласия Колмогорова [Goldstein et al., 2004]. Прекрасный обзор наиболее распространенных на практике методов оценки $b$ представлен в работе [Marzocchi, Sandri, 2003].

В статье [Bengoubou-Valérius, Gibert, 2013] был выполнен сравнительный анализ этих методов. Авторы исследования пришли к выводу о том, что наиболее надежным является метод максимального правдоподобия. С учетом поправки за дискретизацию магнитуды [Utsu, 1966], оценка максимального правдоподобия $b$ будет иметь вид:

где $\bar {M}$ – среднее выборочное значение магнитуды при $m \geqslant {{M}_{c}}$, $\Delta M$ – интервал группировки, т.е. точность определения магнитуды (в современных каталогах обычно $\Delta M = 0.1$).

Недооценка, т.е. занижение величины ${{M}_{c}}$ приводит к получению некорректных оценок параметров сейсмичности и ошибочной интерпретации данных. Так, включение событий с непредставительной регистрацией в оценку наклона графика повторяемости $b$ может привести к существенному занижению этой оценки (см. рис. 2б), что в свою очередь влечет недооценку частоты возникновения слабых землетрясений и переоценку частоты возникновения сильных землетрясений. Переоценка, т.е. завышение ${{M}_{c}}$ может быть не столь критична [Chen et al., 2011], но приводит к снижению объема анализируемой выборки, а значит повышает неопределенность всех дальнейших результатов анализа. Значение ${{M}_{c}}$ изменяется в пространстве и во времени и эти вариации необходимо учитывать при работе с каталогами землетрясений для оценки различных параметров сейсмического режима.

Таким образом, анализ пространственно-временных вариаций ${{M}_{c}}$ становится неотъемлемым элементом практически любого исследования, опирающегося на каталоги землетрясений.

Одна из первых методик оценки ${{M}_{c}}$ была предложена в работе [Stepp, 1972] и была основана на предположении о том, что последовательность основных землетрясений соответствует стационарному Пуассоновскому процессу. В таком случае средняя частота землетрясений $\lambda $ должна быть стабильна в представительной части каталога, а ее стандартное отклонение обратно пропорционально взятому для оценки $\lambda $ интервалу времени: ${{\sigma }_{\lambda }} = \sqrt {{\lambda \mathord{\left/ {\vphantom {\lambda T}} \right. \kern-0em} T}} $. По графику зависимости ${{\sigma }_{\lambda }}\left( T \right)$ для каждого значения магнитуды можно определить момент, когда график начинает отклоняться от степенной зависимости ${1 \mathord{\left/ {\vphantom {1 {\sqrt T }}} \right. \kern-0em} {\sqrt T }}$, что соответствует окончанию периода представительной регистрации.

Позднее был предложен метод оценки ${{M}_{c}}$, основанный на сравнении числа землетрясений, зарегистрированных в дневное и в ночное время [Rydelek, Sacks, 1989]. Этот метод предназначен для анализа локальных каталогов. Предполагается, что вероятность обнаружения землетрясений возрастает в ночное время вследствие снижения уровня индустриальных шумов. Таким образом, если преобладающее число землетрясений в заданном интервале значений магнитуды регистрируется ночью, каталог считается непредставительным в этом интервале.

В дальнейшем проблеме оценки ${{M}_{c}}$ уделялось существенное внимание и со временем число методов оценки ${{M}_{c}}$ значительно возросло. Подробный обзор современных методов оценки ${{M}_{c}}$ приводится в работе [Mignan, Woessner, 2012]. Большинство этих методов можно отнести к одной из двух категорий:

1. Методы, основанные на анализе каталогов землетрясений.

2. Методы, анализирующие возможности регистрации землетрясений станциями сейсмической сети заданной конфигурации [Ringdal, 1975; Gomberg, 1991].

Сопоставление результатов этих двух подходов выполнялось, в частности, в работе [Смирнов, Габсатарова, 2000], где по данным каталога землетрясений Северного Кавказа было показано, что оценки минимального представительного класса, полученные с помощью двух подходов, хорошо согласуются между собой и адекватно отражают изменения сейсмической сети.

Методы из первой категории менее трудоемкие, требуют меньшего объема информации и чаще применяются на практике. Кроме того, информация, необходимая для использования методов из второй категории, не всегда доступна для исследователей, вследствие чего эти методы не всегда применимы. Поэтому в данной работе рассматриваются только методы из первой категории. Существующее разнообразие приводит к тому, что сам по себе выбор предпочтительного метода оценки ${{M}_{c}}$ перестает быть тривиальным.

В настоящей работе рассмотрены шесть современных методов оценки ${{M}_{c}}$, проанализировано поведение оценок этих методов на синтетических каталогах землетрясений в зависимости от объема выборки и формы исходного распределения числа событий по магнитуде, использованного при создании синтетического каталога. Для анализа вариаций оценок ${{M}_{c}}$ использован метод бутстрап анализа.

МЕТОДЫ ОЦЕНКИ МАГНИТУДЫ ПРЕДСТАВИТЕЛЬНОЙ РЕГИСТРАЦИИ

Оценка ${{M}_{c}}$ в точке максимальной кривизны кумулятивного графика повторяемости землетрясений [Wiemer, Wyss, 2000]. Это самая простая и быстрая процедура оценки ${{M}_{c}}$. На практике точка максимальной кривизны кумулятивного графика повторяемости соответствует магнитуде, на которую приходится максимальное число событий в выборке. Несмотря на простоту и надежность, данная процедура имеет тенденцию занижать значение ${{M}_{c}}$, в особенности в тех случаях, когда дискретный график повторяемости не имеет выраженного максимума. Предположительно, такая форма дискретного графика повторяемости, а именно кривизна непредставительной его части, может указывать на неоднородность анализируемой выборки [Mignan, 2012]. Тем не менее, для однородных выборок такая оценка будет адекватной [Mignan et al., 2011]. Процедура будет обозначаться аббревиатурой MAXC (от Maximum curvature), а соответствующая оценка – $M_{c}^{{{\text{MAXC}}}}$. Пример такой оценки показан на рис. 1.

Оценка ${{M}_{c}}$ посредством расчета критерия согласия [Wiemer, Wyss, 2000]. Данный метод основан на сравнении анализируемой выборки значений магнитуды и созданной на ее основе синтетической выборки. На первом этапе выбирается значение магнитуды нижней отсечки ${{M}_{{co}}}$ (это значение не должно превышать $M_{c}^{{{\text{MAXC}}}}$), по части выборки с магнитудами $m \geqslant {{M}_{{co}}}$ методом максимального правдоподобия оцениваются параметры $a$, $b$ из (1). На основе полученных оценок $a$, $b$ создается синтетическая выборка из распределения Гутенберга–Рихтера. Далее сравниваются кумулятивные графики повторяемости для анализируемой и синтетической выборки и расчитывается критерий согласия:

где ${{B}_{i}}$, ${{S}_{i}}$ – кумулятивное число событий в i-ом магнитудном интервале для анализируемой и синтетической выборки соответственно.

Величина $R$ показывает, какой процент от объема анализируемой выборки можно описать соотношением (1) при данном значении ${{M}_{{co}}}$. Значение ${{M}_{{co}}}$ повышается на величину $\Delta M$, соответствующую точности определения магнитуды, и вычисления повторяются. В качестве оценки ${{M}_{c}}$ выбирается наименьшее значение ${{M}_{{co}}}$, для которого величина критерия согласия $R$ достигает заданного уровня (обычно это 90 или 95%). Процедура будет обозначаться аббревиатурой GFT (от Goodness-of-fit test), а соответствующая оценка – $M_{c}^{{{\text{GFT}}}}$. Пример применения данной процедуры показан на рис. 2a.

Оценка ${{M}_{c}}$ по стабилизации значений b-value. Данная процедура была предложена в работах [Арефьев и др, 1989; Cao, Gao, 2002], в ее основе лежит наблюдение, согласно которому величина $b$ возрастает при ${{M}_{{co}}} < {{M}_{c}}$ и при ${{M}_{{co}}} \gg {{M}_{c}}$, но остается практически неизменной при ${{M}_{{co}}} \geqslant {{M}_{c}}$. Соответственно, в качестве оценки ${{M}_{c}}$ можно взять магнитуду, для которой разность двух последовательных значений $b$ достаточно мала (например, меньше 0.03). В работе [Woessner, Wiemer, 2005] процедура была доработана и был предложен формализованный критерий стабилизации $b$. Модифицированная процедуры выполняется следующим образом: для некоторого диапазона значений ${{M}_{{co}}}$ по формуле (3) оцениваются значения $b$ и рассчитываются среднеквадратичные ошибки этих оценок ${{\sigma }_{b}}$ по формуле [Shi, Bolt, 1982]:

где $\bar {M}$ – среднее выборочное значение магнитуды, $N$ – число событий.

Затем в окне шириной $dM = 5\Delta M$ рассчитываются усредненные значения $b$:

В качестве оценки ${{M}_{c}}$ выбирается наименьшее значение ${{M}_{{co}}}$, для которого разность между текущим и усредненным значениями $b$ не превосходит среднеквадратичную ошибку оценки этого параметра, т.е. выполняется условие $\Delta b = \left| {\bar {b} - b} \right| \leqslant {{\sigma }_{b}}$. Процедура обозначается аббревиатурой MBS (${{M}_{c}}$ by $b$-value stability), а соответствующая оценка – $M_{c}^{{{\text{MBS}}}}$. Пример применения метода показан на рис. 2б.

Оценка ${{M}_{c}}$ методом медианного анализа последовательности значений сдвигов дискретного графика повторяемости [Amorèse, 2007]. Метод основан на итерационной процедуре поиска точек, в которых происходит излом дискретного графика повторяемости. Для этого в каждой точке ${{M}_{i}}$, ($i > 1$) рассчитывается величина сдвига:

где $N\left( {{{M}_{i}}} \right)$ – число событий в интервале ${{M}_{i}} \pm {{\Delta M} \mathord{\left/ {\vphantom {{\Delta M} 2}} \right. \kern-0em} 2}$.

Сдвиги $s\left( {{{M}_{i}}} \right)$ характеризуют направленность и выраженность изменений дискретного графика повторяемости. Каждому значению $s\left( {{{M}_{i}}} \right)$ присваивается ранг ${{R}_{i}}$ так, что наименьшее значение получает минимальный ранг – 1. Затем ранги суммируются в каждой точке $S{{R}_{i}} = \sum\nolimits_{k = 1}^i {{{R}_{k}}} $ и рассчитывается модифицированная сумма рангов $S{{A}_{i}}$:

где $n$ – длина последовательности.

Последовательность $S{{A}_{i}}$ достигает максимума в точке ${{n}_{1}}$, это наиболее вероятная точка излома дискретного графика повторяемости. В этой точке последовательность $s\left( {{{M}_{i}}} \right)$ делится на две части: ${{i}_{1}} = 2, \ldots ,{{n}_{1}}$ и ${{i}_{2}} = {{n}_{1}} + 1, \ldots ,n + 1$. Проверяется гипотеза о том, что в точке ${{n}_{1}}$ нет излома дискретного графика повторяемости, против альтернативной гипотезы о том, что излом есть. Решение принимается на основе непараметрического критерия суммы рангов Уилкоксона [Wilcoxon, 1945], также известного как U-критерий Манна–Уитни [Mann, Whitney, 1947]. Если основная гипотеза отвергается при заданном уровне значимости, точка ${{n}_{1}}$ заносится в список и итерации продолжаются. Процедура может обнаружить несколько значимых точек излома (рис. 3), магнитуде ${{M}_{c}}$ соответствует та точка, для которой вероятность ошибки первого рода при проверке гипотезы минимальна. Процедура обозначается аббревиатурой MBASS (Median-based analysis of the segment slope), а соответствующая оценка – $M_{c}^{{{\text{MBASS}}}}$.

Оценка ${{M}_{c}}$ методом моделирования полного диапазона значений магнитуды [Ogata, Katsura, 1993]. Согласно этому методу, для любого значения $m$ наблюдаемую интенсивность потока событий $\lambda \left( m \right)$ можно представить в виде:

где ${{\lambda }_{0}}\left( m \right) = {{e}^{{ - \beta m}}}$ – теоретическая интенсивность, $q\left( m \right)$ – вероятность обнаружения событий, которая задается нормальным распределением: где $\mu $ и $\sigma $ – среднее и стандартное отклонение.

С учетом сделанных предположений в работе [Ogata, Katsura, 1993] была получена функция плотности распределения магнитуды на всем диапазоне ее значений:

Параметры этой модели оцениваются методом максимального правдоподобия. Магнитуда ${{M}_{c}}$ не входит в модель в явном виде и задается через параметры $\mu $ и $\sigma $: ${{M}_{c}}\left( n \right) = \mu + n\sigma $, $n$ – показатель уровня доверия. При $n = 0$ регистрируются лишь 50% событий с магнитудой $m \geqslant {{M}_{c}}$ при $n = \left( {1,\,\,2,\,\,3} \right)$ регистрируются 84, 98 и 99% событий, соответственно. Вопрос о том, как соотносятся уровни представительности ${{M}_{c}}\left( n \right)$ модели (11) с оценками ${{M}_{c}}$ остальных методов, рассматривался в статье [Mignan, Woessner, 2012], где на синтетических и реальных каталогах землетрясений было показано, что оценки ${{M}_{c}}$ преимущественно попадают в интервал $\left( {\mu ,\mu + \sigma } \right)$, в редких случаях могут достигать уровня $\mu + 2\sigma $, то есть 95%-ой вероятности обнаружения этой модели. Модель (11) обозначается аббревиатурой OK (от Ogata and Katsura), пример применения этой модели показан на рис. 4.

В дальнейшем метод был модифицирован для того, чтобы явно включить в модель ${{M}_{c}}$ [Woessner, Wiemer, 2005]. Была предложена составная модель, в которой нормальное распределение описывает вероятность обнаружения непредставительных событий:

Важно отметить, что для этой модели уравнение (9) не выполняется, интенсивность потока событий задается следующей функцией:

Эта модель будет обозначаться аббревиатурой WW (от Woessner and Wiemer), метод будет обозначаться аббревиатурой EMR (от Entire Magnitude Range), а соответствующая оценка – $M_{c}^{{{\text{EMR}}}}$. Параметры ${{\mu }}$ и ${{\sigma }}$ оцениваются методом нелинейной регрессии по данным выборки в диапазоне $m < {{M}_{c}}$, параметры $a$ и $b$ оцениваются методом максимального правдоподобия по данным в диапазоне $m \geqslant {{M}_{c}}$. На основе этих оценок и уравнения (13) рассчитывается теоретический дискретный график повторяемости на всем диапазоне значений магнитуды и оценивается правдоподобие получения такой выборки. Вычисления повторяются для некоторого диапазона значений ${{M}_{c}}$, итоговая оценка берется в максимуме функции правдоподобия. Такой подход дает наиболее всестороннее описание сейсмичности, но является наиболее затратным в плане времени вычислений. Применение метода показано на рис. 5.

Оценка ${{M}_{c}}$ в точке начала линейной части графика повторяемости [Писаренко, 1989; Смирнов, 2009]. Процедура основана на проверке гипотезы о прямолинейности графика повторяемости. На первом этапе выбирается начальное значение магнитуды нижней отсечки ${{M}_{{co}}}$ (это значение не должно превышать $M_{c}^{{{\text{MAXC}}}}$), по части выборки с магнитудами $m \geqslant {{M}_{{co}}}$ методом максимального правдоподобия оцениваются параметры закона повторяемости. Проверяется гипотеза ${{H}_{1}}$ о прямолинейности графика повторяемости в области $m \geqslant {{M}_{{co}}}$ против гипотезы ${{H}_{2}}$ о том, что график повторяемости в этой области нелинейный. Для этого рассчитываются наблюдаемые частоты попадания значений магнитуды в интервалы ${{M}_{i}} \pm {{\Delta M} \mathord{\left/ {\vphantom {{\Delta M} 2}} \right. \kern-0em} 2}$, где ${{M}_{i}} = {{M}_{{co}}} + i\Delta M,$ $i = 0, \ldots ,n$:

где ${{N}_{i}}$ – число событий в интервале ${{M}_{i}} \pm {{\Delta M} \mathord{\left/ {\vphantom {{\Delta M} 2}} \right. \kern-0em} 2}$, ${{Q}_{0}} = \sum\nolimits_{i = 0}^n {{{N}_{i}}} $.

И теоретические частоты, соответствующие закону Гутенберга–Рихтера с параметром $\hat {b}$:

Рассчитывается выборочная статистика, характеризующая отклонение гипотезы ${{H}_{2}}$ от гипотезы ${{H}_{1}}$:

Если верна гипотеза ${{H}_{1}}$, то величина $2{{Q}_{0}}\hat {I}$ имеет асимптотически при ${{Q}_{0}} \to \infty $ распределение ${{{{\chi }}}^{2}}$ с $\left( {n - 1} \right)$ степенью свободы. Для проверки гипотезы ${{H}_{1}}$ определяется уровень значимости $\mu $ величины $2{{Q}_{0}}\hat {I}$. Если $\mu $ меньше заданного уровня (например, 0.01 или 0.05), то гипотеза ${{H}_{1}}$ отвергается. В этом случае на следующем этапе решается вопрос о причинах нелинейности графика повторяемости. Если же гипотеза ${{H}_{1}}$ принимается, то ${{M}_{{co}}}$ – искомая оценка ${{M}_{c}}$ и процедура на этом завершается.

Вероятность регистрации землетрясений с $m = {{M}_{{co}}}$ обозначается как $p$, если $p = 1$ регистрация представительна, если $p < 1$ имеются пропуски землетрясений с $m = {{M}_{{co}}}$. Для решения вопроса о представительности регистрации землетрясений с $m = {{M}_{{co}}}$ строится оценка максимального правдоподобия $\hat {p}$ для $p$ и проверяется гипотеза об отличии этой оценки от 1. Оценка $\hat {p}$ определяется отношением наблюденной величины ${{N}_{0}}$ к оценке числа событий с $m = {{M}_{{co}}}$, построенной по старшим магнитудам с $m \geqslant {{M}_{{co}}} + \Delta M$. Пусть $\tilde {b}$ – оценка максимального правдоподобия, полученная по части выборки с $m \geqslant {{M}_{{co}}} + \Delta M$, ${{\tilde {x}}^{i}} = {{10}^{{ - \tilde {b}\left( {{{M}_{i}} - {{M}_{{co}}}} \right)}}}$. Тогда оценка $\hat {p}$ имеет вид:

а ее асимптотическая дисперсия оценивается как: где ${{{{\psi }}}_{k}} = \sum\nolimits_{i = 1}^n {{{i}^{k}}{{{\tilde {x}}}^{i}}} ,~\,\,\,\,k = 0,\,\,1,\,\,2.$

Для проверки гипотезы ${{H}_{0}}$ о представительной регистрации землетрясений с $m = {{M}_{{co}}}$ рассчитывается величина $t = \frac{{1 - \hat {p}}}{{{{{\left( {Var\left( {\hat {p}} \right)} \right)}}^{{{1 \mathord{\left/ {\vphantom {1 2}} \right. \kern-0em} 2}}}}}}$, имеющая распределение Стьюдента, и определяется ее уровень значимости ${{{{\mu }}}_{0}}$. Гипотеза ${{H}_{0}}$ отвергается если ${{{{\mu }}}_{0}} \geqslant 1 - {{\varepsilon }}$, где ${{\varepsilon }}$ – некоторое малое число (например 0.05 или 0.1). В этом случае значение ${{M}_{{co}}}$ повышается на $\Delta M$ и вся процедура повторяется. Если гипотеза ${{H}_{0}}$ принимается, то ${{M}_{{co}}}$ – искомая оценка ${{M}_{c}}$. Метод будет обозначаться аббревиатурой LLS (от Lower end of linear segment), а соответствующая оценка – $M_{c}^{{{\text{LLS}}}}$. Пример применения метода показан на рис. 6.

Бутстрап [Efron, Tibshirani, 1993] – непараметрический метод исследования распределений искомых параметров, основанный на многократном извлечении методом Монте-Карло повторных выборок из имеющейся выборки. Метод позволяет быстро и просто оценивать разнообразные статистики (например, дисперсию, доверительные интервалы) для сложных моделей. Суть метода состоит в том, что из имеющейся выборки случайным выбором с повторением формируется некоторое множество ${{n}_{b}}$ повторных выборок заданного размера. На множестве повторных выборок оцениваются искомые параметры и по полученным эмпирическим распределениям определяются все необходимые статистики.

В частности, бутстрап применялся для анализа дисперсии оценок ${{M}_{c}}$ в работах [Woessner, Wiemer, 2005; Amorèse, 2007; Mignan et al., 2011; Mignan, Woessner, 2012]. В книге [Chernick, 1999] в качестве рекомендованного числа повторных выборок, достаточного для надежной оценки дисперсии, приводится цифра ${{n}_{b}} = 100$, однако при наличии достаточных вычислительных мощностей автор рекомендует использовать ${{n}_{b}} = 1000$. В работе [Woessner, Wiemer, 2005] авторы отмечают, что оценки дисперсии ${{M}_{c}}$ стабилизируются при ${{n}_{b}} \geqslant 200$. Эту рекомендацию использовали в работах [Mignan et al., 2011; Mignan, Woessner, 2012]. В работе [Amorèse, 2007] автор использовал ${{n}_{b}} = 1000$ повторных выборок. В настоящей работе бутстрап применяется с ${{n}_{b}} = 500$ повторных выборок.

АНАЛИЗ ПРЕДСТАВИТЕЛЬНОСТИ ВЫБОРОК ИЗ ИНСТРУМЕНТАЛЬНЫХ КАТАЛОГОВ ЗЕМЛЕТРЯСЕНИЙ

Все описанные процедуры были программно реализованы В.А. Павленко в MATLAB. В качестве проверки корректности работы процедур была предпринята попытка воспроизвести результаты, полученные в работах [Woessner, Wiemer, 2005; Amorèse, 2007]. Для этого были взяты те же выборки из каталогов землетрясений, которые использовали авторы этих работ в качестве демонстрационного материала:

1. Выборка из регионального каталога Швейцарской сейсмологической службы, the Earthquake Catalog of Switzerland (ECOS; http://ecos09.seismo.ethz.ch/query) за 1992–2002 гг.

2. Выборка из глобального каталога the Global Centroid Moment Tensor Catalog (CMT; https://www.globalcmt.org/CMTsearch) за 1983–2002 гг.

3. Выборка из регионального каталога землетрясений Северной Калифорнии, the Northern California Seismic Network (NCSN; https://www.ncedc.org/ncedc/catalog-search) за 1998–2000 гг.

4. Выборка из каталога землетрясений вулканической зоны сейсмичности Канто, the National Research Institute for Earth Science and Disaster Prevention Earthquake Catalog (NIED; http://evrrss.eri.u-tokyo.ac.jp/db/nied) за 1992–2002 гг.

Параметры получившихся выборок сведены в табл. 1, построенные по выборкам графики повторяемости показаны на рис. 7. Сравнение графиков повторяемости и объемов выборок с представленными в работах [Woessner, Wiemer, 2005; Amorèse, 2007], показывает, что ни одна из полученных выборок не совпадает в полной мере с теми выборками, которые использовали эти авторы. Вероятнее всего, причиной различий стали изменения, внесенные в каталоги со времени публикации этих работ.

По каждой выборке с помощью описанных процедур и метода бутстрап оценивались средние значения и дисперсии оценок ${{M}_{c}}$. Полученные результаты представлены в табл. 1. Результаты для выборок из каталогов ECOS, NCSN и NIED хорошо согласуются с результатами, представленными в работах [Woessner, Wiemer, 2005; Amorèse, 2007]. Оценки ${{M}_{c}}$, полученные по выборке из каталога CMT, существенно отличаются от тех оценок, которые представлены в работе [Woessner, Wiemer, 2005], но это расхождение, по всей видимости, обусловлено различиями выборочных данных.

В этих результатах проявляются некоторые характерные особенности рассматриваемых методов: методы MAXC и GFT обычно дают самые низкие оценки ${{M}_{c}}$, оценки метода MBS обычно оказываются наиболее консервативными. Оценки остальных методов занимают промежуточные позиции.

Из табл. 1 видно, что для всех методов, кроме метода MBS, дисперсии оценок ${{M}_{c}}$, полученных по выборкам из каталогов NCSN и CMT, оказывались минимальными или близкими к минимальным. Форма распределения выборки из каталога NCSN близка к характерной угловой форме с выраженным максимумом, которая присуща выборкам с равномерным уровнем регистрации. Распределение выборки из каталога CMT имеет сглаженный максимум, для этой выборки разброс оценок ${{M}_{c}}$ заметно выше, чем для выборки из каталога NCSN.

Распределения выборок из каталогов ECOS и NIED имеют сглаженную форму с заметной кривизной непредставительной части дискретного графика повторяемости, которая может свидетельствовать о неоднородности выборочных данных. В целом это приводит и к возрастанию дисперсий оценок ${{M}_{c}}$ и к более значительным различиям в оценках ${{M}_{c}}$, полученных разными методами.

Полученные результаты демонстрируют чувствительность методов оценки ${{M}_{c}}$ к форме распределения анализируемой выборки. Этот эффект более детально рассмотрен на синтетических каталогах землетрясений в следующих секциях.

МОДЕЛИ РАСПРЕДЕЛЕНИЯ МАГНИТУДЫ

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

Модель, предложенная в работе [Mignan, 2012], описывает элементарный график повторяемости, построенный по выборке значений магнитуды, в которой пространственно-временные неоднородности ${{M}_{c}}$ сведены к минимуму. Вероятность обнаружения событий при $m < {{M}_{c}}$ изменяется по экспоненте:

где: $\kappa = k\ln \left( {10} \right) > {{\beta }}$ – параметр обнаружения; ${{\beta }}$ – показатель экспоненциального распределения (2). Модель будет обозначаться аббревиатурой AN (от английского Angular).

В статье [García-Hernández et al., 2019] авторы предложили модель, в которой вероятность обнаружения событий задается полиномом второй степени (модель POL):

где ${{M}_{i}}$ – минимальная магнитуда обнаружения.

Такая форма распределения характерна для ранних периодов инструментальных наблюдений, когда ввиду сравнительно невысокой чувствительности приборов слабые события не регистрировались. Обе модели AN и POL подчиняются уравнению (9).

Авторы работы [García-Hernández et al., 2019] проанализировали формы выборочных распределений магнитуды в зависимости от объема выборки и радиуса области, из которой выбирались землетрясения. Для анализа были использованы данные Международного Сейсмологического Центра [ISC; www.isc.ac.uk/iscbulletin/search/catalogue/; Storchak et al., 2013], модели оценивались по критерию согласия Колмогорова и информационному критерию Акаике [Akaike, 1974]. Результаты этого исследования показали, что для всех моделей степень согласия с данными снижалась с возрастанием объема выборки, но практически не зависела от радиуса области выборки событий. Выраженное снижение согласия с данными наблюдалось для выборок, содержащих более 1000 событий. В целом модель OK чаще других оказывалась наилучшей в плане аппроксимации данных; модель AN показывала хорошие результаты для малых выборок объемом 100–300 событий; модель WW показывала хорошие результаты для более объемных выборок, включающих 1000–3000 событий; модель POL оказалась существенно хуже, показав относительно неплохой результат только для выборок объемом 100 событий.

Таким образом, для создания синтетических каталогов были выбраны модели WW, AN и POL. Модель WW описывает наиболее распространенную ситуацию, при которой выборка значений магнитуды неоднородна по ${{M}_{c}}$. Модель AN описывает выборку, в которой неоднородности ${{M}_{c}}$ сведены к минимуму. Модель POL описывает выборку, из которой удалены слабые землетрясения с магнитудами $m \leqslant {{M}_{i}}$. Различие этих трех моделей продемонстрировано на рис. 8.

ТЕСТИРОВАНИЕ МЕТОДОВ ОЦЕНКИ ${{M}_{c}}$ НА СИНТЕТИЧЕСКИХ КАТАЛОГАХ ЗЕМЛЕТРЯСЕНИЙ

При выборе значений параметров модели WW учитывались их оценки, полученные по выборкам из каталогов ECOS, CMT, NCSN и NIED (табл. 1). Из этих результатов видно, что магнитуда $\mu $ располагается примерно посередине между минимальной магнитудой в выборке и оценкой $~M_{c}^{{{\text{EMR}}}}$, а ${{\sigma }}$ имеет величину порядка 0.1–0.3.

Выборки из каталогов ECOS и NIED имеют типичную для модели WW форму распределения (рис. 7а, 7г). При фиксированном ${{\mu }}$, снижение ${{\sigma }}$ приводит к тому, что левый хвост распределения приобретает более крутой наклон, максимум распределения сглаживается, форма распределения приближается к форме выборки из каталога CMT (рис. 7б). При смещении ${{\mu }}$ ближе к ${{M}_{c}}$ при фиксированном ${{\sigma }}$ форма распределения становится близка к форме выборки из каталога NCSN (рис. 7в) с острым максимумом и крутым наклоном левого хвоста распределения.

На основе этих трех вариантов формы распределения WW были выбраны следующие комбинации значений параметров для создания синтетических каталогов: ${{{{\mu }}}_{1}} = 0.5,~\,\,{{{{\sigma }}}_{1}} = 0.25$ (такие значения параметров использовали для создания синтетического каталога авторы [Woessner, Wiemer, 2005]), ${{\mu }_{2}} = 0.5,\,\,~{{\sigma }_{2}} = 0.1$ и ${{\mu }_{3}} = 0.75,~\,\,{{\sigma }_{3}} = 0.25$.

В работе [Mignan, 2012] приводятся оценки параметра $k$ модели AN, полученные по данным каталогов землетрясений Невады и Южной Калифорнии. Для Невады оценки $k$ лежат в пределах от 1.4 до 5.4, максимум плотности распределения приходится на $\hat {k} = 3$, для Южной Калифорнии оценки лежат в пределах от 2 до 5 с максимумом плотности при $\hat {k} = 2.9$. По выборке из каталога NCSN была получена оценка максимального правдоподобия $\hat {k} = 2.95$. С повышением $k$ наклон левого хвоста распределения становится все более крутым, а число событий в непредставительной части распределения снижается. Для создания синтетических каталогов с распределением AN были использованы значения ${{k}_{1}} = 2$, ${{k}_{2}} = 3$, ${{k}_{3}} = 4$.

Значения параметра ${{M}_{i}}$ модели POL были заданы произвольно: ${{M}_{{i1}}} = 0.5$, ${{M}_{{i2}}} = 0.6$, ${{M}_{{i3}}} = 0.7$.

Для представительной части каталогов для всех трех моделей распределения магнитуды были заданы значения параметров ${{M}_{c}} = 1$, $b = 1$. Для каждой модели и для каждой комбинации значений параметров был сгенерирован синтетический каталог объемом 10 000 событий.

Далее по синтетическим каталогам с помощью метода бутстрап оценивались зависимости двух первых моментов распределений ${{M}_{c}}$ от объема выборки $N$. Объем выборки варьировался в диапазоне от ${{N}_{{{\text{min}}}}} = 20$ до ${{N}_{{{\text{max}}}}} = 2000$. Нижняя граница диапазона соответствует минимальному объему выборки, позволяющему применять описанные процедуры оценки ${{M}_{c}}$. Верхняя граница выбиралась как достаточно большое значение объема выборки, позволяющее оценить асимптотические свойства оценок ${{M}_{c}}$ при $N \to \infty $.

Как известно, состоятельной называется оценка, которая сходится по вероятности к оцениваемому параметру. Признаком состоятельной оценки служит асимптотическая несмещенность и убывание дисперсии с ростом объема выборки. Таким образом, от оценок ${{M}_{c}}$ ожидалось, что с ростом объема выборки средние значения будут приближаться к истинному значению ${{M}_{c}}$, а их дисперсии будут снижаться.

РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЕ

Результаты тестирования методов оценки ${{M}_{c}}$ на синтетических каталогах показаны на рис. 9–рис. 11. На этих рисунках слева направо по столбцам показаны результаты для каталогов с номерами 1, 2, 3. Сверху вниз по строкам показаны результаты для методов MAXC, GFT, MBS, MBASS, EMR и LLS. На графиках горизонтальной сплошной чертой показано истинное значение ${{M}_{c}}$, точками показаны средние значения оценок ${{M}_{c}}$, пунктирными линиями обозначены границы одного стандартного отклонения от среднего.

Результаты для модели WW показаны на рис. 9. Оценка $M_{c}^{{{\text{MAXC}}}}$ для каталога WW₁ медленно сходится к истинному значению ${{M}_{c}}$ и, следовательно, является состоятельной. Оценка $M_{c}^{{{\text{MAXC}}}}$ для каталога WW₃ также является состоятельной. Оценка $M_{c}^{{{\text{MAXC}}}}$ для каталога WW₂ недооценивает истинное значение ${{M}_{c}}$ на величину порядка 0.12 и имеет неубывающую дисперсию в силу особенностей формы распределения.

Оценка $M_{c}^{{{\text{GFT}}}}$ для всех трех каталогов WW имеет при $N \to \infty $ асимптотическое смещение размером порядка –0.1, –0.11 и –0.05 соответственно. Методы MBS и LLS дают схожие результаты, оба метода несколько недооценивают ${{M}_{c}}$ на малых выборках, но с ростом $N$ оценки постепенно сходятся к истинному значению ${{M}_{c}}$.

Результаты метода MBASS очень похожи на результаты метода MAXC, однако оценка $M_{c}^{{{\text{MBASS}}}}$ для каталога WW₂ все же оказывается состоятельной. Метод EMR дает состоятельные оценки ${{M}_{c}}$ для каталогов WW₁ и WW₂. Оценка $M_{c}^{{{\text{EMR}}}}$ для каталога WW₃ переоценивает ${{M}_{c}}$ на малых выборках, но с ростом $N$ постепенно сходится к истинному значению ${{M}_{c}}$.

Результаты для модели AN показаны на рис 10. Метод MAXC дает состоятельные оценки ${{M}_{c}}$ для всех трех каталогов AN. Метод GFT занижает значение ${{M}_{c}}$ для каталогов AN₁ и AN₂ на величину порядка 0.08 и 0.02 соответственно. Для каталога AN₃ метод GFT дает состоятельную оценку ${{M}_{c}}$. Методы MBS и LLS вновь демонстрируют схожие результаты, однако оценка $M_{c}^{{{\text{MBS}}}}$ для всех трех каталогов AN имеет при $N \to \infty $ малое асимптотическое смещение размером порядка 0.02.

Оценки $M_{c}^{{{\text{MBASS}}}}$ в среднем соответствуют оценкам $M_{c}^{{{\text{MAXC}}}}$, но имеют более высокую дисперсию. Для всех трех каталогов AN метод EMR переоценивает ${{M}_{c}}$ на малых выборках, но с ростом $N$ оценка $M_{c}^{{{\text{EMR}}}}$ постепенно сходится к истинному значению ${{M}_{c}}$. Результаты слабо зависят от параметра $k$.

Результаты для модели POL показаны на рис. 11. Описанные методы оказываются не в состоянии корректно определить значение ${{M}_{c}}$ для модели POL. Особенно ярко это демонстрируют результаты метода MBASS, который катастрофически завышает значение ${{M}_{c}}$ для всех трех каталогов POL. Ни одна из оценок не является состоятельной. Результаты слабо зависят от параметра ${{M}_{i}}$.

Полученные результаты показывают, что все описанные методы хорошо справляются с задачей определения ${{M}_{c}}$ для модели AN, имеющей резкий переход между непредставительной и представительной частями выборки. Несколько хуже это получается в случае модели WW, со сглаженным максимумом распределения. Наконец, в случае модели POL, у которой этот переход практически никак не выражен, описанные методы оказываются не в состоянии корректно определить ${{M}_{c}}$.

С одной стороны, такой результат подсказывает, что следует с особым вниманием подходить к анализу представительности данных, относящихся к началу инструментальных каталогов. С другой стороны, наблюдения [García-Hernández et al., 2019], согласно которым модель POL в реальных каталогах встречается достаточно редко, дают надежду на то, что, возможно, на практике с таким распределением сталкиваться не придется.

Метод MAXC оказывается незаменимым при анализе представительности выборок из малых пространственно-временных объемов, для которых характерной является форма распределения AN. Это единственный метод, применимый для анализа представительности выборок экстремально малого объема $N \geqslant 4$ [Mignan et al., 2011]. Этот метод применялся для создания карт пространственных вариаций ${{M}_{c}}$ высокого разрешения на Тайване [Mignan et al., 2011]. Для более объемных выборок, включающих сотни и более событий, в которых максимум распределения магнитуды обычно сглажен, метод MAXC занижает ${{M}_{c}}$ и не рекомендуется к использованию.

Методы GFT и MBS – несложные методы, дающие реалистичные оценки ${{M}_{c}}$. Подходят для анализа представительности выборок среднего размера, содержащих от нескольких сотен до нескольких тысяч событий. Метод GFT применялся для картирования ${{M}_{c}}$ на Аляске, в Западных США и Японии [Wiemer, Wyss, 2000]. Учитывая, что метод GFT, как правило, дает более мягкую оценку ${{M}_{c}}$, а метод MBS – более консервативную, следует выбирать один из них, основываясь на тех требованиях к полноте каталога, которые предъявляет поставленная задача. Применять эти методы к более объемным выборкам следует с большой осторожностью.

Метод MBASS не обладает преимуществами перед более простым непараметрическим методом MAXC. В тех случаях, когда метод MBASS корректно определяет ${{M}_{c}}$, оценки оказываются близки к тем, что дает метод MAXC. При этом дисперсия оценок $M_{c}^{{{\text{MBASS}}}}$ как правило оказывается выше. По этой причине метод MBASS не рекомендуется к использованию.

Метод EMR долгое время позиционировался как лучший метод для анализа представительности каталогов землетрясений. В частности, он применялся для картирования значений ${{M}_{c}}$ по мировым каталогам ISC и CMT [Woessner, Wiemer, 2005], а также по каталогам землетрясений Южной Калифорнии [Hutton et al., 2010] и Японии [Nanjo et al., 2010]. Однако в итоге сами авторы метода EMR признали его неудачным [Mignan, Woessner, 2012] и призвали отказаться от его дальнейшего использования.

Метод LLS демонстрирует замечательные результаты: оценки $M_{c}^{{{\text{LLS}}}}$ достаточно быстро сходятся к истинному значению ${{M}_{c}}$ и имеют сравнительно невысокую дисперсию. Этот метод был внедрен в практику в организациях ЕГС РАН и применялся для картирования пространственно-временных неоднородностей ${{M}_{c}}$ по каталогам землетрясений Новой Зеландии, Греции, Камчатки, Кавказа, Киргизии и Северного Китая [Smirnov, 1998]. При этом, судя по публикациям, этот метод остается практически неизвестным за рубежом. Метод LLS можно использовать для анализа представительности выборок среднего и большого размера.

В заключение стоит отметить, что картирование ${{M}_{c}}$ в областях с высокими вариациями уровня регистрации, таких как Камчатка или Япония, должно опираться на анализ локальных распределений магнитуды, как было сделано, например, в работе [Mignan, 2012].

ВЫВОДЫ

В работе выполнен сравнительный анализ шести методов оценки магнитуды представительной регистрации землетрясений ${{M}_{c}}$. Для сравнения методов использовались как выборки из реальных инструментальных каталогов, так и синтетические каталоги землетрясений, сгенерированные на основе трех моделей распределения магнитуды. Для создания синтетических каталогов были использованы модели, соответствующие выборочным распределениям, встречающимся при анализе реальных данных инструментальных каталогов: модель неоднородной выборки, модель однородной выборки, модель выборки, из которой удалены слабые землетрясения. Показано, что результаты применения описанных методов оценки ${{M}_{c}}$ в значительной мере зависят от формы распределения и объема анализируемой выборки. По результатам анализа сформулированы рекомендации по выбору подходящих методов оценки ${{M}_{c}}$.

БЛАГОДАРНОСТИ

Авторы благодарят профессора В.Б. Смирнова за помощь в реализации метода LLS.

Авторы благодарят рецензентов Д.А. Сторчака и И.А. Воробьеву, а также анонимного куратора статьи за их ценные замечания.