Устойчивое оценивание одновременного тренда среднего и дисперсии случайного сигнала Текст научной статьи по специальности «Математика»

Управление, вычислительная техника и информатика

УДК 519.2

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

А.В. Китаева

Томский политехнический университет E-mail: [email protected]

Рассматриваются параметрические робастные оценки одновременного тренда среднего и дисперсии случайного процесса, построенные по дискретным независимым наблюдениям, аналогичные по структуре оценкам квартилей распределений. Квартили распределения шумов предполагаются фиксированными. Показана сильная состоятельность и асимптотическая нормальность предложенных оценок

Ключевые слова:

Робастное оценивание, временные ряды, тренды среднего и дисперсии, сходимость почти наверное, асимптотическая нормальность.

Введение

Классические статистические процедуры дают хорошие результаты, как правило, в случае нормального распределения шумов. Если же распределение помех наблюдений имеет «тяжелые хвосты», например, F(x)~e—м, или представляет собой смесь распределений, к примеру, гауссовских, F(x)=(ls)N(0,al2)+sN(0,a22), s<<1, cjl2<<a22, где второе распределение вносит «аномальные» помехи, то их эффективность может быть очень низка. Устойчивые (или робастные) статистические методы сохраняют работоспособность для достаточно широкого класса распределений и менее подвержены влиянию «аномальных» ошибок наблюдений. Основы современной теории робастности были заложены в 60-х годах прошлого века [1]. В настоящее время робаст-ные методы находят широкое применение в самых разнообразных областях не только техники, где их применение стало обыденным, но и экономики, например, - макроэкономике [2], маркетинге [3, 4], финансовом анализе [5, 6]. Следует отметить также устойчивый интерес к применению робастных статистических методов (в частности, М-оценок) в машинном зрении (computer vision) [7, 8].

Метод наименьших модулей, предложенный и использованный Р.Дж. Босковичем (R.J. Boscovich) в 1757 г. является, по-видимому, исторически пер-

вой робастной процедурой, появившейся, кстати, на 50 лет раньше метода наименьших квадратов (МНК) - самого популярного, хотя и неустойчивого, метода оценивания. В экономике и финансах метод наименьших модулей является старейшей и наиболее распространенной робастной альтернативой МНК - В. Шарп (М 8Иагре) применял его для анализа рынка ценных бумаг еще в 1971 г. [9]. Процедуры оценивания, основанные на норме Ц, обладают хорошими робастными свойствами [10], но сложны для исследования и вычисления по сравнению с МНК. В данной работе предлагается использовать знаковую функцию в качестве меточной, что характерно для оптимизационных критериев по норме Ц.

Параметры масштаба случайного сигнала приходится оценивать наряду с параметрами сдвига при применении М-оценок параметров положения. Эта проблема возникает, например, при обработке сигнала в большинстве процедур машинного зрения [11]. В данной работе предполагается, что среднее и дисперсия меняются с течением времени и искомые функции допускают разложение по некоторым заданным ортогональным системам функций, так что задача сводится к оцениванию параметров этого разложения, что является обычной постановкой в параметрическом анализе временных рядов [12].

1. Постановка задачи

Будем рассматривать следующую модель наблюдений

в т _

X = X (' / N) + п< ехР[Ха Фк (' / N)], ' = 1 N

к=0 к=0

где ц - независимые одинаково распределенные случайные величины с плотностью {%(•), к=0,д} и {фк(.), к=0,т} - системы ортонормирован-ных на [0,1] полиномов порядка от нуля дотсоответственно, причем полиномы {%(•), к=0,д} орто-

т

нормированны с весом ехр-2[ХакФк(х)], аки ак -

к =0

параметры, подлежащие оцениванию. Экспоненциальное представление тренда дисперсии вызвано удобством последующих вычислений.

Оценки параметров трендов среднего и дисперсии соответственно ak, к=0,д и ак, к=0,т предлагается искать из следующей системы уравнений

х sign

n=Ew-(i 1N )exp-1

s

xi -Ea Wk(i 1N)

_ k =0

n.=E фj s-i(i / n ) х

i=1

s

xi "E a Wk (i 1 N) "

k =0

m ^

"exP[Ea^k(i 1 N)]

E ^ (i i n )

= 0, j = 0,s;

(1)

sign

+ 1I2

преобразования у(-) не вырожден в этой точке. Тогда последовательность {хп} такая, что уп(хп)=0, сходится к х0.

Доказательство. Заметим, что ||у (х0)—у„ (хи)|| —>0, поскольку у„(хв)—у(х»)=0, и у„(х„)=0.

Рассмотрим окрестность У точки х0, в которой функция у непрерывна - т. к. у(х) имеет единственный нуль в £ и ёйР(х0)^0, то существует некоторая окрестность УсУ, в которой функция у строго монотонна.

Для любых точек ^ ^еУ для Уе>035>0:||у(^1)-у(52)||<5=о||51-^2||<е'; с другой стороны, в силу сходимости последовательности функций для Уе>03N:n>N:\\yя(sl)-y(sl)||<ел|\yя(s2)-yis2)||<е• Заметим, что для Уп>И

Цх*) - У С^Ц ^ I к^) - У М! +1 \Уп С^) - У (в 2^1 +

41 Уп СО - Уп С^Ц < 2е +||Уп («1) - Уп С^Ц-

Возьмем е<5/4, тогда из условия ||уп(^1)-уп(^2)||<5/2 будет следовать Ц^-^Ке' Таким образом, для любых ^ ^еУ и для любого е0 существует 50 такое, что из условия ||уп(^1)-уп(^2)||<50 при достаточно больших номерах п следует ||51-52||<е0.

Итак, из сходимости ||у(х0)-уп(хп)||—0 будет следовать сходимость ||х0-хп||—0. Лемма доказана.

Введем векторные обозначения

(a > (~ > ас ( a ^ 0 (" > a0

a = a j = \am J , a = a \ s J , a = va s J W (•) =

= 0, у = в +1, в + т +1,

-1, х < 0;

где знаковая функция sign(х) = <! 0, х = 0; . Исполь-

1, х > 0.

зование знаковой функции в силу ее ограниченности обеспечивает устойчивость оценок к «загрязнению» наблюдений, но, с другой стороны, затрудняет исследование их свойств, вследствие ее не диф-ференцируемости в нуле.

Исследуем асимптотические свойства оценок при неограниченном возрастании числа наблюдений, пользуясь методикой, примененной в работе [13].

2. Сходимость оценок почти наверное

Прежде всего, покажем сходимость оценок к истинным значениям параметров. Докажем вспомогательный результат.

Лемма. Пусть {уп} - последовательность функций Вк—Вк, которая поточечно сходится к функции у, причем существует единственная точка х0:

у(х0)=0, и якобиан аег до = аег т(•);У = 1, к

Теорема 1. Пусть плотность распределения вероятностей p(-) симметрична относительно 0 и

i

Jp(x)dx = 1/4, причемp(0)^0,p(1)^0, иp(.) непре-0 6 6 рывна в точках х=0 и х=1. Тогда оценки а и а сильно состоятельны.

Доказательство. Рассмотрим систему (1). Найдем Gi (а, а) = lim M (п)/ N =

j N ^да j

1 m ui(x) _

= -2jtyj(х)ехр-1[Еа*Фк(x)] J p(t)dtdx, j =0,5;

Gi (a, a) = lim M (n. )I N =

1

=i<

0

где

u2 (x)

"s"1(x)[-2 J p(t)dt + 1I2]dx, j = s +1, s + m +1,

Ui( x) = XAaWk (x) exp-1

u2(x) =

EAakWk(x)+exp

k=0

m

;exp-1 EaA(x)

Еаф(x)

_k=0 _

m

Ea4k(x)

k=0 _

Aak = a k - a.

Сходимость почти наверное ц/Ж к Ц при следует непосредственно из усиленного закона больших чисел.

Покажем, что система

G (а, a) = 0, j = 0,5

(2)

Gj (а, a) = 0, j = s + 1,5 + m +1

(3)

da k

a=a ,a=a

-2p(°)¡Wk (x)Wj (x) exp-2[ Xa^k (x)]dx -0 k =0

= -2 p(0)8jk, j, k = Ö~s;

= 0, j = Ö~s, k = 0m;

dGj (a, a)

dak

dGj (a, a)

a=a ,a=a

da k

1

-2 P(1)¡ Wk(x^j - s-1(x) exp-1

1*кФк(x)

dx,

j = s +1, s + m +1, k = 0, s; dGj (a, a)

dak

= -2 P(1)<5

j-s-1k'

где sß =

0, k Ф j,

1, k = j.

имеет единственное решение а =а при любых значениях а и а.

х

Действительно, функция g(x)=\p(t)dt обращается в нуль только в точке х=0, поскольку она монотонна в силу неотрицательности плотности и строго монотонна в некоторой окрестности точки х=0 в силу непрерывности плотности в этой точке, и р(0)^0. Согласно системе (2) функция ¡(и1(х)) должна быть ортогональна с положительным весом ехр-1[|;акфк(х)] системе полиномов {%(.), к=0,д}, следовательно, она должна иметь не менее я+1-го нуля на [0,1] [14]. С другой стороны, уравнение и1(х)=0 имеет не более д корней, т. к. $Дакщ(х) является полиномом не более чем д-го порядка и ехр-1[Ёакфк(х)]>0. Отсюда следует, что для удовлетворения системы (2) необходимо потребовать ¿Дакук(х)=0, т. е. Дак=0, к=0,д.

Аналогично показывается, что система

также имеет единственное решение та=а при а =а. В этом случае и2(х)=ехр[Ё Дакфк(х)], Дак=ак-ак, и функция f(x)=-2¡p(t)dt+f/2 обращается в нуль только в точке х=1. 0

Найдем матрицу производных преобразования, задаваемого системами (2) и (3):

дО/ (а, а)

Привлекая лемму, получаем требуемый результат.

3. Асимптотическая нормальность оценок

На6йде6м асимптотическое распределение оценок а и а.

Теорема 2. В—словиях теоремы 1 случайные векторы NДа и жДа имеют при Ж^ж нормальное распределение с нулевым средним и ковариационными матрицами соответственно Т/(4р2(0)) и

3/41 + p(1)/p(0)BBT(p(1)/ p(0) -1)

где I - еди-

4 p (1)

ничная матрица соответствующей размерности, матрица

B = ,bjk = ¡Ф -s-1( x)Wk (x)exp-1 \ Yak ф (x) f^

j = s +1, s + m +1, k = 0, s

Доказательство. Рассмотрим случайные величины

N m

j (С) = N-1/2 YwWj (i / N) exp-1[ Y^k (i / N)] X

i=1 k= 0 s m

X sign[«i - XN~ll2CkWk (i/N)exp-1[X^k (x)]],

k=0 k=0

j = M,

где с=(c0,...,cs)T - некоторый вектор параметров. Заметим, что цр^Еа )=N-1/2nj(o° ,а). Обозначим

v(i / N) = N-т Су/(i / N) exp-1[ X акфк(i / N)].

k =0

Найдем

lim M(n(N'(c)) = -2 lim Nx

N^tt J N^tt

N m v / N )

xXWj(i/N)expЛТчМ/N)]N1/2 ¡ p(x)dx =

i=1 k=0 0

s 1 [ m I

= -2p(0)Xck¡Wj(x)Wk(x)exp-21 Xakф (x) Idx =

j = s +1, s + m +1, k = 0, m,

k=0 0 \k=0 = -2 p(0)Cj,

lim cov(n(N)(c),n(N )(С)) =

N^tt j k

= ¡Wj(x)Wk(x)exp-2 [X^(x) Idx = 8jk.

0 \ k=0 )

Таким образом, ) сходится в среднеква-дратическом при N^tt к случайной величине ПЛ с )=-2p(0)Cj+gp где M(?i)=0, cov(^,q)=8ip причем распределение вектора q=(g0,...,gs)T совпадает с аоиоптотическим распределением вектора

Г(° )=(п0(Л,(° ),...,пГ(о ))T.

k=0

_( Найдем предельное распределение вектора

ц{М(0). Для любого действительного вектора 7=(/0,...,/5)т сумма

/ ^ (0) =

N s т

= N-т(I/N)ехр-'[X«кФк0/N]э!^)

.=1 к=0 к = 0

по центральной предельной теореме имеет при И—ж нормальное распределение, следовательно [15], вектор д распределен нормально.

Если последовательность векторов сИ удовлетворяет П(И)(сИ)=0, то нетрудно показать, аналогично лемме, что || сИ-С||—>0 по вероятности, где с удовлетворяет П( с)=0. Таким образом, асимптотическое распределение вектора си=лИд« совпадает с распределением вектора д/(2р(0)). Аналогично рассмотрим

П )(с, й) = N-1/2 ХФ; - s Л / N) X

I=1

s ( т

^ - N -1/2 X Ск¥к 0 / N )ехр-1|Х

к=0

т

ехр N-1/2 X ^Фк (I / N)

х-^ sign

(x) -

+1/2 U

= 0, j = s +1, s + m +1,

где d =(d6,...,dm)T - некоторьш^вектор параметров, ,#/2Да)=И-у2щ(аa). Раскрывая неопределенность, имеем

lim M (n( N '(с, d)) =

= -2 p (1)\ф] -i( x)

X cWk(x) exp-1 |X аф (x) |+

k=0

m

+X акфк(x)

dx,

как и вектор

-ШАа, т. е.

2 p(1)

имеет гаус-

совское распределение с нулевым средним. Найдем ковариационную матрицу этого распределения, учитывая, что

М [э1яп(п,. ^п( п,)] = ,

М[ф§п(п1 -1) + 1/2)^п(п; -1) +1/2)] = 3/48,

М [^п( п1 -1) +1/2)э1мп( п;)] =1/28;.

Получаем

lim M

N

n(N)(0,0) Bn(N)(0)

2 p(1)

2 p(0)

Ншсоу(п(" )(с, й),ц(" )(с, й)) =

N—ж J к

1

= |Ф;-s-1(х)фк-s-1(ХУ^Х = 8;к ■

0

Аналогично предыдущему, получаем, что в асимптотике вектор лиДа распределен также,

(0,0) - П (0) 2 р(1) 2р(0)

= 3/41 + р(1)/ р(0)ВВт (р(1)/ р(0) -1)

= 4^ '

что и требовалось показать. Теорема 2 доказана.

Замечание. Требование симметричности распределения шумов, вообще говоря, не является обязательным, необходимо только фиксировать квартили распределения &=-1, б2=0, 03=1.

Выбор весовых функций перед 81§п(.) в системе (1) вызван стремлением повысить асимптотическую точность оценивания. Покажем это для оценок параметров тренда среднего. Рассмотрим следующую систему, определяющую оценки а(Ь)=(а(Ь), к=0,я): N s Л(4)

П = X Ь; (. / ^^[х,.-X а к Щ (. / N)] = 0,

1=1 к= 0

;=м,

причем весовые функции Ь(.) таковы, что решение системы

1 и (х)

I bj (х) | р(/)сНсХ = 0, ; = 0, s,

0 0

s т

и1(х)=XKЩk (х)ехр-1[X«kфk(х)]

к=0 к=0

относительно ащ единственно в условиях теоремы 1 - это обеспечивает сходимость почти наверное рассматриваемых оценок к а. Ковариационная матрица асимптотического распределения оценок а щ

будет иметь вид А РА—, где матрицы 4 р 2(0)

A =

a]k = j b] (x)]Vk (x) exP-11 X ак Фк (x)

], к = 0,5

D = {djk =jbj (x)bk (x)dx, j, k = 0, s}.

0

Решая задачу минимизации Sp(A-1DA-1T) по вектору b (x)=(b0(x),...,bs(x))T, получаем

_ ___I

b(x) =^(x)exp-11 Хакфк(x)

но, т. к. истинные значения параметров тренда дисперсии а неизвестны, а оценки а сильно состоятельны, то следует положить

b(x) =^(x)exp 11 XакФк(x)

V к =0

что и сделано в системе (1).

Разрывность меточной функции в нуле создает проблемы при применении стандартных методов решения системы (1). К настоящему времени разработано много специальных алгоритмов для численного оценивания по норме Ь1 [16]. В системе (1) все оцениваемые параметры присутствуют во всех уравнениях системы, что создает дополнительные трудности в ее решении. Можно упростить систему и оценивать параметры тренда среднего независимо от параметров тренда дисперсии, если исключить из весовых^функций первых (я+1)-го уравнения а, положив Ь(х)=у(х).

Заметим, что тогда оценки а(Ь) удовлетворяют следующему критерию наименьших модулей:

I

xi -Ia t Wk (i1N )

■ min.

- (b) a k

В этом случае будут потери в эффективности оценивания 6. Как показывают результаты расчетов, сделанные для линейных трендов среднего и дисперсии при гауссовском распределении шумов, т. е. у/0(х)=ф0(х)=1, у/1(х)=ф(х)=^2(х-1/2), р(.)=Л(0,о2), эти потери могут достигать 30 % в диапазоне значения а/а от -2,5 до 2,5.

СПИСОК ЛИТЕРАТУРЫ

1. Huber P.J. Robustness: Where are we now? // Student. - 1995. -V. 1. - № 2. - P. 75-86.

2. Atkinson A.C., Koopman S.J., Shephard N. Detecting shocks: outliers and breaks in time series // Journal of Econometrics. - 1997. -V. 80. - P. 387-422.

3. Franses P.H., Kloek T., Lucas A. Outlier robust analysis of longrun marketing effects for weekly scanning data // Journal of Econometrics. - 1999. - V. 89. - P. 293-315.

4. Seung-Hoon Yoo. A robust estimation of hedonic price models: least absolute deviations estimation // Applied Economics Letters. -2001. - V. 8. - P. 55-58.

5. Franses P.H., Van Dijk D., Lucas A. Short patches of outliers, ARCH and volatility modeling // Applied Financial Economics. -2004. - V. 14. - P. 221-231.

6. Martin R.D., Simin T. Outlier-Resistant Estimates of Beta // Financial Analysts Journal. - 2003. - V. 59. - № 5 - P. 56-69.

7. Stewart C.V. Robust Parameter Estimation in Computer Vision // SIAM Review. - 1999. - V. 41. - № 3. - P. 513-537.

8. Robust Statistical Techniques in Image Understanding. Special Issue of Comput. Vision Image Understand. - 2000. - V. 78. - 156 p.

9. Sharpe W.F. Mean-Absolute-Deviation Characteristic Lines for Securities and Portfolios // Management Science. - 1971. - V. 18. -№ 2. - P. B1-B13.

10. Shevlyakov G.L., Vilchevski N.O. Robustness in Data Analysis: Criteria and Methods. (Modern Probability and Statistics Series). - Boston: Vsp International Science Publishers, 2002. - 310 p.

11. Dahyot R., Wilson S. Robust Scale Estimation for the Generalized Gaussian Probability Density Function // Metodoloski zvezki. -2006. - V. 3. - № 1. - P. 21-37.

12. Андерсон Т. Статистический анализ временных рядов. - М.: Мир, 1976. - 758 с.

13. Китаева А.В. Медианные оценки параметров квадратичного тренда временного ряда // Автометрия. - 1990. - № 1. - С. 87-89.

14. Сеге Г. Ортогональные многочлены. - М.: Физматгиз, 1962. -500 с.

15. Андерсен Т. Введение в многомерный статистический анализ. - М.: Физматгиз, 1963. - 500 с.

16. Koenker R. L1 computation: An interior monologue. In: L1-Stati-stical Procedures and Related Topics. (IMS Lecture Notes Monograph Series. V. 31). - Hayward, California, 1997. - P. 15-32.

Поступила 10.10.2008 г.