УДК 519.22, 517.521, 53.088
А. А. Ерошенко1, А. А. Кудрявцев2, О. В. Шестаков3
ПРЕДЕЛЬНОЕ РАСПРЕДЕЛЕНИЕ ОЦЕНКИ РИСКА МЕТОДА ВЕЙГЛЕТ-ВЕЙВЛЕТ-РАЗЛОЖЕНИЯ СИГНАЛА В МОДЕЛИ С КОРРЕЛИРОВАННЫМ ШУМОМ*
В работе рассматривается задача оценки функции сигнала после применения линейного однородного оператора в модели с коррелированным шумом. Исследуются асимптотические свойства оценки риска при пороговой обработке коэффициентов вейглет-вейвлет-разложения сигнала. Приводятся условия, при которых имеет место асимптотическая нормальность оценки риска.
Ключевые слова: вейвлеты, пороговая обработка, несмещенная оценка риска, коррелированный шум, асимптотическая нормальность.
1. Введение. Во многих прикладных задачах анализа и обработки сигналов и изображений данные измеряются не напрямую, а после применения некоторого линейного оператора. Кроме того, в измерениях всегда присутствует шум, обусловленный несовершенством оборудования и различными случайными помехами.
В последние десятилетия значительно возросла популярность нелинейных методов обработки сигналов и изображений с помощью аппарата вейвлет-анализа. В задачах удаления шума обычно используется процедура пороговой обработки вейвлет-коэффициентов, которая обнуляет коэффициенты, не превышающие заданного порога. В получаемой таким образом оценке сигнала или изображения неизбежно содержатся погрешности. Асимптотические свойства оценки этих погрешностей (риска) исследовались при различных условиях измерения сигнала в работах [1-14]. В частности, в работе [14] доказывается состоятельность и асимптотическая нормальность оценки риска при использовании так называемого вейвлет-вейглет-разложения сигнала (Wavelet-Vaguelette Décomposition) в модели с коррелированным шумом. В данной работе мы предполагаем, что линейный оператор К, который применяется к сигналу, является однородным с показателем /3, т. е.
K[f(a(x - ж0))] = а~Р(К/)[а(х - ж0)]
для любого xq и любого а > 0, и рассматриваем альтернативный метод представления сигнала, предложенный в работе [15] и получивший название вейглет-вейвлет-разложение (Vaguelette-Wavelet Décomposition). Примерами однородных линейных операторов служат оператор интегрирования, преобразование Абеля и некоторые виды операторов свертки [15].
Оказывается, что при выполнении определенных условий гладкости наблюдаемого сигнала и анализирующего вейвлета оценка риска в методе вейглет-вейвлет-разложения также является асимптотически нормальной. Это обстоятельство позволяет строить асимптотические доверительные интервалы для оценки риска, не используя значений истинных ненаблюдаемых отсчетов сигнала.
2. Модель данных и коэффициенты разложения. Вейвлет-разложение функции К/ € Ь2(Ш), описывающей линейное преобразование сигнала, представляет собой ряд
Kf=Yl (кМк)Ф,к, (1)
j,kez
где = 2^2ф(2Н — к), a 'ip(t) — некоторая материнская вейвлет-функция (семейство {^jkjjkez
образует ортонормированный базис в Ь2(Ш)). Индекс j в (1) называется масштабом, а индекс к — сдвигом. Функция ф должна удовлетворять определенным требованиям, однако ее можно выбрать таким
1 Факультет ВМК МГУ, асп., e-mail: aeroshikQgmail.com
2 Факультет ВМК МГУ, доц., к.ф.-м.н., e-mail: nubigenaQhotmail.com
3 Факультет ВМК МГУ, доц., д.ф.-м.н.; ИПИ РАН, ст. науч. сотр., e-mail: oshestakovQcs.msu.su
* Работа выполнена при финансовой поддержке РНФ (проект № 14-11-00364).
образом, чтобы она обладала некоторыми полезными свойствами, например была дифференцируемой нужное число раз и имела заданное число М нулевых моментов [16], т. е.
ею
J гкф(г)йг = о, к = о,...,м- 1.
В дальнейшем будут рассматриваться функции линейно преобразованного сигнала К/ € Ь2(Ж) на конечном отрезке [а, Ь], равномерно регулярные по Липшицу с некоторым параметром 7 > 0. Для таких функций известно [17], что если вейвлет-функция М раз непрерывно дифференцируема (М ^ 7), имеет М нулевых моментов и достаточно быстро убывает на бесконечности, т. е. существует такая константа С а > 0, что
сю
г (1 + \^)\ф(1)\ ль^Са,
то найдется такая константа А > 0, что
2*^172)- (2)
На практике функции линейно преобразованного сигнала всегда заданы в дискретных отсчетах на конечном отрезке. Не ограничивая общности, будем считать, что это отрезок [0,1] и функция К/ задана в точках г/2-7, г = 1,... ,2,/: (К/)^ = К/ (г/2-7). Дискретное вейвлет-преобразование представляет собой умножение вектора значений функции К/ (обозначим его через К/) на ортогональную
матрицу Ш, определяемую вейвлет-функцией ф: К/ = ШК/ [17]. При этом дискретные вейвлет-коэффициенты связаны с непрерывными следующим образом: (К/)У| и 2•,!2{К ¡^ф^к). Это приближение тем точнее, чем больше </. Здесь не будут обсуждаться методы борьбы с краевыми эффектами, связанными с использованием вейвлет-разложения на конечном отрезке. Познакомиться с этими методами можно, например, в [18].
В реальных наблюдениях всегда присутствует шум. Пусть {е^, ! е 2} — стационарный гауссов-ский процесс с ковариационной последовательностью = соу(вг, е^). Будем полагать, что е* имеют нулевое среднее и единичную дисперсию. Рассмотрим следующую модель данных:
¥г = (К/)г + вг, г^,...^. (3)
Для £ € [0,1] определим наблюдаемый процесс
! [з-7*] г Р-7*]
ЪЮ = з7 £ = + £ е"
г=1 г=1
[2-4]
где FJ(t) = 1/2-7 #/(»/2Г) — "суммарный сигнал".
г=1
Предположим, что автоковариационная функция шума убывает медленно согласно модели Гк ~ Ак~а, где 0 < а < 1. Положим г2 = 2А/(1 — а)(2 — а) (без ограничения общности далее будем полагать, что г=1)иЯ = 1- а/2 € (1/2,1).
Определим дробное броуновское движение В #(£) как гауссовский процесс на Ж с нулевым средним и ковариационной функцией
= ^(\з\2Н + щ2Н з\2Н), ж,
где
тг пт ^Г(2 — 2Н) соз(ттН) Н = °(Вя(1)) =-к II {211 - 1)-'
Из леммы 5.1 в [19] следует, что
2»^(Г,^) ^>тВя(г), ¿6 [0,1].
Таким образом, полагая е = г1/"2 •7</2, можно аппроксимировать наблюдаемый процесс с помо-
щью У (г) для г € [0,1]:
у(г) = Р(г) + £авн(г), (4)
г
где ¿^(г) = / К/(з) йз. о
Применяя к (4) вейвлет-разложение и аппроксимируя его дискретным вейвлет-преобразованием, приходим к следующей модели дискретных вейвлет-коэффициентов [20]:
лмм Л]к — Нк + 1 2 Цк,
где /х^тг = (К/)У| и 2,112{К и = 2Л ^ ) J 'ф^к(Шн- Шумовые переменные х^к имеют стан-
дартное нормальное распределение, но не являются независимыми.
Заметим, что если преобразование К однородно с показателем /3, то К~1 тоже однородно с показателем —/3. Пусть Хук = \\К~1ф^\\. Легко показать, что А^ =
Поскольку оператор К линеен и однороден, существуют такие функции и-^к,,, что {/, «¿/г) = = {К/, ф:-]к)- При этом функция / представляется в виде ряда
/ = Е ^к(К/,'фук)иук,
з,кег
где и:-]к = К~1ф:]к/\:1к (функции и^ь называются вейглетами). Последовательность {ujk} не образует ортонормированную систему, однако если выполнены некоторые условия гладкости на К*ф и К~1'ф [21], то последовательность {ujk} образует устойчивый базис.
3. Пороговая обработка и оценка риска. Смысл пороговой обработки вейвлет-коэффициентов заключается в удалении достаточно маленьких коэффициентов, которые считаются шумом. Будем использовать так называемую мягкую пороговую обработку с порогом Т^-, зависящим от уровня
К каждому вейвлет-коэффициенту применяется функция рт}(%) = sgn(ж) (|ж| т.е. при такой
пороговой обработке коэффициенты, которые по модулю меньше порога обнуляются, а абсолютные величины остальных коэффициентов уменьшаются на величину порога.
Погрешность (или риск) мягкой пороговой обработки определяется следующим образом:
,7-1 2^-1
для = Е Е - рт, ■ (6)
3 = 0 к=О
В работах [2] и [3] было предложено использовать порог Т) = а¿л/2Ы2-1. Было показано, что при таком пороге риск близок к минимальному [2]. Этот порог получил название "универсальный". В дальнейшем будет использоваться именно такой вид порога. В выражении (6) присутствуют неизвестные величины рь^к, поэтому вычислить значение нельзя. Однако его можно оценить. В качестве оценки риска используется следующая величина [1]:
,7-1 2^-1
ыг) = Е Е (?)
3 = 0 к=О
где .Р[ж,Т,сг] = (х — о2)1(\х\ < Т2) + (а2 +Т2)1(|ж| > Т2). Можно показать [17], что Е,т(/) является несмещенной оценкой для _??,/(/).
4. Вспомогательные результаты. В этом разделе будут получены некоторые результаты, касающиеся характера зависимости эмпирических вейвлет-коэффициентов. Всюду далее предполагается, что используются вейвлеты Мейера [17], обладающие нужным количеством М непрерывных производных. Для вейвлет-функции Мейера ф(х) при любом натуральном М0 существует константа СМо > 0, такая, что \ф(0\ < СМо |С|М° ЦевиРР$у
Рассмотрим ковариацию случайных величин в модели (5) [20]:
cov
cov
cov
r(Xjk,Xit) = 2J^E I ф:)к dB и I фа dBH = 2J(-1~a)-L J ФМОФиШГ^
Для преобразования Фурье вейвлет-функции справедливо Рассмотрим ковариацию в пределах одного уровня j:
Поскольку вейвлет-функция ф имеет М непрерывных производных, найдется такая константа См > 0, что
\cov(XjQ,Xjk)\ sC--. (8)
Теперь рассмотрим ковариацию на разных уровнях сох(Х^,Хц), предполагая, что j > i:
= 2J(l-a)J_2-±-i(l-a) f ^
27Г J
где А = j — i. Повторяя рассуждения из работы [13], можно получить следующую оценку:
(Сме2J(1-"b»-i(i-«)2-AMo 1 — > 1,
со v(Xik,Xü)^l . {к2 0 9
где Сме — некоторая положительная константа.
Рассмотрим теперь структуру дисперсии оценки риска. Введем обозначение для последовательностей: a,j ~ bj, если lim —^ = 1 при J оо.
bj
Лемма 1. Пусть а >
1/2 и 7 > (4а - 2)"1, тогда Dj = DRj(f) ~ ,
где константа С
не зависит от функции сигнала f.
Доказательство. При выполнении условий леммы (27 + 1)-1 < 1^(2а)-1. Выберем р", такое, что (27 +1)-1 < р" < 1 — (2а)-1 и р" J — целое число. Тогда в силу (2) ßjk —> 0 для всех j: р" J ^ j < J при J ^ оо. Разобьем выражение (7) на две суммы:
p"J2i-l J-1 2j-l
Rj(f) = Е E tfkFlXjk<,Tj<-a:j] + £ ^ikFl^k^j^j]-
j=0 k=0 j=p"j+1 fc=o
Так как существует такая константа С'р > 0, что F[X?k, Tj, aj] ^ CpTf = CpJдля первой суммы имеем
Е Е < < (10)
j=0 k=0 j=0 k=0
с некоторой константой C'F > 0. Далее
J-12j-l J-l 2*-1 J-l 2j-l
DRj(f) = Е Е b%DF[Xik>Tj,°j] + Е Е Е Е Ti,(Ji}2, F[Xjk, Tj,(jj}2), (И)
j=0 k=0 i=0 1=0 j=О k=о
где во второй сумме предполагается, что (i,l) ф (j,k). Рассмотрим сумму дисперсий:
J-12i-l р" J 2j-1 J-l 2j-l
E E ^jkDFlXjk'.T3'-a:il = E E ^jkDFlXjk'.T3'-a:il + E E
j=0 fc=0 i=0 fc=0 j=p"J+1 fc=o
В силу (10) первая сумма не превосходит Ср,]222(1_а+(а!+2'3)р")-7; где Ср — некоторая положительная константа, и так как р" < 1 — (2а)-1, имеем 2(1 — а + (а + 2/3)р") <4/3 + 1. Учитывая вид порога Т) и принимая во внимание, что а > 1/2, для второй суммы имеем
7-1 2^—1 7-1 2^—1 7-1 2^—1
Е Е ~ Е Е = 2 Е Х%2°2з(°2з + -
3=р"7+\ к=0 з=р"7+1 к=0 з=р"7+1 к=0
/-1 2'-1 /-1 92(1 —а)/ -7"1
- Е Е Х%2а* = Е Е = 22(1"а)/+1 Е (2«Аоо)42-+2^ С"2<^>',
(12)
где С — положительная константа.
Рассмотрим теперь сумму ковариаций в (11). Аналогично сумме дисперсий имеем
7-1 2i — 1 7-1 2^—1
Е Е Е Е а^Л^соу^^2,,^,^],^^,^,^]) ~
i=Q 1=0 3 = 0 к=0
7-1 2* — 1 7-1 2'-1
* Е Е Е
¿=р"7+1 г=0 З=р"7+1 к=о Известно, что если вектор (X, У ) имеет двумерное нормальное распределение, то
СОУ(Х2, У2) = 4ЕХЕГсоу(Х, ¥) + 2СОУ'(.V. Г). (13)
Используя (8), (9), (13) и действуя, как в работе [13], получаем
Е Е Е Е1^2зкСОУ2(Хи,Х]к)^С"2^\ (14)
1=р"7+1 1=0 З=р"7+1 к=о
где (7" — некоторая константа. Аналогично с учетом (2)
7-1 21-1 7-1 2^—1
Е Е Е Е ^utfkfi'iiHkCOviXinXjk)
г=р"7+1 1=0 j=p"7+1 к=0
sC C*2(-1-a(-1-p"^4^J, (15)
где С* — некоторая константа.
Объединяя (12), (14) и (15), получаем, что DRj(f) ~ . Лемма доказана.
Сформулируем еще одно свойство эмпирических вейвлет-коэффициентов. Говорят, что последовательность случайных величин {У^}?^ обладает свойством />перемешивания, если для функции
р(т) = sup |corr(Yj, Yj)\
i,j:\i—j\>m
справедливо p(m) ^ 0 при m ^ oo.
Лемма 2. Последовательность ^X2kF[X2k, Tj, j = 0,..., J — 1, к = обладает
свойством р-перемешивания. Причем для некоторой положительной константы С\
р
р{ш) < (т + 1)2М°
С
—-- для элементов на одном уровне (i = j),
р для элементов на разных уровнях.
2(то+1 )а
Доказательство этой леммы полностью аналогично доказательству соответствующего утверждения в работе [13].
5. Основная теорема. Докажем асимптотическую нормальность оценки риска. Теорема. Пусть а > 1/2 и функция / регулярна с параметром 7 > (4а — 2)-1. Тогда при пороговой обработке с "универсальным" порогом Т^ имеет место сходимость по распределению:
ЯЛ/)-ЯЛ/) ^(ОД), ./^оо, (16)
где = ^ а константа С не зависит от функции сигнала /.
Доказательство. Из леммы 1 следует, что 01?ЛЯ ~ I)2 = . Разобьем выражение
в (16) на две суммы, как это было сделано в лемме 1:
для - для 4
DJ ¿ъ
7-1 2^ — 1 ,
X) X Л/Л Т3 > - > Т3 > ^.У
3=:р"/+1 к=0 4 _
где р" > (27 + I)-1. Поскольку 7 > (4а — 2)-1, имеем (27 + I)-1 < 1 — (2а)-1. Следовательно, можно выбрать такое р", что (27 + I)"1 < р" < 1 - (2а)"1. Отсюда 1 - а + (а + 2/3)р" < 2/3 + 1/2, и из (10) следует, что первая сумма стремится к нулю п.в.
Далее, действуя, как в лемме 1, нетрудно показать, что
1J 1 1 2
V ^ V ^ , /1 тт-^ / -ж—. Г ^ О 1 -ГТ-* -ж—. Г О Л
Т LJ т
J j=p"J+1 к=о
= зир^ Е (17)
•т 3 к=О
Из леммы 2 следует, что последовательность ] = 0,..., I — 1, к = 1,...,2:г, обла-
дает свойством р-перемешивания и, следовательно [22], обладает свойством а-перемешивания. Наконец, выполнено условие Линдеберга: для любого е > 0
¿2 Е Е
3 к=о
х 1(\)к\Р[Х?к,Т5,о5] - ЕР[Х?к,Т5,о5]\ > еБ^ 0, J ^ оо. (18)
Действительно, так как 1-[Х]к.Т).а)\ < СР3(см. лемму 1), а = то
при а > 1/2, начиная с некоторого </, все индикаторы в (18) обращаются в нуль.
Таким образом, выполнены все условия теоремы 2.1 из работы [23] и справедлива сходимость по распределению (16). Теорема доказана.
СПИСОК ЛИТЕРАТУРЫ
1. Donoho D., Johnstone I. М. Adapting to unknown smoothness via wavelet shrinkage // J. Amer. Stat. Assoc. 1995. 90. P. 1200-1224.
2. Donoho D., Johnstone I. M. Ideal spatial adaptation via wavelet shrinkage // Biometrika. 1994. 81. N 3. P. 425-455.
3. Donoho D., Johnstone I.M., Kerkyacharian G., PieardD. Wavelet shrinkage: asymptopia?//J. R. Statist. Soc. Ser. B. 1995. 57. N 2. P. 301-369.
4. Marron J.S., Adak S., Johnstone I.M., Neumann M.H., Patil P. Exact risk analysis of wavelet regression //J. Comput. Graph. Stat. 1998. 7. P. 278-309.
5. Antoniadis A., Fan J. Regularization of wavelet approximations // J. Amer. Statist. Assoc. 2001. 96. N 455. P. 939-967.
6. Маркин А.В. Предельное распределение оценки риска при пороговой обработке вейвлет-коэффициентов // Информатика и ее применения. 2009. 3. № 4. С. 57-63.
7. Маркин А. В., ШестаковО.В. О состоятельности оценки риска при пороговой обработке вейвлет-ко-эффициентов // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2010. № 1. С. 26-34. (Markin А. V., Shestakov O.V. Consistency of risk estimation with thresholding of wavelet coefficients / / Moscow University Computational Mathematics and Cybernetics. 2010. N 1. P. 22-30.)
8. Шее так о в О. В. Аппроксимация распределения оценки риска пороговой обработки вейвлет-коэффициентов нормальным распределением при использовании выборочной дисперсии // Информатика и ее применения. 2010. 4. № 4. С. 73-81.
9. Шестаков О.В. О точности приближения распределения оценки риска пороговой обработки вейвлет-коэффициентов сигнала нормальным законом при неизвестном уровне шума // Системы и средства информатики. 2012. 22. № 1. С. 142-152.
10. Шестаков О. В. Асимптотическая нормальность оценки риска пороговой обработки вейвлет-коэффициентов при выборе адаптивного порога // Докл. АН. 2012. 445. № 5. С. 513-515. (Shestakov O.V. Asymptotic normality of adaptive wavelet thresholding risk estimation // Doklady Mathematics. 2012. 86. N 1. P. 556-558.)
11. Шестаков О.В. Зависимость предельного распределения оценки риска пороговой обработки вейвлет-коэффициентов сигнала от вида оценки дисперсии шума при выборе адаптивного порога // T-Comm — Телекоммуникации и Транспорт. 2012. № 1. С. 46-51.
12. Шестаков О. В. Центральная предельная теорема для функции обобщенной кросс-валидации при пороговой обработке вейвлет-коэффициентов // Информатика и ее применения. 2013. 7. № 2. С. 40-49.
13. Ерошен ко А. А., Шестаков О. В. Асимптотические свойства оценки риска при пороговой обработке вейвлет-коэффициентов в модели с коррелированным шумом // Информатика и ее применения. 2014. 8. № 1. С. 36-44.
14. Ерошенко А. А., Шестаков О. В. Асимптотическая нормальность оценки риска при вейвлет-вейглет-разложении функции сигнала в модели с коррелированным шумом // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2014. № 3. С. 110-117.
15. Abramovich F., Silverman B.W. Wavelet decomposition approaches to statistical inverse problems // Biometrika. 1998. 85. N 1. P. 115-129.
16. Добеши И. Десять лекций по вейвлетам. Ижевск: НИЦ Регулярная и хаотическая динамика, 2001.
17. Mallat S. A Wavelet Tour of Signal Processing. N.Y.: Academic Press, 1999.
18. Boggess A., Narkowich F. A First Course in Wavelets with Fourier Analysis. Upper Saddle River: Prentice Hall, 2001.
19. T a q q u M.S. Weak convergence to fractional Brownian motion and to the Rosenblatt process // Z. Wahrsehein-lichkeitsth. Verw. Geb. 1975. 31. P. 287-302.
20. Johnstone I.M., Silverman B.W. Wavelet threshold estimates for data with correlated noise // J. R. Statist. Soc. Ser. B. 1997. 59. P. 319-351.
21. Lee N. Wavelet-vaguelette decompositions and homogenous equations. PhD Dissertation. Purdue University, 1997.
22. Bradley R. C. Basic properties of strong mixing conditions. A survey and some open questions // Probab. Surveys. 2005. 2. P. 107-144.
23. Peligrad M. On the asymptotic normality of sequences of weak dependent random variables // J. Theor. Probab. 1996. 9. N 3. P. 703-715.
Поступила в редакцию 17.09.14
LIMIT DISTRIBUTION OF RISK ESTIMATE FOR VAGUELETTE-WAVELET DECOMPOSITION OF A SIGNAL IN THE MODEL OF DATA WITH CORRELATED NOISE
Eroshenko A. A., Kudryavtsev A. A., Shestakov O. V.
In the present paper we consider the problem of estimating a signal function after applying linear homogeneous operator in the model of data with correlated noise. We study asymptotical properties of risk estimate of the thresholding method for vaguelette-wavelet decomposition of a signal. We find the conditions under which the risk estimate is asymptotically normal.
Keywords: wavelets, thresholding, unbiased risk estimate, correlated noise, asymptotic normality.