Министерство образования и науки Российской Федерации

УДК

ГРНТИ

Инв. №

УТВЕРЖДЕНО:

Исполнитель:

Государственное учебно-научное учреждение

Факультет вычислительной математики и

кибернетики Московского государственного

университета имени М.В.Ломоносова

От имени Руководителя организации

/ Моисеев Е.И. /

М.П.

НАУЧНО-ТЕХНИЧЕСКИЙ

ОТЧЕТ о выполнении 3 этапа Государственного контракта № 14.740.11.0996 от 23 мая 2011 г.

Исполнитель: Государственное учебно-научное учреждение Факультет вычислительной математики и кибернетики Московского государственного университета имени М.В.Ломоносова Программа (мероприятие): Федеральная целевая программа «Научные и научнопедагогические кадры инновационной России» на 2009-2013 гг., в рамках реализации мероприятия № 1.2.2 Проведение научных исследований научными группами под руководством кандидатов наук.

Проект: Математические методы анализа и обработки стохастических изображений.

Руководитель проекта:

/Шестаков Олег Владимирович (подпись) Москва 2012 г.

СПИСОК ОСНОВНЫХ ИСПОЛНИТЕЛЕЙ

по Государственному контракту 14.740.11.0996 от 23 мая 2011 на выполнение поисковых научно-исследовательских работ для государственных нужд Организация-Исполнитель: Государственное учебно-научное учреждение Факультет вычислительной математики и кибернетики Московского государственного университета имени М.В.Ломоносова Руководитель темы:

Шестаков О. В.

кандидат физикоматематических наук, подпись, дата доцент Исполнители темы:

Кудрявцев А. А.

кандидат физикоматематических наук, без подпись, дата ученого звания Мизотин М. М.

кандидат физикоматематических наук, без подпись, дата ученого звания Насонова А. А.

без ученой степени, без ученого звания подпись, дата Насонов А. В.

кандидат физикоматематических наук, без подпись, дата ученого звания Семашко А. С.

без ученой степени, без ученого звания подпись, дата Сторожилова М. В.

без ученой степени, без ученого звания подпись, дата Сергеев В. В.

без ученой степени, без без ученой степени, без без ученой степени, без Отчет 65 с., 1 ч., 14 рис., 2 табл., 36 источн., 0 прил.

стохастические изображения, томографические изображения, локальные дескрипторы, ортогональные многочлены, разреженные предстваления В отчете представлены результаты исследований, выполненных по этапу Государственного контракта № 14.740.11.0996 "Математические методы анализа и обработки стохастических изображений." (шифр "2011от 23 мая 2011 по направлению "Проведение научных исследований научными группами под руководством кандидатов наук в следующих областях:- математика;- механика" в рамках мероприятия 1.2. "Проведение научных исследований научными группами под руководством кандидатов наук.", мероприятия 1.2 "Проведение научных исследований научными группами под руководством докторов наук и кандидатов наук", направления 1 "Стимулирование закрепления молодежи в сфере науки, образования и высоких технологий." федеральной целевой программы "Научные и научно-педагогические кадры инновационной России" на 2009годы.

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

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

Для проведения вычислительных экспериментов на втором была задействована вычислительная и программная база кафедры математической статистики и лаборатории математических методов обработки изображений факультета ВМК МГУ имени М.В. Ломоносова.

На третьем этапе работ получены следующие результаты:

1. Методы нелинейной пороговой обработки сигналов и изображений, направленных на подавление шума.

2. Метод обработки и анализа изображений на основе разреженных представлений.

3. Метод индексации и поиска изображений на основе систем ортогональных многочленов.

СОДЕРЖАНИЕ

Введение

Проведение III этапа исследований по проблеме «Математические методы анализа и обработки стохастических изображений»

Разработка методов нелинейной пороговой обработки сигналов и изображений, направленных на подавление шума

Разработка методов обработки и анализа изображений с помощью разреженных представлений

Создание методов обработки и анализа изображений на основе систем ортогональных многочленов.

Заключение

Список использованных источников

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

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

Математический аппарат разреженных представлений динамично развивается и все больше и больше применяется для решения задач обработки изображений, видео и др. сигналов. В то время как в некоторых типах задач (например, шумоподавлении, сжатии, восстановлении потерянных данных) методы разреженных представлений применяются использования – например, для классификации изображений или их частей.

классификации.

Среди различных применений томографических методов особенно важным является их применение в медицине. Магнитно-резонансная томография (МРТ) является незаменимым современным диагностическим средством и уже получила широкое распространение. Число медицинских изображений быстро растет в связи с появлением новых методов оборудования для проведения обследований, и ручной анализ такого огромного объема данных очень трудоемок. В литературе можно найти результаты исследований, которые показали, что эффективность работы докторов увеличивается при использовании поисковых систем по изображениям. В связи с этим, задача поиска МРТ является актуальной.

В разделе 3 рассматривается применение дескрипторов на основе ортогональных многочленов к поиску медицинских изображений, в частности для поиска МРТ изображений мозга и обнаружения болезни Альцгеймера (БА). Наличие у радиологов и нейрологов «под рукой» МРТ изображений похожих случаев позволяет вынести более обоснованное решения о состоянии здоровья и диагнозе.

Проведение III этапа исследований по проблеме «Математические методы анализа и обработки стохастических сигналов и изображений, направленных на подавление шума.

1.1 Метод реконструкции изображений В данном разделе на примере интегрального преобразования, нелинейной пороговой обработки, который также применим к широкому кругу других приложений, связанных с подавлением шума в регистрируемых данных. Томографические методы реконструкции изображений получили широкое распространение в медицине, биологии, физике плазмы, газовой динамике, геофизике, астрономии и радиолокации (см., например, [1–3]). Во реконструкции изображений играет преобразование (или оператор) Радона:

где f ( x, y ) L2 ( 2 ) - функция с компактным носителем, описывающая интегрируется функция f ( x, y ) (здесь s - расстояние от начала координат до прямой, а - угол, образованный с осью x перпендикуляром, опущенным из начала координат на эту прямую). При фиксированном значении функцию Rf (, s) также называют проекцией функции f ( x, y ) под углом. Для преобразования Радона справедлива так называемая теорема о центральном сечении, играющая ключевую роль во многих методах реконструкции изображений по проекциям:

где Rf (, ) - одномерное преобразование Фурье функции Rf (, s) по второй переменной, а f (1, 2 ) - двумерное преобразование Фурье функции На практике проекционные данные регистрируются с некоторой погрешностью (шумом), возникающей из-за несовершенства оборудования и различных случайных помех. И поскольку задача обращения преобразования Радона относится к классу так называемых некорректно поставленных задач, при реконструкции изображений приходится использовать методы регуляризации. В работах [4–7] рассматривается метод регуляризации, заключающийся в разложении функции изображения по двумерному ортонормированному вейвлет-базису с последующей пороговой обработкой коэффициентов разложения. При этом функция Rf раскладывается по так ортонормированным. Однако в силу «почти ортогональности» вейглетов среднеквадратичный риск метода оценивается при дополнительных предположениях, эквивалентных условию ортонормированности вейглетбазиса.

Методы анализа среднеквадратичного риска и его оценки в задаче реконструкции томографических изображений во многом опираются на тот факт, что преобразование Радона обладает многими свойствами линейных однородных операторов, которые рассматривались в работах [8] и [9]. В коэффициентов вейглет-вейвлет разложения сигнала, и показывается, что при некоторых условиях регулярности оценка среднеквадратичного риска этого метода асимптотически нормальна. Далее мы представим метод обращения преобразования Радона, альтернативный методу вейвлет-вейглет разложения и исследуем свойства оценки среднеквадратичного риска этого метода.

Разложим функцию Rf (, s) в ряд Фурье по переменной на [0,2 ].

вещественный базис Фурье Для сокращения записи обозначим элементы этого базиса через {Z n }=0, перенумеровав их соответствующим образом. Имеем (x) - некоторая вейвлет-функция. Индекс j называется масштабом, а индекс k - сдвигом. Семейство { j,k } j,kZ образует ортонормированный базис в L2 ( ) (см. [10]). Далее разложим функции Rf n (s ) по базису После применения к разложению (1.3) обратного преобразования R 1, получаем представление для функции f, которое является аналогом вейглетвейвлет разложения, рассмотренного в работах [8], [9] и [11]:

Последовательность {u n, j,k } не образует ортонормированную систему, гладкости, то последовательность {u n, j,k } образует устойчивый базис, т.е.

существуют такие константы 0 < A B <, что для всех квадратично суммируемых последовательностей {cn, j,k }.

Иногда свойство (1.5) называют «почти ортогональностью» (см. [7]).

Функции u n, j,k по аналогии с терминологией работ [5] и [11] назовем «вейглетами». Справедливо следующее утверждение, являющееся полным аналогом леммы 1.1 из работы [8].

Лемма 1.1 Пусть существуют такие константы A > 0, a > 1/2 и b > 3/2, что для всех, тогда последовательность {u n, j,k } образует устойчивый базис в L2 ( ).

Доказательство. Требуется показать, что существуют константы 0 < A B < такие, что для произвольной функции f L2 ( ) выполнено (это соотношение эквивалентно условию (1.5)).

Определим оператор обратного проецирования, действующий на функцию g (, s ) :

Оператор R* является сопряженным к оператору R (см. [13]).

Обозначим Таким же образом, как в доказательстве леммы 1.1 из работы [8], положительных констант C1 и C2, что Найдем константы n, j,k. Используя теорему о центральном сечении (1.2), имеем где через обозначена функция, для которой пропорциональны 2 j/2. Докажем теперь (1.9). Имеем периодична с периодом 2 j 1. Следовательно, Далее, поступая, как в доказательстве леммы 1.1 из работы [8], получаем, что при выполнении условия (1.6) для некоторой положительной константы C1 выполнено о центральном сечении Следовательно, (1.9) выполнено. Докажем теперь (1.10).

Здесь через обозначена функция, для которой ( ) = ( ) Далее, поступая, как при доказательстве (1.9), получаем, что найдется положительная константа C2 такая, что Лемма доказана.

среднеквадратичного риска изображения, имеет компактный носитель и равномерно регулярна по Липшицу с некоторым показателем > 0. Оказывается, что при этом функция Rf равномерно регулярна по Липшицу с показателем 1/2 (см.

[13] и [14]). Также будем считать, что вейвлет-функция удовлетворяет условиям теоремы 6.3 из [10] и леммы 1.1. Таким образом, найдется такая константа A > 0, что при всех неотрицательных j и всех k.

коэффициенты получаются умножением матрицы значений функции Rf на ортогональную матрицу F и ортогональную матрицу W, определяемые базисом Фурье {Z n } и вейвлет-базисом { j,k }. Причем, поскольку в измеряемых данных:

где J - некоторое положительное число (“разрешение” изображения), X u,v - наблюдаемые данные, R - преобразование Радона, f - истинная (незашумленная) функция изображения, а u,v - независимые случайные погрешности, которые имеют одинаковое нормальное распределение с нулевым средним и дисперсией 2. Таким образом, в силу ортогональности описываются следующей моделью:

независимы и имеют такое же распределение, как и u,v.

Поскольку в измерениях присутствует шум, необходимо использовать некоторые процедуры для его удаления. Мы рассмотрим процедуру пороговой обработки (см. [10]). Смысл пороговой обработки вейвлеткоэффициентов измеряемого сигнала заключается в удалении достаточно маленьких коэффициентов, которые считаются шумом. Будем использовать так называемую мягкую пороговую обработку с порогом T j, зависящим от масштаба j. Для подавления шума к каждому коэффициенту Yn, j,k применяется функция T (x ), определяемая соотношением При такой пороговой обработке коэффициенты, которые по модулю меньше порога T j, обнуляются, а абсолютные величины остальных коэффициентов уменьшаются на величину порога. В результате функция сигнала f в разложении (1.4) (точнее, ее нормированная версия) оценивается следующим образом:

Среднеквадратичная ошибка (риск) оценки f определяется как Так как система функций {u n, j,k } не является ортонормированной, риск нельзя представить в виде суммы вкладов погрешностей отдельных коэффициентов, однако в силу устойчивости базиса {u n, j,k } найдутся такие константы A и B, зависящие от, что ArJ r BrJ, где (здесь T = (T0,, TJ 1 ) - вектор порогов). Таким образом, с точностью до множителя величина r асимптотически ведет себя так же, как rJ, поэтому в дальнейшем мы будем рассматривать асимптотическое поведение В выражении (1.15) присутствуют неизвестные величины n, j,k, поэтому вычислить значение rJ нельзя. Однако его можно оценить. В каждом слагаемом если Yn, j,k > T j, то вклад этого слагаемого в риск составляет n, j,k ( 2 T j2 ), а если Yn, j,k T j, то вклад составляет n, j,k n, j,k.

Поскольку EYn2, j,k = 2 n, j,k, значение n, j,k можно оценить разностью Таким образом, в качестве оценки риска можно использовать следующую величину:

причем эта оценка является несмещенной (см. [10]).

1.3 Предельные теоремы для оценки риска По аналогии с работами [6], [8] и [9] мы рассмотрим ситуацию, когда для каждого масштаба j выбирается “универсальный” порог. Такой порог масштабе количество коэффициентов равно 2 j J, и получил название «универсальный», так как он не зависит от наблюдаемых данных и при выборе этого порога риск rJ близок к минимальному (см. [12]).

Для оценки риска (1.16) справедлива следующая теорема, являющаяся аналогом теоремы 1.1 из [8].

теоремы 6.3 из [10] и леммы 1.1, а функция f L2 ( ) имеет компактный носитель и равномерно регулярна по Липшицу с показателем > 0, тогда Зачастую дисперсия шума 2 не известна и ее также необходимо оценивать с помощью некоторой оценки 2. При этом выражение (??) принимает вид (обозначим вектор этих порогов через TU ).

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

Далее мы рассмотрим оба случая.

В случае, когда 2 строится по независимой выборке, для оценки риска (1.18) справедливо следующее утверждение.

теоремы 6.3 из [10] и леммы 1.1, а функция f L2 ( ) имеет компактный носитель и равномерно регулярна по Липшицу с показателем > 0. Пусть 2 - не зависящая от Yn, j,k оценка дисперсии 2, для которой выполнено где (x ) - функция распределения нормального закона с нулевым средним и дисперсией 2. Тогда где (x ) - функция распределения нормального закона с нулевым средним и дисперсией Если дисперсия оценивается по коэффициентам разложения Yn, j,k и функция f удовлетворяет требуемым условиям регулярности, то в силу (1.12) ее оценивают по половине всех коэффициентов на уровне j = J 1, так как эти коэффициенты фактически содержат только шум. Мы рассмотрим популярные оценки неизвестной дисперсии: S - выборочная дисперсия, R - интерквартильный размах и M - абсолютное медианное отклонение от медианы. Эти оценки определяются следующими соотношениями:

где Y2 J 1,(1/4) и Y2 J 1,(3/4) - выборочные квантили порядка 1/4 и 3/4, построенные по набору коэффициентов Yn, J 1,k, 0 n 2 J 1, 0 k 2 J 1 1, 3/4 - теоретическая квантиль порядка 3/4 стандартного нормального распределения ( 3/4 0.6745 ), а « med » обозначает выборочную медиану.

теоремы 6.3 из [10] и леммы 1.1, а функция f L2 ( ) имеет компактный носитель и равномерно регулярна по Липшицу с показателем > 1/2. Пусть оценка дисперсии 2 определяется выражением (1.21). Тогда где (x ) - функция распределения нормального закона с нулевым средним и дисперсией 2 = 2/9.

теоремы 6.3 из [10] и леммы 1.1, а функция f L2 ( ) имеет компактный носитель и равномерно регулярна по Липшицу с показателем > 1. Пусть оценка дисперсии 2 определяется выражением (1.22) или (1.23). Тогда где (x ) - функция распределения нормального закона с нулевым средним и дисперсией Доказательство теорем 1-4 полностью аналогично доказательству соответствующих теорем из работ [8] и [9].

Замечание 1.2 Вектор порогов TU фактически используется для построения оценки функции Rf. В работе [11] вместо вектора порогов TU рассматривается единый для всех масштабов порог (в данном случае он показывается, что при построении оценки функции f этот порог является асимптотически близким к оптимальному в минимаксном смысле. При рассмотрении порога TR все представленные результаты останутся неизменными.

2 Разработка методов обработки и анализа изображений с помощью разреженных представлений Методы разреженных представлений применимы во многих задачах обработки изображений – сжатии, подавлении шума, ретушировании и др.

При этом изображение, закодированное в виде вектор-столбца y из n элементов, представляется в виде линейной комбинации нескольких изображений, называемых атомами, из определенного набора изображений – разреженным, если число ненулевых элементов вектора x гораздо меньше числа атомов (столбцов матрицы A); при этом необходимо учитывать величину ошибки представления ||y – Ax||2.

Будем обозначать ||x||0 число ненулевых элементов в x, а supp x – множество индексов ненулевых элементов x (носитель вектора х), таким образом, ||x||0 = |supp x|.

2.1 Построение разреженных представлений представления данного изображения y с данным словарем A возможна в двух вариантах:

а) с ограничением количества k0 атомов, используемых для построения представления:

б) с ограничением ошибки представления :

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

приближенно оптимального разреженного представления, среди которых – так называемые жадные методы. Базовый алгоритм OMP (orthogonal matching pursuit), использующийся в этих методах, основан на итеративном расширении носителя решения. На каждой итерации выбирается элемент, добавление которого в носитель позволяет максимально уменьшить ошибку представления.

Введем обозначения: k – номер итерации, x k – решение, r k – невязка, S k = supp x k – носитель решения, ai – i-й столбец матрицы A. Тогда алгоритм OMP можно записать следующим образом:

2) на каждой итерации увеличиваем k на 1 и выполняем следующие шаги:

i S k 1 вычисляем минимальную норму невязки, которую можно достичь, добавляя в носитель решения i-й элемент:

добавляем элемент, позволяющий достичь минимальную норму невязки:

вычисляем новое решение:

где AS обозначает матрицу, состоящую из столбцов матрицы A с номерами, входящими в множество S, а xS обозначает вектор из элементов вектора x с номерами из этого множества.

проверяем критерии останова:

k = k0, если задана максимальная мощность носителя y Ax k, если задана максимальная ошибка представления.

В данной работе применялся вариант вышеописанного алгоритма, позволяющий в некоторых случаях получить более точное решение. Его отличие от базового алгоритма в том, что для выбора наилучшего элемента для добавления при вычислении остаточных ошибок (i) представления варьируются значения всех элементов решения, входящих в носитель.

2.2 Построение словаря Для применения разреженных представлений крайне важно выбрать способ построения словарей, используемых в решении задачи. Можно выделить два подхода к выбору словарей:

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

представления более быстрыми, чем стандартные, алгоритмами, но в большинстве случаев для реальных изображений представления, получаемые с такими словарями, не могут быть сильно разреженными (т.е. мощность обрабатываемые изображения принадлежат к достаточно узкому классу.

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

Для построения словарей в данной работе применялся алгоритм KSVD, позволяющий построить словарь заданного размера, наилучшим образом описывающий заданный набор изображений в расчете на заданную мощность носителя представлений для изображений из набора.

Обозначим изображения обучающей выборки (y1, y 2,, y d ) = Y Rnd.

изображений обучающей выборки (x1, x2,, xd ) = X Rmd.

Формальные требования можно записать следующим образом:

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

K-SVD – эвристический алгоритм, основанный на итеративном улучшении словаря. Изначально словарь задается произвольным образом (например, в качестве атомов берется часть изображений из обучающей выборки), и по нему каким-либо алгоритмом вычисляются начальные представления xi для yi. Под улучшением словаря понимается минимизация суммы квадратов ошибок представлений изображений обучающей выборки:

оптимальное значение для каждого из атомов словаря одновременно с соответствующими элементами представлений, использующих этот атом.

Перепишем минимизируемую функцию в виде Здесь и далее x j обозначает не разреженное представление для y j, а строку из -ых элементов разреженных представлений всех изображений обучающей выборки. Задачу минимизации можно теперь рассматривать, как задачу поиска наилучшего приближения для матрицы E j Y матрицей ранга 1. При этом, варьируя x j, необходимо сохранять не только количество нулей в нем, но и их расположение. Для этого зададим матрицу проекции P j, такую, что вектор-строка x j P j содержит только ненулевые приближение матрицей ранга 1 можно найти, вычислив сингулярное разложение:

и использовав наибольшее сингулярное число:

Вместо ресурсоемкого вычисления сингулярного разложения иногда этот шаг заменяют на несколько итераций поочередной оптимизации по x R и a j с помощью метода наименьших квадратов:

Критерий остановки алгоритма K-SVD – изменение минимизируемой функции становится меньше заданного порога.

2.3 Применение разреженных представлений для задачи обнаружения хориоретинальной атрофии В данной работе методы разреженного представления применялись для решения задачи выявления областей хориоретинальной атрофии (ХРА) сетчатки на изображениях человеческого глазного дна. При автоматическом анализе фотографий сетчатки информация об областях ХРА позволяет точнее определять границу диска зрительного нерва. Это важно, так как характеристики диска зрительного нерва (ДЗН) составляют основные критерии, по которым диагностируется глаукома, и, кроме того, размеры ДЗН часто служат опорой при определении размеров других важных составляющих сетчатки. Также, выраженность ХРА или ее отсутствие позволяют делать более точные предположения о наличии у пациента специфических болезней, таких как диабетическая ретинопатия и глаукома.

здорового глазного дна. атрофированного участка зрительного нерва, также показана граница области lutea) – желтое пятно.

Fovea – центральная ямка макулы, содержащая только фоторецепторы.

Рисунок 2.1 – Пример изображения глазного дна и Наиболее распространен в литературе подход для детектирования областей атрофии, основанный на классификации по текстурным признакам.

Изображение разбивается на блоки небольшого размера, для каждого из которых вычисляется некоторое количество признаков (чаще всего эти признаки - характеристики матрицы смежности), после чего к полученным данным применяется какой-либо классификатор. Можно либо говорить на основании результатов о наличии или отсутствии в целом заметной патологии, либо использовать результаты классификации отдельных блоков в дальнейшем анализе изображения. В первом случае, хотя и достижима точность классификации порядка 80%, но практическая ценность такого результата крайне мала. Во втором же случае достаточно надежный алгоритм позволял бы получить упомянутые выше преимущества от детектирования ХРА, но точность классификации отдельных блоков слишком мала (порядка 60%). При такой низкой точности использование результатов детектирования в других этапах анализа изображения может привести к ухудшению их результатов, поэтому на практике в системах автоматического анализа изображений глазного дна детектирование областей ХРА обычно не производится.

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

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

где n - количество классов, ei - ошибка представления изображения по словарю i-го класса, pi - вероятность того, что изображение принадлежит к iму классу.

Для обучения классификаторов была использована база из изображений, на которых специалист-офтальмолог вручную отметил области ХРА, для тестирования - база из 130 изображений. Все изображения были разбиты на блоки 8х8 пикселей со всеми возможными перекрытиями. Были отобраны блоки, более чем на половину принадлежащие областям атрофии, для построения словаря класса А. Для построения словаря класса Б были отобраны блоки, целиком лежащие внутри области видимости камеры и не имеющие пересечений с областями атрофии. Для усиления значимости величины ошибки представления как показателя принадлежности к классу были вычислены усредненные изображения всех блоков, принадлежащих к каждому из классов. Эти усредненные изображения были вычтены из всех блоков, входящих в обучающие выборки для соответствующих классов.

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

Для каждого класса методом K-SVD были построены словари из 100, 200 и 400 атомов со значениями параметра k0 3, 4, 6 и 8. Для полученных пар словарей была оценена точность классификации с использованием параметры составляют: размер словаря – 200 атомов, k0 = 4; точность классификации при этом составила 59,3%.

Рисунок 2.2 – Результаты оценки точности классификации для словарей размером в 100, 200 и 400 атомов в зависимости от параметра k0.

Увеличение размера словаря практически не дает прироста точности, а уменьшению точности из-за того, что классы А и Б не слишком сильно отличаются (при достаточно большом k0 для блоков обучающей выборки класа А возможно достаточно точное представление по словарю класса Б, и наоборот).

2.3.2 Классификатор на основе текстурных признаков Текстурные признаки также вычислялись для всех блоков 8х8 пикселей со всеми перекрытиями; критерии принадлежности блоков к классам при обучении были те же, что описаны в предыдущем разделе.

Для каждого блока вычислялись матрицы смежности уровней яркости С1,1, С1,-1, С-1,1, С-1,-1 для синего канала I изображения с квантованием на уровня по формуле:

Для каждой из четырех матриц вычислялись следующие значения:

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

Такой набор признаков был выбран с помощью поэтапного исключения малозначащих и нерелевантных признаков. Изначально характеристики для перечисленных свойств матриц смежности, рассматривался также второй угловой момент. В первую очередь была проверена гипотеза, что уменьшение числа анализируемых каналов не сильно ухудшит качество классификации, так как при такой редукции сильнее всего уменьшается объем вычислений. Проверка производилась с использованием наибольших, наименьших и усредненных значений каждого из признаков по четырем матрицам. При сравнении четырех групп признаков: RGB (63 признака), R, G и B (по 21 признаку для каждого из каналов по отдельности) были получены следующие значения точности классификации: 65,6% (RGB), 63,2% (R), 65,1% (G), 65,4% (B). Поэтому решено было использовать признаки только синего канала изображения. Дальнейший переход к шести наиболее значащим признакам при увеличении скорости работы в 2,4 раза привел к падению точности классификации всего на 0,5%. Таким образом, точность классификатора на основе текстурных признаков составила 64,9%.

2.3.3 Фильтрация результатов классификации Описанные в предыдущих разделах классификаторы оперируют перекрывающимися блоками размером 8х8 пикселей. Чтобы получить область, классифицированную как атрофическую, необходимо для каждого пикселя учесть результаты классификации всех блоков, включающих данный пиксель (их может быть меньше, чем 64, для пикселей, находящихся рядом с границей поля зрения и совместить результаты классификации для соседних пикселей (атрофированный регион должен быть сплошным и с гладкими границами; он не может, например, состоять из пикселей, расположенных в шахматном порядке).

Таким образом, для каждого пикселя имеется несколько ответов классификатора, обозначающих правдоподобность того, что соответствующие блоки изображения принадлежат области атрофии.

Среднее значение этих ответов приписывалось каждому пикселю, затем полученная «карта вероятности атрофии» сглаживалась гауссовским фильтром с = 5. Результат бинаризовался с порогом, равным 0,5, и подвергался морфологическому закрытию со структурным элементом в виде круга, радиус которого равен 1/50 диаметра поля зрения камеры. Полученное бинарное изображение и считалось маской области атрофии – итоговым результатом работы классификатора. Именно эти результаты использовались для оценки качества работы классификаторов; несколько таких результатов приведено ниже в качестве примеров.

2.3.4 Составной классификатор Для совместного использования двух базовых классификаторов была использована несимметричная схема, в которой результат работы классификатора, основанного на разреженных представлениях – ошибки представления блока по двум словарям – использовался в качестве пары дополнительных параметров для классификатора, основанного на текстурных признаках. Таким образом, для блоков, имеющих значительные ошибки представления по обоим словарям, соотношение этих ошибок автоматически учитывалось с малым весом, и наоборот, при ненадежности ответа, полученного по текстурным признакам, результат классификации по ошибкам разреженных представлений имел больший вес. Такой подход позволил достичь точности классификации в 69,7%. Примеры результатов, полученных с помощью описанного алгоритма, приведены на рисунке 2.3.

Рисунок 2.3 – Примеры результатов классификации. Синим цветом обозначены границы детектированных областей атрофии, зеленым – эталонные границы областей атрофии, черным – границы диска зрительного нерва. Перед фильтрацией результатов классификации исключались пиксели, принадлежащие диску зрительного нерва.

3 Создание методов обработки и анализа изображений на основе систем ортогональных многочленов.

особенностей Извлечение ключевых точек и их параметризация - одна из основных низкоуровневых задач обработки изображения. Например, это - начальный шаг в задачах сопоставления стереоизображений [15], распознавания объектов [16], индексация видео и изображений [17], склейки панорам и других. Существует большое количество подходов к задаче выделения ключевых точек, например, детектор уголков Харриса [18], детектор пятен на основе разности гауссиан (DoG), предложенный Lowe [19], подход, основанный на круговых гармонических функциях [20, 21], и другие. Задача построения дескриптора ключевой точки также широко представлена в литературе [19, 20, 22]. В построения дескриптора основным требованием фотометрических преобразований. Это требование крайне важно для того чтобы получить высокий процент правильных соответствий ключевых точек на изображениях, полученных с разных ракурсов и при различных условиях освещения.

выделению ключевых точек и построению дескрипторов на основе круговых гармонических функциях Гаусса-Лагерра (КГФ ГЛ) [20] и применяется к проблеме поиска МРТ изображений мозга. Дескрипторы на основе КГФ ГЛ имеют несколько преимуществ по сравнению с другими дескрипторами, особенно для рассматриваемой задачи. Например, они не содержат основанную на гистограмме часть, как, например, в дескрипторе SIFT, которая делает дескриптор инвариантным к небольшим локальным сдвигам определенных прикладных задачах, такая инвариантность нежелательна, практически отсутствуют. Как показано в работах [23, 24], дескрипторы без применения гистограммы в некоторых случаях превосходят SIFT, который на данный момент является де-факто текущим стандартом качества.

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

Другие преимущества и свойства дескрипторов Гаусса-Лагерра рассматриваются в разделе 3.2. Применение дескрипторов и постановка основной задачи приводится в разделе 3.3. В разделе 3.4 описан предложенный алгоритм для поиска МРТ мозга. В разделе 3.5 приводится алгоритм предварительной классификации, призванный улучшить общую эффективность работы системы. Наконец, результаты представлены в разделе Ошибка! Источник ссылки не найден..

3.2 Ключевые точки и дескрипторы на основе круговых гармонических функций Гаусса-Лагерра (КГФ ГЛ) 3.2.1 Обнаружение ключевых точек полярных координатах функций:

Их радиальные профили - функции Лагерра:

где n 0,1,; 0, 1, 2... и L (x) полиномы Лагерра:

Функции Лагерра n (x) может быть вычислены, используя следующие рекуррентные соотношения:

Применение рекуррентных формул позволяет избежать проблем потери точности при умножении экспонент на степенную функцию.

Эти функции n, называются круговыми гармоническими функциями Гаусса-Лагерра (КГФ ГЛ) и имеют 2 индекса, определяющих радиальный ( n 0,1,4; 1, 2,..., 5 ) приведены на рис. 3.1.

Рисунок 3.1 – Действительная часть n ( n 0,1,4; 1, 2,..., 5 ).

Функции КГФ ГЛ являются самоуправляемыми, то есть они могут быть повернуты на угол путем умножения на ei. Они также сохраняют свою форму при применении преобразования Фурье. Также, они подходят для многомасштабного и многокомпонентного анализа изображения [20, 25].

действительной плоскости R2. Вследствие ортогональности системы функций n, изображение I ( x, y ) может быть представлено в окрестности точки x0, y0 для любого фиксированного масштаба в декартовой системе координат следующим образом:

[ 2 smax,2 smax ] дискретизовано в (2smax + 1) октавы, каждая октава в свою очередь содержит Ns масштабов с равным шагом. Таким образом, набор весов определяется как { j }, где j 0,1,..., 2 N s (2smax 1) 1. Принимая во внимание свойства функций КГФ ГЛ, оказыватся, что при значениях параметра n 0, 3,4 они являются детекторами особенностей типа «развилка» и «крест».. Набор 2Ns(2smax + 1) откликов детектора определяется как:

и назевается скалограммой. Скалограмма просматривается скользящим 3D окном (5x5x3 пикселей). Кандидаты ключевых точек K ( x, y, ) определяются как локальные максимумы скалограммы в пределах окна.

Здесь ( x, y ) координата ключевой точки и соответствующий масштаб ключевой точки. Таким образом, получается набор ключевых точек изображения {K }. Этот набор прореживается, выбрасывая те ключевые точки K у которых одни и те же координаты ( x, y ), но на разных (два и более) масштабах. И, наконец, отбрасываются ключевые точки K, в которых отклик детектора S ( x, y; ) по величине меньше, чем выбранный порог:

3.2.2 Дескрипторы на основе КГФ ГЛ Здесь мы кратко рассматриваем алгоритм построения дескрипторов ключевых точек, предложенный в [6]. Любая ключевая точка K ( x, y, ) связана с локальным дескриптором { (n,, j )}. В данном случае дескриптор – это комплекснозначный вектор, состоявший из локальной проекции изображения на систему функций КГФ ГЛ n на 2 jmax масштабах, соседних с масштабом определяются как: