Среднеквадратичный риск пороговой обработки в моделях с негауссовым распределением шума Текст научной статьи по специальности «Математика»

удк 519.22

0. В. Шестаков1

МИНИМАКСНЫЙ СРЕДНЕКВАДРАТИЧНЫЙ РИСК ПОРОГОВОЙ ОБРАБОТКИ В МОДЕЛЯХ С НЕГАУССОВЫМ РАСПРЕДЕЛЕНИЕМ ШУМА*

В работе рассматривается задача непараметрического оценивания функции сигнала по наблюдениям с шумом с помощью пороговой обработки ее вейвлет-коэффициентов. Вычисляются порядки среднеквадратичного риска и асимптотически оптимальных порогов при общих предположениях о распределении шума.

Ключевые слова: вейвлеты, пороговая обработка, среднеквадратичный риск.

1. Введение. Вейвлет-преобразование разлагает функцию сигнала в линейную комбинацию сдвигов и сжатий некоторой вейвлет-функции и использует свойство "гладкости" сигнала для "экономного" представления данных. Этот принцип лежит в основе популярных методов пороговой обработки подавления шума: небольшие по абсолютной величине коэффициенты вейвлет-разложения считаются шумом, и их обнуление подавляет большую часть шумов, не сильно затрагивая при этом полезный сигнал. В условиях, когда шум предполагается гауссовым, данные процедуры хорошо изучены, и разработаны методы выбора оптимальных пороговых значений при различных предположениях об объеме доступных данных и регулярности функции сигнала [1-4]. В данной работе рассматривается модель вейвлет-коэффициентов с аддитивным шумом, относительно распределения которого делаются лишь самые общие предположения. Для класса регулярных по Липшицу функций приводятся соотношения, позволяющие вычислить асимптотически оптимальный порог и оценить порядок среднеквадратичного риска. Приводятся примеры вычисления порога для различных распределений шумовых коэффициентов.

2. Пороговая обработка вейвлет-коэффициентов. Предположим, что функция сигнала / € L2(R) задана на конечном отрезке и равномерно регулярна по Липшицу с некоторым показателем 7 > 0 и константой Липшица L > 0: / € Lip(7, L). Наблюдатель может регистрировать сигнал в дискретных отсчетах. При этом значения сигнала искажены шумами, которые могут возникать в силу несовершенства оборудования, помех и других причин. Если шум аддитивен, то применение дискретного разложения приводит к тому, что полезный сигнал "упаковывается" в небольшое число относительно больших по модулю вейвлет-коэффициентов, и в силу линейности вейвлет-разложения шум остается аддитивным. Таким образом, в данной работе рассматривается следующая модель "загрязненных" вейвлет-коэффициентов:

хз,к = Hj,k + Zj,k, j = 0,. ..,J-1, к = 0,..., 2J — 1, (1)

где fj,jtk — дискретные вейвлет-коэффициенты "чистого" сигнала, a Zj^ — "шумовые" коэффициенты, относительно которых предполагается, что они независимы и имеют симметричную дифференцируемую функцию распределения Р(zj^ < х) = 1 — Н(х) с конечной дисперсией а2. Индекс j в (1) называется масштабом (или уровнем), а к — сдвигом (см. [5]), J — максимальный уровень разложения, который определяется числом отсчетов сигнала (предполагается, что это число равно 2J). Обозначим через h(x) производную (плотность) функции распределения и дополнительно предположим, что 0 < Н(х) < 1 и h(x) не имеет разрывов второго рода.

Помимо распределения Н(х) модель (1) определяется видом вейвлет-функции ф, участвующей в разложении (см. [5]). Если вейвлет-функция М раз непрерывно дифференцируема (М ^ 7), имеет М нулевых моментов и достаточно быстро убывает на бесконечности вместе со своими

1 Факультет ВМК МГУ, доц.; Институт проблем информатики Федерального исследовательского центра "Информатика и управление" Российской академии наук, ст. науч. сотр., д.ф.-м.н., e-mail: oshestakovQcs.msu.su

* Работа выполнена при частичной финансовой поддержке РФФИ (проект № 16-07-00736).

производными, то найдется такая константа С/ > 0, что [5]

С/2-7/2

\Ы < 2Л7+1/2) • V2)

Далее предполагается, что ф удовлетворяет этим требованиям. Неравенство (2) дает численную оценку "экономности" представления данных при вейвлет-разложении.

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

Обозначим через Xj^ оценку вейвлет-коэффициента, которая вычисляется с помощью некоторой пороговой функции рт(х) с порогом Т: Х^ = рт(Ху^). В данной работе рассматриваются наиболее популярные функции жесткой пороговой обработки р^\х) = ж1(|ж| > Т) и мягкой пороговой обработки Рт\х) = sign(a;)(|ж| —Т)+.

Среднеквадратичный риск пороговой обработки определяется следующим образом:

j=0 к=О

Целью данной работы является поиск асимптотически оптимального в минимаксном смысле порога для обработки наблюдаемого сигнала в классе функций / € Lip(7, L), т. е. порога, обеспечивающего оптимальный порядок среднеквадратичного риска

Rj(T)= sup rj(f, Т). (3)

/GLip(7,L)

Подробное исследование поведения асимптотически оптимального порога для среднеквадратичного риска в модели с аддитивным гауссовым шумом можно найти в работах [2, 3]. Также в [1] был предложен метод поиска адаптивного оптимального порога, с помощью которого можно оценить риск пороговой обработки конкретной функции. Этот метод основан на построении несмещенной оценки риска, статистические свойства которой подробно исследованы в работах [6, 7]. Кроме среднеквадратичного риска также исследовались другие функции потерь, основанные на вероятностях ошибок вычисления вейвлет-коэффициентов, для которых получены асимптотически оптимальные значения порога (см. [4, 8 ,9]). Заметим, что "разумный" порог должен возрастать с ростом J (см. [3]). Однако для упрощения записи далее в явном виде зависимость порога от J указываться не будет.

Далее символом х обозначается порядок рассматриваемой величины по J, т. е. aj х bj, если, начиная с некоторого J, выполнено C\bj ^ a,j ^ C2bj для некоторых положительных констант С\ и С2- Обозначение a,j ~ bj будем использовать в том случае, если lim aj/bj = 1.

J—>СХ)

3. Оптимизация значения порога при жесткой пороговой обработке. Рассмотрим процедуру жесткой пороговой обработки: Х^ =

Пусть функция (¿1(7) > 0 сколь угодно медленно неограниченно возрастает, а функция («Л > 0 сколь угодно медленно стремится к нулю по </. Определим

g1(J) = T + d1(J), g2(J) = mmU2(J), inf {C/(i)}

L te[T-d2(J),T+d2(J)]

где С — некоторая положительная константа, а

оо

= о-! f x2h(x) dx. t h(t) J

Неравенство (2) дает возможность разбить все множество индексов {0,..., 7 — 1} на три класса в зависимости от величины Пусть индексы ^ и ] 1 < ]21 таковы, что

I< ЫА к < 3 < 32 - 1, \н,к\ < ЫА 32 < 3 < 3 - !• (4)

При этом в силу (2)

= ¿ = 1,2, (5)

где для удобства введено обозначение А = 1/(1 + 27). Разобьем Т) на три суммы:

= Е (X,-* -/х,-*) =5!+52 + 53, (6)

3=0 к=о

где в 51 суммирование по ] ведется от 0 до ,7*1 — 1, в ¿>2 — от до з2 — 1, а в 53 — от до 7 — 1. Рассмотрим 53. Для каждого слагаемого, начиная с некоторого </, имеем

оо

Е — = J х2к{х) йх + J х2к{х) йх + J ^2 кк{х) йх.

В силу (4) имеем

-т-шл

ОС -1

! х2к{х) йх + J х2к{х) (¿ж х 2 У х2к{х) йх, (7)

Т-Щ,к -ос Г

т.е. при ^ ^ — 1 величины не влияют на порядок правой части (7). Учитывая, что

число слагаемых в 53 имеет порядок 2,т, получаем

ОО .7—1 Т —

53 х 2'7 ! х2]г(х)ёх+^2 ^ [12кк{х)йх. (8)

г 3=32_т_

1*3,к

Найдем верхнюю оценку для (3) при жесткой пороговой обработке. Для этого оценим слагаемые в суммах 51 и 52 из (6):

т-щ.

к

Е - < я2 + J ^,кНх) <1х.

Далее

т-щл

I ^кЦх) йх = (Щ-Т - Нгк) - /1(1' - Нгк)) = (Т + у)2Н(у) + (Т - у')2Н(у% (9)

-т-шл

где V = —Т—/х^й, у' = Т—• Поскольку второй момент у Н(х) конечен, то х2Н(х) —> 0 при х ^ оо. Следовательно, в силу симметричности распределения Н(х), выражение (9) ограничено величиной

С^Т2, где С$ — некоторая положительная константа. Поэтому, учитывая (5), получаем, что

вг + Б2 ^ С^Т22Х" (д2^))~2Х . (10)

Кроме того, учитывая (2) и (5), нетрудно показать, что для некоторой константы С > О

.1-1 ,1-1

3=32_т_^к 3=32

что имеет меньший порядок роста по </, чем правая часть (10). Следовательно, принимая во внимание (8) и (10), заключаем, что порог Тт\ удовлетворяющий соотношению

ОО

Т~2(д2(</))2А J х2Цх) йх х (И)

является нижней границей для асимптотически оптимального порога. При этом верхняя граница для К,г(Т) задается выражением (10).

Теперь найдем нижнюю границу для Т?/(Т). Имеем

т-щл

Е (х^ь — > о2 ~ J х2к(х) ёх.

-т-Ш,

к

Заметим, что найдется такая функция / € 1лр(7, Ь), что в неравенстве (2) будет достигаться равенство для 0 ^ ] ^ — 1 с некоторой константой С/ (см. [5]). Значит, в силу конечности второго момента у Н(х) для любого е > 0 найдется такое «То, что при J > Jo для всех 0 ^ ] ^ — 1

т-щл

J х2к(х) ёх < е.

-Т-ШЛ

Следовательно, в этом случае все слагаемые в 51 из (6) отделены от нуля и

^ С^Ф = С^{91{1))-2\ (12)

Мь)

где 6Аг — некоторая положительная константа. Пусть порог Т^ удовлетворяет соотношению

сю

(51 (</))2А J ®2М®) Ах х (13)

г

Заметим, что сумма ¿>2 и величина М/ в данных рассуждениях не присутствуют. Это означает, что истинное значение 1?/(Т) имеет порядок не ниже, чем правая часть (12), т.е. рассматриваемый порядок является нижней оценкой для истинного порядка Т?/(Т), а — верхней границей для асимптотически оптимального порога Т.

Таким образом, справедливо следующее утверждение.

Теорема 1. При выборе асимптотически оптимального порога для жесткой пороговой обработки среднеквадратичный риск 1?/(Т) удовлетворяет неравенствам

С^Ы!))-^ < 1Ь(Т) < С^Т22^{д2{1))~2\

где Ст ^ и С$ — некоторые положительные константы. Для асимптотически оптимального значения порога, минимизирующего порядок 1?/(Т) при жесткой пороговой обработке, начиная с некоторого ,1 справедливо неравенство

тС1) <■ т <■ т1^)

-м '

где Т4Л) и Т^ находятся из соотношений (11) и (13) соответственно.

4. Мягкая пороговая обработка. Пусть теперь оценки вейвлет-коэффициентов получаются с помощью мягкой пороговой обработки: Х^ = Определим функции д\{,1) и д2(<Т) так

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

сю

т = (¿_т1)2^) /(ж - т)2/^)йх-

г

Повторяя рассуждения из предыдущего раздела, получаем, что нижней границей для асимптотически оптимального порога является значение Тт \ удовлетворяющее соотношению

сю

Т"2(52(«Т))2А У(я - Т)2Цх) йх х (14)

При этом верхняя граница для -Й/(Т) так же, как и в (10), задается соотношением

НЛТ) ^ С^Т^-7 (д2(,7))~2Х ,

которое отличается от (10) лишь константой и видом функции д2(</).

Для получения нижней границы для -Й/(Т) снова выберем такую функцию / € 1лр(7, Ь), чтобы в неравенстве (2) достигалось равенство для 0 ^ ] ^ — 1 с некоторой константой С/. Имеем

— ос

Е - % У (ж + Т)2Цх) йх+ У {х- Т)2к{х) ёх.

-оо

Заметим, что ^ Т + (¿1(7) при 0 ^ ] ^ — 1. Пусть > 0, тогда

сю

(ж - T)2h{ж) dx ^ У (ж - T)2h{ж) (¿ж ^ С^Т2,

где С™ — некоторая положительная константа. Аналогично рассматривается случай < 0. Таким образом, имеет место оценка

51 ^ С$Т22^ = С$Т22Х7(д1^))-2Х. (15)

Повторяя рассуждения предыдущего раздела, получаем, что порог удовлетворяющий со-

отношению

оо

T~2(gi(J))2X у (ж - T)2h( ж) dx х (16)

г

является верхней границей для асимптотически оптимального порога, а правая часть (15) — нижней границей для Rj(T). Заметим, что нижняя граница в случае мягкой пороговой обработки получилась точнее, чем в случае жесткой.

Теорема 2. При выборе асимптотически оптимального порога для мягкой пороговой обработки среднеквадратичный риск Rj(T) удовлетворяет неравенствам

C^T22XJ(9l(J))-2X < Rj(T) < С$Т22XJ (g2(J)r2X ,

где Cm^ и — некоторые положительные константы. Для асимптотически оптимального значения порога, минимизирующего порядок Rj(T) при мягкой пороговой обработке, начиная с некоторого J справедливо неравенство

t(s) < Т < т^

ш ^ А ^ А м '

где и Т^ удовлетворяют соотношениям (14) и (16) соответственно.

5. Примеры. Рассмотрим конкретные примеры распределения шума. Пусть h(x) х xae~0xß, х —У оо, а G R, в ^ 0, ß > 0.

В случае, когда в ф 0, для мягкой и жесткой пороговых обработок справедливо g2(J) ~ если ß > 1, и <72(«Л = d{J), если 0 < ß ^ 1. Таким образом, поскольку в обоих случаях выполняется Тт ^ ~ Т^ и ~ для асимптотически оптимального порога справедливо соотношение

/1 _ л

которое совпадает с асимптотически оптимальным порогом при минимизации функции потерь, основанной на вероятностях ошибок (см. [9]). В частности, при а = 0, в = 1/(2а2), ß = 2 шум

имеет центрированное нормальное распределение с дисперсией с2, для которого, как показано в [3], асимптотически оптимальный порог имеет порядок

2(1 ^А)1п2-7.

Для нормального распределения можно вычислить более точный порядок асимптотически оптимального порога (см. [8]).

Рассмотрим случай 0 = 0, для которого хвост распределения шума убывает полиномиально и при а < —3 существует конечный второй момент. Асимптотические оценки оптимального порога имеют вид:

If) х ((d2(J))2A, Tff x 2^4

T^ x ((d2(J))2A2(1"A)-7)^T , x

Это означает, что с помощью описанного метода не удается точно определить порядок Т. При этом, если |ск| близко к 3, а показатель у достаточно велик (при этом А близко к нулю), то, как видно из неравенства (2), практически весь полезный сигнал может быть удален при пороговой обработке.

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

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. Minimax estimation via wavelet shrinkage // Ann. Statist. 1998. 26. N 3. P. 879-921.

3. Jans en M. Noise Reduction by Wavelet Thresholding. Lecture Notes in Statistics. Vol. 161. New York: Springer Verlag, 2001.

4. Sadasivan J., Mukherjee S., Seelamantula C. S. An optimum shrinkage estimator based on minimum-probability-of-error criterion and application to signal denoising // Proc. IEEE ICASSP 2014. Florence: IEEE, 2014. P. 4249-4253.

5. Ma 11 at S. A Wavelet Tour of Signal Processing. New York: Academic Press, 1999.

6. Маркин А. В., Шестаков O.B. О состоятельности оценки риска при пороговой обработке вейвлет-коэффициентов // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2010. № 1. С. 26-34. (Mar kin А. V., Shestakov O.V. Consistency of risk estimation with thresholding of wavelet coefficients // Moscow Univ. Comput. Math, and Cybern. 2010. 34. N 1. P. 22-30.)

7. Шестаков O.B. Асимптотическая нормальность оценки риска пороговой обработки вейвлет-коэффи-циентов при выборе адаптивного порога // Докл. АН. 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.)

8. Кудрявцев А. А., Шестаков O.B. Асимптотическое поведение порога, минимизирующего усредненную вероятность ошибки вычисления вейвлет-коэффициентов // Докл. АН. 2016. 468. № 5. С. 487491. (Kudryavtsev A. A., Shestakov O.V. Asymptotic behavior of the threshold minimizing the average probability of error in calculation of wavelet coefficients // Doklady Mathematics. 2016. 93. N 3. P. 295-299.)

Э.Кудрявцев А. А., Шестаков O.B. Асимптотически оптимальная пороговая обработка вейвлет-коэффициентов в моделях с негауссовым распределением шума // Докл. АН. 2016. 471. № 1. С. 11-15. (Kudryavtsev A. A., Shestakov O.V. Asymptotically optimal wavelet thresholding in the models with non-Gaussian noise distributions // Doklady Mathematics. 2016. 94. N 3. P. 615-619.)

Поступила в редакцию 16.05.17