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

Область техники
Устройство индикации изменения риска сильных землетрясений по результатам многоканального наблюдения с перерывами (УИР) представляет собой прибор, предназначенный для обнаружения реперных точек при анализе временных рядов (BP) концентрации радона [11], с целью краткосрочного прогнозирования сильных землетрясений.
Устройство может быть использовано в областях практической деятельности, где возникает задача наблюдения за некоторой совокупностью BP с целью выявления аномального поведения объекта (объектов).
УИР может применяться:
- при выявлении предвестников опасных природных явлений и аномального поведения технических объектов с целью предотвращения техногенных катастроф;
- для диагностики состояния операторов в критически важных областях ручного управления,
а также может повысить качество и расширить функциональные возможности:
- различных охранных систем;
- систем медицинской техники, используемых для наблюдения за состоянием пациентов и т.д.
Уровень техники. Аналоги и их недостатки
Для анализа BP, среди других, применяются анализаторы, которые при представлении BP используют собственные векторы (СВ) ковариационных матриц (КМ) и главные компоненты (ПС) [6, 7]. Полезная модель [1] представляет собой такой анализатор, получивший название «Анализатор собственных векторов и компонент сигнала» или «Айгеноскоп» (от англ. «eigenvector» - «собственный вектор»). Конструкция айгеноскопа приведена на фиг. 1. В этой конструкции из исходного BP, подвергнутого масштабированию и аналого-цифровому преобразованию в блоке 18, в блоке 19 формируется матрица смешанных моментов (в простейшем случае КМ), затем в блоке 20 вычисляются ее СВ и собственные значения (СЗ); полученные таким образом СВ используются для представления исходного BP в блоке 21; в этом же блоке с использованием различных статистических методов анализируются признаки BP, СВ и СЗ.
Вычисление КМ в текущем времени требует накопления данных о BP. После того как данные накоплены, их необходимо превратить в оценку КМ. Получение СВ и СЗ из КМ также требует определенного времени на вычисления. Таким образом, имеется некоторая задержка с момента начала анализа BP до момента получения СВ и СЗ. В случае, если BP нестационарен, КМ меняется во времени, и ее оценку периодически необходимо обновлять.
После того, как айгеноскоп вошел в установившийся режим работы (т.е. сформирована текущая оценка КМ), на выходе блока 20 в каждый момент времени имеется набор СВ, которые соответствуют текущей оценке КМ. Именно эти СВ используются для представления BP; для этого в айгеноскопе в каждый текущий момент времени вычисляется матрица-столбец коэффициентов разложения BP в ортонормированном базисе СВ, то есть матрица-столбец:
M×M
ψ(T) - квадратная матрица размера М×М, столбцы которой - СВ для КМ в момент времени Т,
Х(Т) - сегмент BP (матрица-столбец) на интервале анализа длины М, непосредственно предшествующем текущему моменту Т,
α(Т) - матрица-столбец ГК в момент времени Т,
М - длина интервала анализа, (совпадает с размерностью КМ и числом СВ),
Т - текущее время,
(…)' - транспонирование.
Из ГК, полученных с использованием соотношения (1), могут быть сформированы матрицы  и , каждый столбец которых представляет собой столбец ГК некоторых сегментов BP.
M×n
M×m
Матрицы A(Т) и B(T) могут быть использованы в каждый момент времени Т для выявления изменений с помощью методов дискриминантного анализа (ДА).
М×n
T×m
M×n
T×m
Впервые задача о выявлении отличий между двумя множествами многомерных наблюдений была предложена и решена Фишером [5] и получила в дальнейшем [6] название задачи ДА. Для выявления отличий с использованием ДА строится так называемая дискриминантная функция (ДФ) [5]. Качество решения задачи ДА оценивается критерием дискриминации. В случае, когда для выявления изменений BP используется зависимость значения критерия дискриминации [6] от времени, дискриминантная функция (ДФ) обобщается до т.н. дискриминантного функционала или дискриминантного коллектора (ДК) - [3, 4]. ДК использует при своей работе две функционально связанные с текущим временем Т последовательности времен и и формирует матрицы наблюдений А(Т) и B(T) с использованием значений ГК для этих последовательностей времен. Значение ДК совпадает со значением критерия дискриминации, определенного на основании матриц A(Т) и B(Т) для каждого текущего момента времени Г.
A
B
У ДК различают: тип, структуру и емкость. Тип ДК определяется критерием, который используется при построении ДФ. Структура определяется взаимным расположением моментов времени в последовательностях t и t. Емкости определяются величинами m и n.
Такой подход реализуется в полезной модели [2], получившей название «Сигнализатор значимых отличий», конструкция которого приведена на фиг. 2.
М×n
T×n
Недостатком этой конструкции являются задержки, связанные с необходимостью формировать для каждого момента текущего дискретного времени Т две матрицы наблюдения А(Т) и B(T) из ГК, соответствующих предшествующим положениям сегментов BP, причем качество выявления изменений может сильно зависеть от m и n. Увеличение этих величин, с одной стороны, повышает достоверность выявления изменений, а с другой, может недопустимо уменьшить горизонт прогнозирования.
Этот недостаток устраняется в полезной модели "Сигнализатор изменений главных компонент" (СИГК) [9], в которой в качестве статистики для принятия решения без задержки предлагается использовать произведение ГК (или монотонно связанное с этим произведением значение среднего геометрического (СГ) - [7]).
M×M
M×M
На фиг. 3 представлена конструкция СИГК. Временной ряд со входа СИГК поступает на вход блока масштабирования и АЦП - 18, с выхода которого поступает на блок вычисления ковариационной матрицы - 32. В блоке 32 для некоторого интервала времени, предшествующего текущему моменту, как правило, в несколько десятков раз превышающего размер КМ М, с использованием любого известного алгоритма [7] строится оценка КМ K(T). С выхода блока 32 эта оценка поступает на вход блока вычисления собственных векторов 33. В блоке 33 решается задача вычисления СВ для K(T). Матрица СВ передается с выхода блока 33 на первый вход блока вычислителя скалярных произведений и перемножения 34. На второй вход блока 34 с выхода блока 1 масштабирования и АЦП поступает оцифрованный BP. В блоке 34 для каждого текущего момента времени формируется матрица столбец сегмента BP длительностью М, состоящая из М отсчетов, предшествующих текущему моменту времени. Далее осуществляется вычисление ГК по формуле (1). Таким образом, в каждый момент времени Т получается матрица-столбец ГК α(T), в которую входят М значений. Затем эти ГК перемножают, и полученное значение поступает на выход устройства. Величина произведения сравнивается с заданным пользователем порогом принятия решения. Превышение порога рассматривается как реперная точка (предвестник прогнозируемого события).
Прототип
В качестве прототипа выбрано устройство сигнализатора изменений главных компонент - фиг. 3. Акцентируем, что недостатком прототипа является частичная потеря свойств при выходе из строя отдельных датчиков, подключенных к его входу.
Действительно, при использовании в реальных условиях аналогов и прототипа для обработки геофизических данных возникает проблема, связанная с неидеальностью используемых многомерных временных рядов (МВР); поскольку МВР получаются с датчиков, которые имеют конечную надежность, в одномерных временных рядах (BP) (вместе составляющих МВР) появляются перерывы, что приводит к частичной потере рабочих характеристик.
Раскрытие полезной модели
Цель и преимущества: обеспечение минимальной потери рабочих характеристик и автоматическая перестройка конструкции - в условиях перерывов наблюдений, поступающих от части датчиков.
Результат достигается благодаря тому, что в предлагаемой конструкции каждый канал обработки порождает значительное число главных компонент, что позволяет сохранять рабочие характеристики конструкции даже в условиях отказа нескольких датчиков.
В прототипе изменение числа работающих датчиков требует перезапуска устройства, в то время как в предлагаемой полезной модели за счет особенностей конструкции обеспечивается бесперебойная работа при сохранении потенциально возможных рабочих характеристик.
На фиг. 4, 5, 6, 7 приведены конструкции устройств, предлагаемых в полезной модели. Эти устройства отличаются друг от друга только конструкцией своей выходной части, включенной после блока 4. Поэтому описание начнем с общей части этих конструкций.
К каждому входу УИР (количество входов УИР равно U) подключен свой блок 1 вычислителя главных компонент, конструкция которого имеет вид, приведенный на фиг. 5.
Блок 1 состоит из блока вычисления матрицы смешанных моментов 7, блока вычисления и отбора собственных векторов 8, вычислителя скалярных произведений 6 и работает как специализированный айгеноскоп с конструкцией, аналогичной [1].
Все U выходов блоков вычислителя главных компонент 1 соединены с 1-м - U-м входами блока хранения и выборки наборов значений главных компонент 2, в котором накапливаются наборы значений ГК, поступающих со всех U выходов блоков 1, для последовательных моментов дискретного времени. В каждый из дискретных моментов времени на выходе 1 блока 2 появляется набор ГК, который поступает на вход блока вычисления синдрома 3. Этот набор ГК соответствует последнему анализируемому моменту времени.
Для объяснения того, как работает блок вычисления синдрома 3, необходимо пояснить следующее. ГК в блоках 1 могут быть вычислены не всегда. Причиной невозможности вычисления ГК являются отказы датчиков, подключенных ко входам УИР, что приводит к перерывам во BP, поступающих на входы УИР. Эту ситуацию мы рассмотрим далее более подробно; сейчас достаточно отметить, что в случае, когда какой-нибудь из U блоков 1 не может вычислить ГК, он записывает вместо числового значения ГК в память блока 2 некоторый специальный код, который фиксирует (идентифицирует) ситуацию «значение ГК отсутствует». Будем сокращенно называть этот идентификатор «нет ГК». Перенумеруем по порядку все ГК, поступающие на входы блока 2. Синдром - это вектор, который имеет размерность, совпадающую с числом ГК, поступающим на вход блока 2. Координата синдрома принимает значение 1, если соответствующая номеру координаты синдрома ГК имеет числовое значение; координата синдрома принимает значение 0, если соответствующая номеру координаты синдрома ГК имеет идентификатор «нет ГК». Таким образом, в блоке 3 синдром вычисляется для каждого текущего момента времени, и его значение поступает на 2-ой вход блока вычисления сигнала коллективного поведения главных компонент 4.
С выхода 2 блока 2 на 1-ый вход блока 4 поступают наборы ГК для всех последовательных моментов времени, включая последний. Каждый из наборов ГК представляет собой вектор, размерность которого соответствует с размерностью синдромов. В блоке вычисления сигнала коллективного поведения главных компонент 4 происходит сравнение вектора ГК со значениями синдромов. Если в векторе ГК на позициях, соответствующих единичным позициям синдрома, имеются числовые значения, то по таким наборам ГК вычисляется СКП, который в дальнейшем поступает на вход блока оценки сверху для отношения правдоподобия возникновения землетрясения в текущий момент времени и принятия решения о прогнозе возникновения землетрясения 5.
Блок оценки сверху для отношения правдоподобия возникновения землетрясения в текущий момент времени и принятия решения о прогнозе возникновения землетрясения 5, представленный на фиг. 6, состоит из блока отбора последнего значения сигнала коллективного поведения главных компонент 9, блока вычисления квантиль 10 и блока визуального анализа и сравнения 11. В блоке 10 на основе выборки значений СКП, поступающей на вход блока с выхода блока 4, вычисляется шкала квантиль, которая передается на второй вход блока 11, в котором происходит оценка и сравнение последнего значения СКП, поступившего с выхода блока 9 на первый вход блока 11, с величиной квантиль и выявление аномалий СКП.
В конструкции УИР, представленной на фиг. 7 (п. 2 ФПМ), в качестве блока оценки сверху для отношения правдоподобия возникновения землетрясения в текущий момент времени и принятия решения о прогнозе возникновения землетрясения 5 используется индикатор риска сильного землетрясения, имеющий конструкцию индикатора риска по п. 1 ФПМ [15]. В этом случае на первый вход блока 5 индикатора риска сильного землетрясения подается оцифрованный ряд произошедших аномальных событий (ОРАС) - землетрясений, а второй его вход соединен с выходом блока 4. На выходе блока 5 получается оценка сверху для функции правдоподобия - отношения апостериорной (полученной в результате наблюдений) вероятности возникновения землетрясения к ее априорной (до наблюдения) вероятности и отношение правдоподобия для времени возникновения землетрясений, полученная путем свертки оценка сверху для функции правдоподобия адаптивным фильтром, настроенным по значениям ОРАС - [15]. Конструкция индикатора риска сильного землетрясения, выполненная в виде индикатора риска по п. 1 ФПМ [15], приведена на фиг. 8.
Остановимся на работе блоков 1-4 конструкции фиг. 4 и фиг. 7 и их взаимосвязях подробно. Начнем с блока 1 вычислителя главных компонент, конструкция которого представлена на фиг. 5. В блоке 7 вычисляется текущая оценка КМ
(Z)
M,M
K(T) - оценка КМ для Z-ого BP,
(Z)
M,(T-M+1)
X - траекторная матрица размера М×(Т-М+1) для Z-ого BP,
(…)' - транспонирование матрицы,
i,j
(Z)k - коэффициент ковариации для i-го и j-го отсчетов Z-ого BP на интервале анализа,
м - интервал анализа, Т - текущий момент дискретного времени Т>М.
Коэффициенты ковариации и траекторная матрица определяются соотношениями:
где
(Z)
i
S - i-ый отсчет Z-ого BP
d - номер элемента ансамбля, принимающий все целые значения от 1 до Т-М+1. В общем случае индекс d принимает значение из некоторого множества целых чисел в диапазоне от 1 до Т-М+1 (например, с фиксированным не единичным шагом).
В блоке вычисления и отбора собственных векторов 8 для каждого момента времени по оценке КМ (2) с использованием стандартных средств вычисляются СВ и СЗ, удовлетворяющие соотношению [7]
(Z)
(Z)
i
i
где ψ(T) и λ(T) - i-ый СВ (матрица-столбец) и i-ое СЗ (число), соответствующие KM Z-ого BP, определенной в момент дискретного времени Т.
Поскольку в работе УИР используются только СВ, а СЗ не используются, то при вычислении КМ по формулам (2) или (3) достаточно использовать только числитель.
При отборе СВ, используемых далее при вычислении ГК, может использоваться визуальный анализ спектра собственных значений. При таком анализе используется нормированный спектр собственных значений, на который не влияет знаменатель (2) и (3). Подробно вопрос анализа спектра СЗ рассмотрен в описании [1].
Будем обозначать числитель, представляющий собой матрицу, пропорциональную ковариационной, и ее элементы, соответственно  и .
Как следует из вышеприведенных соотношений, оценка КМ и определенные для нее СВ меняются (обновляются) для каждого момента времени. В случае, если в работе датчиков, соединенных со входами УИР, возникают перерывы, то для каждого j-го сегмента BP без перерывов оценивается парциальная величина , а результирующая оценка определяется как их сумма
где
D - число непрерывных (без перерывов) сегментов BP, длительность которых не меньше чем М, предшествующих текущему моменту дискретного времени Т,
(Z)
j
T - момент окончания j-го непрерывного сегмента Z-го BP.
(Z)
i
Собственные векторы ψ(T) определяются для оценки матрицы (6). Именно из этих СВ и отбираются L СВ, которые подаются на вход 2 блока вычислителя скалярных произведений 6. В общем случае в каждом из каналов обработки может использоваться свое число СВ.
Каждый из комплектов отобранных СВ, вычисленных для матрицы (6) конкретных BP, состоит из L СВ, которые и используются далее для анализа - каждый своего BP. ГК BP в дискретный момент времени Т определяются в блоке 6 соотношением, аналогичным (1)
где
(Z)
(Z)
L
- матрица-столбец размера L, в которой каждый j-ый элемент представляет собой скалярное произведение Х(Т) и j-го столбца Ψ(T),
- матрица-столбец, состоящая из последних М отсчетов Z-ого BP, предшествующих текущему моменту времени Т,
- матрица размера М×L, в которую в качестве столбцов записан Z-ый комплект отобранных СВ.
(Z)
Естественно, соотношение (7) используется только тогда, когда X(T) не имеет перерывов.
К числу наиболее эффективных способов формирования СКП в блоке 4 следует отнести среднее геометрическое (СГ) для модулей ГК [14]. Значения СГ для моментов дискретного времени, начиная с M-го момента после начала BP или окончания перерыва BP, вычисляется по формуле
где
Н(T) - СГ, определенное для момента дискретного времени T>М,
- число работающих входов УИР, используемых для вычисления СКП в текущий момент времени,
L - число отобранных в блоке 8 СВ.
При использовании других способов формирования СКП, H(t) вычисляется по одной из формул таблицы 2. Способы формирования СКП в блоке 4 УИР описаны в таблице 2.
Технический результат от использования конструкции, предлагаемой в прототипе, состоит в значительном сокращении задержки - за счет замены статистической обработки (например ДА) на вычисление произведения коэффициентов разложения сегмента BP по специально отобранным СВ. Вычисление такого произведения возможно практически без задержки, что было показано в описании [9].
Наиболее эффективно значимые отличия (аномалии) коллективного поведения главных компонент выявляются, если в качестве порога для их фиксации во времени использовать шкалу квантиль. Напомним, что квантиль
q
H(Т) порядка q определяется соотношением [7]
где
q
H(T) - изменяющаяся во времени квантиль,
q - порядок квантили,
T
F(…) - изменяющаяся во времени оценка функции распределения СКП Н(Т).
Сравнение последнего значения СКП, отобранного в блоке 9, с текущим значением величины квантили, определенное по формуле (9) в блоке 10 конструкции, представленной на фиг. 6, производится в блоке визуального анализа и сравнения 11, на первый вход которого поступает текущее значение квантиль, а на второй - последнее значение СКП.
Как уже отмечалось, при работе устройств, реализующих полезные модели [2] и [9] в реальных устройствах, возникает проблема перерывов в работе датчиков. Это приводит к тому, что в каждый момент времени СКП (в частном случае - среднее геометрическое) вычисляется по разному числу BP. Естественно, в каждый момент времени в качестве порога должна использоваться квантиль, определенная по всей предыстории, когда датчики, работающие в момент времени Т, также работали. Обозначим  матрицу строку, имеющую размер М, элементы которой равны единице (если датчик с соответствующим номером работает) и нулю (если датчик с соответствующим номером не работает). Матрица-строка g(T) и есть ранее упоминаемый синдром.
Синдром формируется в блоке 3. Отметим, что выборка на входе блока 4 получается путем отбора тех значений из хранимых в блоке 2, для которых наблюдается единичное значение ключа
где
А - множество номеров единичных координат синдрома, поступающего по (U+1)-му входу блока 2,
k
z(t) - в дискретный момент времени t принимает единичное значение, если k-ая компонента ГК имеет числовой код, и принимает нулевое значение, если k-ой компоненте ГК соответствует код «нет ГК».
Пусть выборка на входе блока 10 к текущему дискретному моменту времени Т записывается матрицей строкой
где
o
H(i) - значения СКП, вычисленные по отобранным в соответствии с (10) ГК,
W - число отобранных наборов ГК.
q
Упорядочим  по возрастанию и обозначим упорядоченную матрицу-строку . Тогда оценка квантили, H(T) определенная по выборке (11) будет равна
где
floor (…) - целая часть.
Вычисленное с использованием (12) значение квантили и используется в блоке 11.
Приведем более развернутое описание работы блоков УИР (фиг. 5 - фиг. 8), раскрывающих причинно-следственные связи технического результата, заключающего в автоматической перестройки конструкции - в условиях перерывов наблюдений, поступающих от части датчиков.
На фиг. 7 представлена конструкция устройства индикации изменения риска сильных землетрясений по результатам многоканального наблюдения с перерывами (к пункту 2 формулы полезной модели).
Блок 1 преобразует исходный BP в несколько BP главных компонент. Число ГК зависит от числа СВ, отобранных в блоке 10, входящем в состав блока 1. В случае, если в исходном BP появляется перерыв, то на выходе блока появляется идентификатор «нет ГК», который передается на соответствующий вход блока 2 вместо числового значения.
На вход блока 2 из блоков 1 поступают числовые значения ГК или идентификаторы «нет ГК». Все каналы поступления ГК пронумерованы. Все поступающие на вход блока 2 значения хранятся в памяти блока в виде наборов значений ГК. Каждому значению ГК соответствуют уникальные номер ГК и номер дискретного момента времени.
Набор значений ГК, соответствующий последнему анализируемому (текущему) моменту дискретного времени, поступает на вход блока 3, в котором определяется синдром (вектор, единичные и нулевые компоненты которого - в порядке нумерации каналов ГК - соответствуют числовым значениям ГК и идентификторам «нет ГК», соответственно). Синдром используется далее для определения шкалы квантили, которая отвечает набору каналов поступления числовых значений ГК в текущий момент времени.
Наборы ГК для всех последовательных моментов времени, включая последний, передаются из блока 2 для дальнейшей обработки.
В блоке вычисления синдрома 3 вычисление синдрома производится по формуле
i
где  - цикл, V - общее число ГК, F - координата синдрома. Результат передается на выход блока.
В блоке вычисления сигнала коллективного поведения главных компонент 4, для каждой отобранной матрицы столбца ГК, отвечающей соотношению (7), вычисляется значение СКП - по тем координатам, которым соответствуют единичные координаты синдрома. Так на выходе блока 4 формируется выборка значений СКП.
В блоке оценки сверху для отношения правдоподобия возникновения землетрясения в текущий момент времени и принятия решения о прогнозе возникновения землетрясения 5 индикатор риска сильного землетрясения преобразует СКП, поступающий с выхода блока 4, в величину, пропорциональную оценке сверху для отношения правдоподобия возникновения землетрясения в текущий момент времени. Осуществляет прогноз возникновения землетрясений путем сглаживания оценки сверху для отношения правдоподобия при помощи адаптивной фильтрации.
На фиг. 5 представлена конструкция устройства индикации изменения риска сильных землетрясений по результатам многоканального наблюдения с перерывами (по пункту 3 формулы полезной модели).
(Z)
В блоке вычислителя скалярных произведений 6, реализуется соотношение (7) для отобранных в блоке 8 СВ. В случае, если матрица столбец X(T) имеет перерывы (на любой позиции), то на выходе блока 6 выдается идентификатор «нет ГК».
В блоке 7 вычисляется по формуле (6) сумма компонент, представляющих собой матрицы, пропорциональные ковариационным, определенным для каждого из отрезков исходного BP, имеющих длительность, большую интервала анализа.
Блок 8 вычисляет и отбирает собственные векторы для ковариационной матрицы, поступающей на вход блока - для каждого текущего момента времени.
На фиг. 6 представлена конструкция устройства индикации изменения риска сильных землетрясений по результатам многоканального наблюдения с перерывами (по п. 5 ФПМ).
В блоке 9 происходит отбор последнего значения СКП ГК для дальнейшего сравнения с величиной квантиль.
В блоке 10 на основании выборки значений сигналов коллективного поведения главных компонент строится шкала квантиль заданного уровня.
В блоке 11 происходит визуализация последовательности текущих значений сигнала коллективного поведения в интерьере шкалы квантиль, вычисленной для каждого текущего момента времени.
На фиг. 8 представлена конструкция индикатора риска сильного землетрясения, выполненная в виде индикатора риска по п. 1 ФПМ [15].
В блоке 12 для каждого дискретного момента времени T:
T
1. По выборке СКП, предшествующей этому текущему моменту времени, строится оценка выборочной функции распределения (ВФР) - F(x).
T
2. Для текущего значения СКП - X, поступающего на второй вход ИР, по ВФР определяется вероятность превышения текущего значения по формуле .
В блоке 13 формируется величина  - оценка сверху для отношения правдоподобия возникновения землетрясения, которая подается на вход блока фильтров 14.
t
В блоке 14 формируется группа сигналов , каждый из которых является сверткой входного сигнала Y, с ИХ адаптивного фильтра.
В блоке 15 определяется подобие ОРАС - землетрясений и сигналов с выхода блока фильтров 14.
Блок 16 выдает на своем выходе сигнал с того фильтра, номер которого поступил второй вход блока 16.
В блоке 17 происходит визуализация значений оценки сверху функции правдоподобия возникновения землетрясения и принятие решений.
Доказательства технических результатов. Причинно-следственные связи
Основные технические результаты от применения заявляемых конструкций это:
устойчивость работы конструкции УИР при перерывах в работе датчиков, приводящей к перерывам во временных рядах.
На фиг. 9-12 приведены иллюстрации работы конструкции фиг. 4 и фиг. 7 в условиях, когда все датчики работают. Эти иллюстрации практически идентичны приведенным в [9]. На фиг. 9 - исходные BP, получаемые от датчиков концентрации радона, на фиг. 10 - значения ГК, на фиг. 11 - среднее геометрическое ГК, на фиг. 12 - среднее геометрическое ГК (полулогарифмический масштаб).
Для периода с 1 ноября 2013 года по май 2014 года были отобраны землетрясения на Камчатке с магнитудой Ml, превышающей значение 5.6 (сведены в таблицу 3).
Проиллюстрируем достижение заявленных технических результатов, возникающих при использовании конструкций, приведенных на фиг. 4 и 7, и причинно-следственные связи, приводящие к получению этих технических результатов.
Иллюстрацию удобнее начать с конструкции, приведенной на фиг. 7. В качестве отправной точки возьмем ситуацию, представленную на фиг. 13.
Фиг. 13а иллюстрирует разворачивание во времени вычисления СГ и его квантили порядка 0.75 и 0.95. Как видно из верхнего графика, по мере продолжения наблюдения, величины квантили изменяются, что является следствием изменения во времени оценки ФР для СГ. Из нижнего графика фиг. 13а видно, что синдромы всех датчиков единичны, что свидетельствует о том, что перерывов в их работе нет. Отметим, что на нижнем графике синдромы визуализируются следующим образом. К номеру датчика прибавляется величина синдрома, умноженная на . Такая визуализация позволяет на одном графике отслеживать поведение синдромов всех датчиков одновременно.
На фиг. 13б представлен сигнал на выходе блока 5, представленной на фиг. 7. Как видно из графика, имеются два доминирующих максимума (реперные точки), которые соответствуют 30 и 100 суткам с момента начала наблюдения (01.11.2013). Эти реперные точки соответствуют с высоким правдоподобием аномалиям СКП.
Для преобразования оценки сверху функции правдоподобия аномалий СКП в величину, пропорциональную оценке сверху функции правдоподобия землетрясения, использовался гауссовский фильтр с характеристикой, которая была оптимизирована на данных 2013 года (изображена на фиг. 16). Как видно из графика фиг. 13в после наблюдаемого роста значений оценки сверху функции правдоподобия возникновения землетрясений действительно наблюдались две серии из трех и двух землетрясений с магнитудой Ml>5.6, соответственно.
Описанную выше ситуацию будем рассматривать как базовую и смоделируем вариант отказа датчиков. На фиг. 14 и 15 представлены два варианта: фиг. 14 иллюстрирует ситуацию перерыва в работе 4 датчика, который начался на 35 сутки и закончился на 125 сутки. На фиг. 15 представлена ситуация полного отказа 4-го датчика. В этих двух вариантах в течение всего периода анализа не работал второй датчик. Описанное выше состояние датчиков иллюстрируется графиками синдромов для рассматриваемых вариантов (нижние графики фиг. 14а и 15а, соответственно).
Как видно из графиков фиг. 14б, 15б, 14в и 15в, оценки отношений правдоподобия (как для аномалий СКП, так и для момента возникновения землетрясения) позволяют уверенно прогнозировать локализацию во времени предстоящих землетрясений.
Это является следствием того, что каждый канал обработки порождает значительное число ГК (в нашем случае это ГК, которые соответствуют СВ с номерами от 2 до 27, то есть их 26 на каждый датчик); поведение ГК в случае аномалии СКП в значительной степени синхронное (см. на резкое практически одновременное увеличение всех ГК на фиг. 12); это и сохраняет работоспособность конструкции даже в условиях отказа двух датчиков из шести.
В прототипе выход из строя датчика приводил к тому, что канал изымался из обработки на неопределенное время. Восстановление датчика требовало перезапуска устройства, что приводило к выводу его из режима наблюдения на срок, требуемый на его адаптацию (от месяца до двух месяцев), что в условиях локализованных во времени аномалий СКП приводит к повышению вероятности пропуска аномалий.
В рассмотренных выше примерах легко себе представить ситуацию, когда перерыв в работе четвертого датчика возник бы на 25-е сутки, что означало бы, что мы для первой серии (из трех землетрясений) получали бы результаты представленные на фиг. 15 - до момента полной адаптации устройства после устранения неисправности датчика. Предположим, что работоспособность была бы восстановлена за 20 суток, а период адаптации занял 60 суток, таким образом, полностью работоспособность прототипа восстановилась бы на 100 сутки. Таким образом, результаты работы прототипа были бы близки к ситуации, представленной на фиг. 15. И хотя работоспособность устройства (прототипа) сохраняется, она ниже, чем в случае заявляемой конструкции.
Таким образом, технический результат для конструкции фиг. 7 состоит в сохранении эффективности за счет отсутствия ее понижения вследствие перерывов в работе, необходимых для перенастройки устройства после восстановления работоспособности вышедших из строя датчиков. Технический результат причинно обусловлен конструкцией устройства.
В случае конструкции фиг. 4 ситуация аналогична. Отличие состоит лишь в том, что в качестве аномалии СКП рассматривается факт превышения квантили (не ниже порядка 0.95), после чего землетрясения считаются вероятными в течение 50-60 суток с вероятностью 1/(1-q), где q - уровень квантили.
Описание фигур - чертежей
Фиг. 1. Айгеноскоп - анализатор собственных векторов и компонент сигнала (Полезная модель 116242 RU): 18 - блок масштабирования и АЦП, 19 - блок вычислителя матрицы смешанных моментов, 20 - блок вычисления собственных векторов и собственных значений, 21 - блок вычислителя скалярных произведений и анализатора признаков.
Фиг. 2. Сигнализатор значимых отличий (Полезная модель №133642 RU): 22 - блок измерения, предварительной обработки и формирования матриц наблюдений; 23 - блок вычисления дискриминантной функции; 24 - блок памяти, визуализации, управления и принятия решений; 25 - блок индексирования матриц сравниваемых множеств; 26 - блок измерения, хранения и выборки; 27 - блок вычисления скалярных произведений; 28 - блок формирования матриц наблюдений; 29 - блок вычислителя матрицы смешанных моментов; 20 - блок вычисления собственных векторов и собственных значений; 31 - блок селектор собственных векторов.
Фиг. 3. Сигнализатор изменений главных компонент (Полезная модель №141416): 18 - блок масштабирования и АЦП; 32 - блок вычислителя ковариационной матрицы; 33 - блок вычисления собственных векторов; 34 - блок вычислителя скалярных произведений и перемножения.
Фиг. 4. Устройство индикации изменения риска сильных землетрясений по результатам многоканального наблюдения с перерывами (к пункту 1 формулы полезной модели): блок 1 - вычислитель главных компонент; блок 2 - блок хранения и выборки наборов значений главных компонент; блок 3 - блок вычисления синдрома; блок 4 - блок вычисления сигнала коллективного поведения главных компонент; блок 5 - блок оценки сверху для отношения правдоподобия возникновения землетрясения в текущий момент времени и принятия решения о прогнозе возникновения землетрясения.
Фиг. 5. Устройство индикации изменения риска сильных землетрясений по результатам многоканального наблюдения с перерывами - конструкция блока 1 (по пункту 3 формулы полезной модели): блок 6 - вычислитель скалярных произведений; 7 - блок вычислителя матрицы смешанных моментов; 8 - блок вычисления и отбора собственных векторов.
Фиг. 6. Устройство индикации изменения риска сильных землетрясений по результатам многоканального наблюдения с перерывами (по п. 5 ФПМ): блок 9 - блок отбора последнего значения СКП ГК; 10 - блок вычисления квантиль; 11 - блок визуального анализа и сравнения.
Фиг. 7. Устройство индикации изменения риска сильных землетрясений по результатам многоканального наблюдения с перерывами (к пункту 2 формулы полезной модели): блок 1 - вычислитель главных компонент; блок 2 - блок хранения и выборки наборов значений главных компонент; блок 3 - блок вычисления синдрома; блок 4 - блок вычисления сигнала коллективного поведения главных компонент; 5 - блок оценки сверху для отношения правдоподобия возникновения землетрясения в текущий момент времени и принятия решения о прогнозе возникновения землетрясения, выполненный в виде индикатора риска сильного землетрясения по конструкции индикатора риска по п. 1 ФПМ [15].
Фиг. 8. Конструкция индикатора риска сильного землетрясения, выполненная в виде индикатора риска по п. 1 ФПМ [15]: 12 - блок оценки вероятности превышения; 13 - блок нелинейного преобразования; 14 - блок фильтров; 15 - блок оценки подобия; 16 - блок коммутации; 17 - блок визуализации и принятия решений.
Фиг. 9. Шесть исходных BP получаемых с датчиков объемной концентрации радона на геодинамическом полигоне в районе Петропавловска - Камчатского в период с 1 ноября 2013 года по 12 февраля 2014 года. По оси абсцисс время от 1 ноября 2013 года в сутках. Время дискретизации 30 минут.
Фиг. 10. Главные компоненты (ГК) исходных BP по 2-14 СВ ковариационной матрицы размера 96 (интервал анализа 2 суток при времени дискретизации 30 минут).
Фиг. 11. Среднее геометрическое модулей ГК, определенное по данным фиг. 9. Показана шкала квантиль (уровня 0.95 и 0.99).
Фиг. 12. Среднее геометрическое модулей ГК, определенное по данным фиг. 9 в полулогарифмическом масштабе.
Фиг. 13. Пример вычисления среднего геометрического и его квантиль порядка 0.75 (пунктирная линия) и 0.95 (штрихпунктирная линия) (верхний график фиг. 13а). На нижнем графике фиг. 13а показаны синдромы датчиков. Из графиков синдромов следует, что все датчики на интервале анализа функционируют. На фиг. 13б показана оценка сверху для отношения правдоподобия в зависимости от времени при обнаружении аномалий коллективного поведения. Как видно из графика, наблюдаются значительные аномалии в точках близких к 30 и 100 суткам с начала обработки. На фиг. 13в показана сглаженная верхняя оценка отношения правдоподобия (сглаживание проводилось с использованием гауссовского фильтра с задержкой 50 суток и величиной σ=10 суток, импульсная характеристика которого приведена на фиг. 16). Вертикальными линиями на фиг. 13б и 13в показаны моменты землетрясений с магнитудами, превышающими 5.6 (таблица 3). Использованы СВ с номерами 2-27.
Фиг. 14. Пример вычисления среднего геометрического и его квантиль порядка 0.75 (пунктирная линия) и 0.95 (штрихпунктирная линия) (верхний график фиг. 14а). На нижнем графике фиг. 14а показаны синдромы датчиков. Из графиков синдромов следует, что в работе 4-го датчика имеется перерыв с 35 по 125 сутки, а 2-ой датчик не работал в течение всего периода анализа; также начиная с 165 суток прекратил работу 1-ый датчик. На фиг. 14б показана оценка сверху для отношения правдоподобия в зависимости от времени при обнаружении аномалий коллективного поведения. Как видно из графика, наблюдаются значительные аномалии в точках близких к 30 и 100 суткам с начала обработки. На фиг. 14в показана сглаженная верхняя оценка отношения правдоподобия (сглаживание проводилось с использованием гауссовского фильтра с задержкой 50 суток и величиной σ=10 суток, импульсная характеристика которого приведена на фиг. 16). Вертикальными линиями на фиг. 14б и 14в показаны моменты землетрясений с магнитудами, превышающими 5.6. Сравнение фиг. 14в и фиг. 14в показывает, что перерывы привели к снижению сглаженной оценки правдоподобия с 7 до 6 (второй максимум), а также к снижению первого максимума с 35 до 7.
Фиг. 15. Пример вычисления среднего геометрического и его квантиль порядка 0.75 (пунктирная линия) и 0.95 (штрихпунктирная линия) (верхний график фиг. 15а). На нижнем графике фиг. 15а показаны синдромы датчиков. Из графиков синдромов следует, что 2-ой и 4-й датчики не функционировали на всем интервале анализа, а 1-ый прекратил работу, начиная со 165 суток. На фиг. 15б показана оценка сверху для отношения правдоподобия в зависимости от времени при обнаружении аномалий коллективного поведения. Как видно из графика, наблюдаются значительные аномалии в точках близких к 30 и 100 суткам с начала обработки. На фиг. 15в показана сглаженная верхняя оценка отношения правдоподобия (сглаживание проводилось с использованием гауссовского фильтра с задержкой 50 суток и величиной σ=10 суток, импульсная характеристика которого приведена на фиг. 16). Вертикальными линиями на фиг. 15б и 15в показаны моменты землетрясений с магнитудами, превышающими 5.6. Сравнение фиг. 14в и фиг. 15в показывает, что перерывы привели к снижению сглаженной оценки правдоподобия с 7 до 6 (второй максимум), а также к снижению первого максимума с 35 до 7.
Фиг. 16. Импульсная характеристика гауссовского фильтра, с помощью которого проводилось сглаживание оценки отношения правдоподобия. Фильтр имеет задержку 50 суток и параметр σ=10 суток.
Осуществление изобретения
Предлагаемая конструкция используется при построении комплексов краткосрочного прогнозирования и оперативной оценки рисков сильных землетрясений.
Конструкция реализуется как программно-аппаратный комплекс.
Источники информации
1. Грунская Л.В., Исакевич Д.В., Исакевич В.В. Анализатор собственных векторов и компонент сигнала. Полезная модель RU 116242. Патентообладатели ООО «БизнесСофтСервис», Исакевич Д.В., Исакевич В.В. http://bankpatentov.ru/node/207042
2. Исакевич Д.В., Исакевич В.В. Сигнализатор значимых отличий. Полезная модель №133642 RU. Патентообладатели ООО «БизнесСофтСервис», Исакевич Д.В., Исакевич В.В., Балакирев А.Н. http://bankpatentov.ru/node/414199
3. Исакевич В.В. и др. Выявление нестационарных участков при помощи нелинейной модели процесса // Радиотехника и электроника. 1995, том 40, №2, с. 255-260.
4. Грунская Л.В., Исакевич Д.В., Исакевич В.В. и др. Каскады дискриминантных функционалов в задачах анализа временных рядов в базисах собственных векторов ковариационных матриц. Нелинейный мир, No 4, 2012.
5. Fisher R.A. The Use of Multiple Measurements in Taxonomic Problems//Annals of Eugenics. - 1936 T.7. - p. 179-188.
6. Кендалл M., Стьюарт А. Статистические выводы и связи. Главная редакция физико-математической литературы изд-ва "Наука", 1973.
7. Корн Г., Корн Т. Справочник по математике для научных работников и инженеров: Определения, теоремы, формулы: Пер. с англ. М.: Наука, 1970. 720 с.
8. Исакевич В.В., Фирстов П.П., Исакевич Д.В., Грунская Л.В., Макаров Е.О. Использование методики дискриминантных функционалов для выявления реперных точек во временных рядах концентрации почвенного радона на сети пунктов Петропавловск-Камчатского геодинамического полигона / Проблемы комплексного геофизического мониторинга Дальнего Востока России. Труды Четвертой научно-технической конференции. Петропавловск-Камчатский. 29 сентября - 5 октября 2013 г. / Отв. ред. В.Н. Чебров. - Обнинск: ГС РАН, 2013. С.59-63.
9. Исакевич В.В., Исакевич Д.В., Грунская Л.В., Фирстов П.П. Сигнализатор изменений главных компонент. Полезная модель №141416. Патентообладатель ФГБОУ ВО Владимирский государственный университет.
10. Любушин А.А. (мл). Агрегированный сигнал систем низкочастотного геофизического мониторинга // Физика Земли. 1998. №1. С. 69-74.
11. Макаров Е.О., Фирстов П.П., Волошин В.Н. Аппаратурный комплекс для регистрации концентрации подпочвенных газов с целью поиска предвестниковых аномалий сильных землетрясений Южной Камчатки // Сейсмические приборы. 2012. Том 48. №1.
12. http://www.emsd.ru/ Официальный сайт Камчатского филиала геофизической службы РАН. Главная страница. Раздел «Последние 10 землетрясений Камчатки с Ml>4,5»
13. Акофф Р., Эмери Ф. О целеустремленных системах. Пер. с англ. Под ред. И.А. Ушаков. М.: Сов. радио, 1974, 272 с.
222
14. Фирстов ПЛ., Исакевич В.В., Макаров Е.О., Исакевич Д.В., Грунская Л.В. Апробация методики айгеноскопии для поиска предвестников сильных землетрясений в поле почвенного радона (Rn) на Камчатке (август 2012 - август 2013 гг.), Сейсмические приборы, 2014, Т.50, №3, стр. 63-75.
15. Исакевич Д.В., Исакевич В.В., Лукьянов В.Е. Индикатор риска. Полезная модель RU 159604. Патентообладатели ООО «Собственный вектор», Исакевич Д.В., Исакевич В.В., Балакирев А.Н., Черников СВ.