ФИЗИКА ЗЕМЛИ, 2004, № 7, с. 50-56
УДК 550.34.113
ИСПОЛЬЗОВАНИЕ ВЕЙВЛЕТ-АНАЛИЗА ДЛЯ АВТОМАТИЧЕСКОЙ КЛАССИФИКАЦИИ ТРЕХКОМПОНЕНТНЫХ СЕЙСМИЧЕСКИХ ЗАПИСЕЙ
© 2004 г. А. А. Лшбушин (мл.)1, 3. Калаб2, Н. Частова3
1Институт физики Земли им. Г.А. Гамбурцева РАН, г. Москва 2Институт геоники Чешской Академии наук, г. Острава, Чехия 3Технический университет Остравы, Чехия Поступила в редакцию 30.01.2003 г.
Предлагается метод автоматической классификации 3-компонентных сейсмических записей, основанный на использовании информационного признака - уровня вейвлет-сжатия Донохо-Джонсто-на. Метод иллюстрируется на примере анализа массива 3-компонентных записей сейсмических событий в шахтах Силезского угольного бассейна (Чехия). Выделенные кластеры событий предположительно различаются по механизму подвижки в эпицентре события.
ВВЕДЕНИЕ
Проблема классификации элементов большого массива сейсмических записей является актуальной задачей в различных областях сейсмологии и сейсморазведки. Для ее решения традиционно применяют широкий набор статистических подходов, опирающихся на идеи спектрального анализа, распознавания образов, синтаксических процедур и т.д. [Seismic signal..., 1982]. В задачах классификации эффективность в решающей степени определяется удачным выбором вектора признаков небольшой размерности, характеризующих классифицируемый объект [Айвазян и др., 1989; Вапник, Червоненкис, 1974; Duda, Hart, 1973]. Поэтому вопрос о выборе признаков, характеризующих сейсмическую запись, является ключевым при попытках разделить множество сигналов на отдельные подмножества - кластеры.
В данной статье в качестве признака, характеризующего скалярную сейсмическую запись, предлагается рассматривать так называемый уровень вейвлет-сжатия Донохо-Джонстона [Dono-ho, Johnstone, 1994] - величину а, изменяющуюся в пределах от 0 до 1 и задающую ту часть коэффициентов разложения сигнала по базису ортогональных финитных функций (вейвлетов), которая может быть отброшена без существенной потери информации о сигнале. Величина уровня а зависит от выбора базиса вейвлет-функций. Поэтому вычислению признака а должна предшествовать процедура определения базиса, по которому будет раскладываться сигнал. В статье в качестве такой процедуры использовался метод последовательного когерентного вейвлет-сжатия (coherent basis thresholding) [Berger, Coifman, Goldberg, 1994; Mallat, 1998], основанный на критерии
минимума энтропии распределения квадратов вейв-лет-коэффициентов в остаточном сигнале при последовательном удалении из него наиболее информативных составляющих.
При наличии примерно одинакового уровня шума в анализируемом множестве сигналов, признак а характеризует "насыщенность" сигнала разнообразными элементами поведения: чем ближе величина а к 1, тем сигнал "проще" и может быть адекватно описан меньшим количеством коэффициентов оптимального вейвлет-бази-са. Таким образом, от каждой трехкомпонентной сейсмической записи получается 3-мерный вектор безразмерных признаков - уровней сжатия каждой компоненты. Далее полученное облако 3-мерных векторов подвергается стандартной процедуре кластерного анализа методом итерационной минимизации функционала компактности разбиения множества векторов на заданное пробное число кластеров [Duda, Hart, 1973]. Оптимальное число кластеров определяется из условия максимума псевдо-Б-статистики [Vogel, Wong, 1978]. В качестве примера применения рассмотрена задача классификации массива из 111 трех-компонентных записей сейсмических событий в шахтах Силезского угольного бассейна в Чехии.
Отметим, что вейвлет-анализ [Chui, 1992; Daubechies, 1992; Mallat, 1998] по своему смыслу является более адекватным инструментом для анализа и классификации нестационарных сигналов, чем традиционно применяемый Фурье-анализ. Это связано с основным свойством используемых базисных функций - их финитностью, что позволяет анализировать резко нестационарные и негауссовы сигналы. Для решения многих задач геофизики использование вейвлетов дает новый
и свежии взгляд на свойства и структуру данных [Любушин (мл.), 2000; 2001; 2002; Cästovaä, Kaläb, 1999; Cästovä, Kaläb, KuCera, 1999].
ОПИСАНИЕ МЕТОДА КЛАССИФИКАЦИИ
Пусть с<к - веИвлет-коэффициенты анализируемого сигнала x(t) (t = 1, ..., N - дискретное время) при разложении по системе ортогональных финитных базисных функциИ. ВерхниИ индекс k соответствует номеру уровня детальности веИв-лет-разложения, а нижниИ индекс j - центру временной окрестности. Максимально возможное значение m номера детальности зависит от объема анализируемоИ выборки. Подробно об используемых обозначениях можно узнать из [Любушин (мл.), 2000; 2001; 2002]. Здесь мы использовали словарь из 17 веИвлетов: 10 обычных ортогональных веИвлетов Добеши с порядками от 2 до 20 (использование более высоких порядков сопряжено с численноИ неустоИчивостью) и 7 так называемых "симлетов" - модификациИ веИвлетов Добеши, в которых форма базисных функциИ является более симметричноИ, чем для обычных веИвлетов [Chui, 1992; Daubechies, 1992; Mal-lat, 1998]. Симлеты обладают теми же своИствами компактности, ортогональности, полноты и гладкости, что и веИвлеты, но для порядков от 2-го до 6-го они совпадают с обычным ортогональным базисом Добеши, а затем, для порядков от 8-го до 20-го, появляются различия в форме базисноИ функции. Вследствие этого общее число используемых вариантов ортогональных компактных базисных функциИ равно 17.
При выборе оптимального веИвлет-базиса наиболее часто используется критериИ минимума энтропии распределения квадратов модулеИ веИв-лет-коэффициентов:
п( m - k) m 2
E(x) = -11 PTln(Pfmin'
k = ij = i (1)
j = j 2/Ц
„( l)l2
l, i
Метод (1) подбирает для сигнала х(0 такой базис, в котором распределение его вейвлет-коэффици-ентов максимально отличается от равномерного. Тем самым максимум информации сосредотачивается в минимальном количестве коэффициентов разложения. Чаще всего простое применение критерия (1) дает вполне удовлетворительные результаты. Однако ниже мы использовали более изощренный метод подбора оптимального базиса, который имеет вид итерационной процедуры, многократно использующей критерий (1). Это связано с желанием выделить, путем применения того или иного базиса, как можно более тонкие разли-
чия в структуре сигналов. Метод был предложен в работе [Berger, Coifman, Goldberg, 1994] для решения задачи очищения старых вокальных записей оперных классиков от характерных шумов (шипение, трески, щелчки) и был назван методом последовательного когерентного отсечения (coherent basis thresholding). Кратко изложим метод в виде последовательности операций:
1) инициализация: в рабочий буфер x(t) помещаем исходный сигнал x0(t);
2) определяем порядок вейвлета из критерия (1) для сигнала x(t): E(x) —► min;
3) отсортируем вейвлет-коэффициенты cj) сигнала x(t) для базиса, определенного в пункте 2, в порядке убывания их абсолютных величин и отсортированные коэффициенты обозначим через dj, j = 0, 1, ..., (N - 1). Таким образом, в отсортированной последовательности коэффициент d0 является максимальным по модулю;
4) определим минимальное целое число M = 0, 1, ..., (N - 1) из условия выполнения неравенства:
\d
м
N-1
<
I \dj
2 In ( N - M); (N - M) '
(2)
] = м +1
5) если условие (2) выполняется сразу, для значения М = 0, то считать, что оптимальный порядок найден, и выйти из алгоритма;
6) если ни для какого М = 0, 1, ..., N - 1) условие (2) не выполняется, то положить значение оптимального порядка вейвлета равным найденному в пункте 2 из условия минимума энтропии сразу после инициализации и выйти из алгоритма;
7) обнулить все коэффициенты с^), для которых \с) ^ > \йм\, совершить обратное вейвлет-пре-образование с оставшимися коэффициентами, получившийся остаточный сигнал поместить в рабочий буфер х(0 и перейти к пункту 2.
Смысл этой процедуры состоит в следующем. Считается, что сигнал состоит из "полезного сигнала", вариации которого отражены в значениях вейвлет-коэффициентов, достаточно больших по модулю, и из "шума", которому соответствуют все прочие коэффициенты. Задача состоит в выборе порога значений модулей коэффициентов, выше которого они отвечают "полезному сигналу", а ниже - "шуму". Неравенство (2) как раз и призвано определять такой порог. Это условие взято из формулы для вероятности асимптотических максимальных уклонений значений гауссовского белого шума Б(г) (например, [Королюк и др., 1985]):
N-1
lim Prob\ max |B(t)|2/1 |B(j)|
N 0 < t <(N-1)
2 < 2 ln N < N I
= 1. (3)
2
2
52
ЛЮБУШИН (мл.) и др.
Кроме того, в дальнейшем нам будет полезна еще одна формула, непосредственно вытекающая из (3):
lim Prob { max \B (t )\<aj 2ln N } = 1, (4)
N 0 < t < (N-1)
где а - значение стандартного отклонения гаус-совского белого шума B(t).
Таким образом, смысл условия (2) состоит в разделении вейвлет-коэффициентов на "шумовые" и "полезные". Шумообразующими считаются те достаточно малые по модулю (нижний предел суммирования в (2) равен M + 1) коэффициенты, максимальные абсолютные значения которых лежат в асимптотических пределах для белого шума (формула (3)). Но такое выделение "шума" из сигнала зависит от используемого базиса (то, что для одного базиса является "шумом", для другого может не удовлетворять критерию (2)). Поэтому, выбрав базис из условия минимума энтропии (пункт 2), далее выделяем "шум" по отношению к выбранному базису (пункт 7) и уже для этого остаточного (обедненного информацией) сигнала опять ищется оптимальный базис (пункт 2) и так далее, до тех пор, пока остаточный сигнал не будет являться шумом даже по отношению к своему оптимальному базису (пункт 5). Последний определенный оптимальный базис считается оптимальным, поскольку он смог хоть что-то определить в самом обедненном остаточном сигнале. В этом есть смысл, поскольку из начального сигнала x0(t) можно извлечь информацию с помощью любого базиса, тогда как для обедненного остаточного сигнала надо поискать наилучший. Заметим, что для обработанных ниже сейсмических записей в 90% случаях оптимальный базис, найденный в результате такой итерационной
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.