CC BY
26
Вычислительные технологии
Том 6, № 3, 2001
ПОГРЕШНОСТИ ОЦЕНИВАНИЯ ПАРАМЕТРОВ
ПРИЛИВНОГО ОТКЛИКА В ИМПУЛЬСНОМ ЭЛЕКТРОМАГНИТНОМ ИЗЛУЧЕНИИ ЗЕМЛИ
Г. М. ВодинчАР Камчатский государственный педагогический университет Петропавловск-Камчатский, Россия
The algorithm of the interval assessment of the harmonic component parameters with known frequencies in a noise-affected signal is described. The theoretical justification is based on the fact, that the noise component of a signal is of Poisson character. The result of the implementation of algorithm to assessing the tide's response phase in the intensiveness of electromagnetic pulsed terrestrial emission is described. It has been concluded that this method is efficient on assessing errors of pointwise estimates of the tide response phase.
В работе [1] выполнено выделение гармонических компонент, модулированных приливными волнами, из сигнала, характеризующего импульсное электромагнитное излучение Земли (ЭМИ), которое оценивается интенсивностью — числом превышающих шумовой порог импульсов в минуту. Проведено точечное оценивание фаз гармонических компонент, вопрос о погрешностях полученных значений остается открытым.
В настоящей работе предложен алгоритм интервального оценивания, позволяющий определить погрешности, допускаемые при точечном оценивании. При оценке применяется прямоугольное временное окно. В [1] для выделения различных приливных волн используются временные окна разной длительности. Выделяются только одна волна и ее высшие гармоники, остальные рассматриваются как шум. Алгоритм, предлагаемый в настоящей работе, одновременно учитывает разные приливные волны.
Предполагается, что заполняющая окно вырезка x(n) является суммой детерминиро-
p
ванного полигармонического сигнала A cos uin + Bi sin uin и случайной последователь-
i= 1
ности s(n). Полигармоническая составляющая образована гармониками приливного отклика. Поэтому частотами являются частоты гармоник из спектра приливного потенциала [2], выделяющихся на шумовом фоне. В интенсивности ЭМИ такими являются группы суточных и полусуточных волн [1].
Превышение электромагнитным излучением Земли шумового порога образует пуассо-новский поток событий (пуассоновский процесс), значит, число регистрируемых в минуту импульсов распределено по закону Пуассона [3]. Типичными средними значениями являются от 1 до 40 импульсов в минуту. Малое среднее число импульсов не позволяет заменить пуассоновский закон его гауссовским приближением. После выделения гармонических составляющих среднее значение шума в пределах окна постоянно. Следовательно, элементы
© Г. М. Водинчар, 2001.
последовательности s(n) являются независимыми случайными величинами, распределенными по закону Пуассона с одинаковым средним значением Л. Требуется по наблюдаемым значениям x(n) и известным ui построить доверительные интервалы для параметров Ai и Bi. Подобная задача рассматривается в [4] в предположении нормального распределения шума. Последнее обстоятельство не позволяет непосредственно применять результаты [4] для анализа данных типа полученных в [1].
Обоснование методики. Сначала аналогично тому, как это сделано в [4], получим точечные оценки параметров Ai и Bi методом наименьших квадратов, и на основе этих оценок построим доверительные интервалы для параметров.
Положим, что целочисленное время n меняется от — N до N. Обозначим ai(n) =
cosUin— 0S1An(N + в'5^ и Pi(n) = sinUin. Тогда x(n) = Ao + ^ (A^(n) + В^(п)) + s(n), (2N + 1) sinü'5ui i=i
1 A . sin(N + 0'5)ui
где A0 = - y Ai-. Непосредственным вычислением можно установить,
2N +1^ sin0'5ui
i=1
N N N N
что Y1 ai(n) = Y1 ei(n) = Y1 ai(n)^i(n) = 0. Суммы вида Y1 ^(n)^(n) можно рас-
n=-N n=-N n=-N n=-N
сматривать как скалярные произведения (2N + 1)-мерных векторов, поэтому будем обозначать их (ai,ej'). Оптимальные в смысле метода наименьших квадратов оценки A* и В* параметров Ai и Bi соответственно являются решениями двух независимых линейных систем:
p
(ai,aj)A* = (ai,x), i = 1, 2, .„, p, (1)
ijj
j=i
p
^(А,вз)В* = (А,х), г = 1, 2, ..., р. (2)
3=1
Введем в рассмотрение р х р матрицы БА = ||(а^а3-)|| и БВ = ||(в^вз)||. Тогда А* =
р , р Л Л
^2 ¿А (а3-, х) и В* = ^ (¿В(вз, х), где через (¿А и (¿В обозначены соответственно элемен-3=1 3=1
ты матриц, обратных матрицам ДА и ОВ. Для получения интервальных оценок А и
^ рассмотрим случайные величины — А* и Bi — В*. Умножив почленно равенство р
х(п) = А0 + (А3-а3-(п) + В3-вз(п)) + з(п) на аДп) и просуммировав по п от — Ж до N,
3=1
р
получим ^, х) = ^2 (ат,,аз)А3- + (а^з). С учетом выражения для (а^х) из системы (1) 3=1
величины А — А* являются решениями системы
p
)(Aj — A*) = — (ai,s), i = 1, -, p.
j=i
р ^ р
Тогда Ai — А* = — ^2 ¿А (аз, з). Аналогично можно установить, что Bi — В* = — ^ (¿В(вз, з).
3=1 з=1
Рассмотрим последовательность аi (п)з(п) и покажем, что она удовлетворяет условиям
Ляпунова [5]. Предварительно заметим, что, специализируя формулу для центрального
к-2
момента к-го порядка распределения Пуассона ^к = А ^2 С-фг (см., например [3], с. 114)
1=0
ПОГРЕШНОСТИ ОЦЕНИВАНИЯ ПАРАМЕТРОВ ПРИЛИВНОГО ОТКЛИКА... 5
на случай к = 3, получим, что = Л. Тогда непосредственным вычислением можно уста-
N М ( ()() М ( ( \ ( \х3 \( 81п(3^ +1-5)^ 3в1п(# + 0.5)^
новить, что м (аДп)з(п) — м (аДп)з(п))) = Л
n=-N V 4sin1.5^j 4sinü.5wi
3sin(N + ü.5)w¿ sin(2N + 1)w¿ 2sin3(N + ü.5)w¿ . ^
+ --———3- . Видно, что каждое слагаемое во
2(2N + 1) sin 0.5w¿ sin (2N + 1)2 sin3 0.5w.
втором множителе при N ^ то либо ограничено, либо бесконечно мало. Тогда в целом все выражение ограничено. Величина
^ , , , „ Л Лт sin(2N + 1W 2sin2(N + 0.5W\
Bn = > D(ai(n)s(n)) = A N + 0.5 + —--— — --Ц-^ то
N ^ ^ ^ V 2 sin иг (2N + 1) sin2 0.5wJ
n=-N
N
^ 3 3/2
при N ^ то. Значит, lim { Е M (а(n)s(n) — M(аДп)з(п))) }/B^/ = 0, и условия Ля-
n=-N
пунова выполняются. Выполнение условий Ляпунова для последовательности ßj(n)s(n)
N3
непосредственно следует из соотношения Е M (ßi(n)s(n) — M(ßi(n)s(n))) =0.
n=-N
Таким образом, в соответствии с центральной предельной теоремой случайные ве-
NN
личины У] аДп)з(п) и Е ßi(n)s(n) распределены асимптотически нормально. При
n=-N n=-N
определении параметров Aj и B^ в сигнале используются временные окна длительностью 40 000 отсчетов и более, поэтому можно с большой степенью точности считать величины
NN
У] аДп)з(п) и Е ßi(n)s(n) распределенными нормально. Тогда и величины Aj — A* и
n=-N n=-N
B^ — B*, являющиеся их линейными комбинациями, тоже распределены по нормальному закону. Математическое ожидание и дисперсия величины Aj — A* определяются выра-
Р N N / p „ \
жениями M(Ai — A*) = — Ё dAА Е aj(n) = 0 и D(Aj — A*) = А Е Е ^а(n) ,
j= 1 n=-N n=-N У j=1 У
в чем можно убедиться непосредственным вычислением. Аналогично M(Bj — B*) = 0 и
/ \ 2
n i p
D(Bi — В*) = A Е ( Ё dBjßj(n) ) . Таким образом, Аг — A* и Вг — В* — центрированные
n=-N \j=1
нормальные величины со стандартными отклонениями a a
\
N I Р . \2
A Е (Е^ aj (n)| и
n=-N \j=1
ав
\
A Д (b ¿Bfí ( Л Об A* 1 ^ sin(N + 0.5Уг A Е Е d^BBßj(n) . обозначим через Ao величину т^тгт У^ A* —-
n=-N \j=1 / 2N + 1 -г=1 sin 0.5^г
0 2 N + 1
г=1
N
N А
£ x(n) — A0
n=-N У
Можно проверить, чтоМ ( ^^-1 У у х(п) — Ад ] = Л, следовательно, несмещенной оцен-
1 М *
кой параметра Л будет ^^-1 ^ ] х(п) — А*. Для построения доверительных интервалов
+ и=-М
выберем доверительную вероятность 7. Тогда с вероятностью 7 параметры А^ и В^ лежат соответственно в интервалах [А* — а а; А* + и [В* — ; В* + ], где опре-
деляется с помощью функции Лапласа Ф(^) из условия 2Ф(г7) = 7 [3]. Большое число используемых отсчетов позволяет использовать функцию Лапласа, несмотря на то, что
среднеквадратические отклонения точно неизвестны и вычисляются с помощью оцениваемого параметра Л.
Апробирование методики. Апробирование предложенной методики проводилось на реальном ряде данных интенсивности ЭМИ за три года, предоставленном В. К. Павлю-ковым. Определялись разности между теоретически рассчитанными фазами приливных волн и фазами гармоник, выделенными из сигнала. Программа для вычисления теоретических фаз волн приливных потенциалов была составлена и предоставлена А. Н. Кролевцом.
Доверительные интервалы для выделенных значений начальных фаз ^ строились следующим образом: составлялись четыре возможные комбинации границ интервалов для Ai и В^ затем из получаемых для каждой комбинации значений начальных фаз выбирались наименьшее и наибольшее в качестве границ доверительного интервала для ^. При выделении гармоник в сигнале использовалось временное окно длительностью 40 279 отсчетов. Определение фаз проводилось с шагом 180 отсчетов (3 часа). На рисунке представлена диаграмма вычисленных разностей фаз для приливной волны 0\ [2] в период с 14.03.2000 по 28.03.2000. Доверительные интервалы считались при доверительной вероятности 7 = 90 %. Величины построенных доверительных интервалов дают погрешность точечного оценивания разностей фаз.
Выводы
1. Разработан алгоритм методики интервального оценивания параметров гармонических компонент на фоне пуассоновского шума и проведено его теоретическое обоснование.
2. Апробирование алгоритма на реальных данных показало, что его можно применять для нахождения погрешностей точечного оценивания фазы приливного отклика в ЭМИ.
Автор выражает признательность А. Н. Кролевцу за постановку задачи.
Список литературы
[1] Кролевец А. Н., Плвлюков В. К. Приливной отклик импульсного электромагнитного излучения и краткосрочный прогноз сильных землетрясений // Проблемы сейсмичности Дальнего Востока / Под ред. А. В. Викулина. Петропавловск-Камчатский, 2000. С. 171-181.
[2] Мельхиор П. Земные приливы. М., 1968. 244 с.
[3] Справочник по теории вероятностей и математической статистике / В. С. Королюк и др. М.: Наука, 1985.
[4] Серебренников М. Г., Первозванский А. А. Выявление скрытых периодичностей. М.: Наука, 1965.
[5] Ширяев А.Н. Вероятность. М., 1989.
Поступила в редакцию 19 июля 2000 г., в переработанном виде —19 февраля 2001 г.
