Программирование, 2021, № 3, стр. 64-72

ГИБРИДНЫЙ МЕТОД ПОДАВЛЕНИЯ ОСЦИЛЛЯЦИЙ ГИББСА НА ИЗОБРАЖЕНИЯХ МАГНИТНО-РЕЗОНАНСНОЙ ТОМОГРАФИИ

М. А. Пенкин a*, А. С. Крылов a**, А. В. Хвостиков a***

a Московский государственный университет имени М.В. Ломоносова, Факультет вычислительной математики и кибернетики, Лаборатория математических методов обработки изображений
119991 Москва, Ленинские горы, д. 1, Россия

* E-mail: penkin97@gmail.com
** E-mail: kryl@cs.msu.ru
*** E-mail: khvostikov@cs.msu.ru

Поступила в редакцию 10.10.2020
После доработки 20.10.2020
Принята к публикации 12.01.2021

Полный текст (PDF)

Аннотация

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

1. ВВЕДЕНИЕ

Подавление осцилляций Гиббса (ложного оконтуривания, англ. ringing) – актуальная задача математических методов обработки изображений. Они часто наблюдаются при изменении разрешения изображений, при повышении резкости изображений, а также при визуализации данных магнитно-резонансной томографии (МРТ) (см. рис. 1). Повышение качества МРТ изображений является важной для медицинской диагностики задачей, например, при работе с атласами мозга [1].

Рис. 1.

Примеры осцилляций Гиббса на изображениях МРТ (указаны стрелками).

Феномен Гиббса был открыт Генри Уилбрахамом в 1848 г., а затем повторно открыт Дж. Уиллардом Гиббсом в 1898 г. Интересным является факт, что в отличие от других математических областей, исследование эффекта Гиббса было не очень активным вплоть до 1977 года [2]. Затем работы на эту тему стали появляться чаще [35]. Отчасти такая динамика объясняется развитием магнитно-резонансной томографии (МРТ), основанной в 1973 году, когда профессор химии Пол Лотербург опубликовал в журнале Nature статью “Создание изображения с помощью индуцированного локального взаимодействия; примеры на основе магнитного резонанса”. Математические методы получения изображений МРТ были затем усовершенствованы Питером Мэнсфилдом, и в 2003 году обоим исследователям была присуждена Нобелевская премия по физиологии и медицине.

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

$\xi (t) = \left\{ \begin{gathered} a,\quad t \in [ - \tau {\text{/}}2,\tau {\text{/}}2] \hfill \\ 0,\quad t \in [ - T{\text{/}}2,T{\text{/}}2]{{\backslash }}[ - \tau {\text{/}}2,\tau {\text{/}}2] \hfill \\ \end{gathered} \right.,$
$\xi (t) = \xi (t + T),$
где $\xi (t)$ – модельный пример конечного разрыва; a – величина скачка; $T$ – период модельной функции; τ – ширина импульса.

Используя преобразование Фурье в комплексной форме и предполагая что $T = 2\tau $, (1.1) может быть записано в виде:

(1.2)
$\xi (t) = \sum\limits_{k = - \infty }^{ + \infty } {{{d}_{k}}} {{e}^{{i{{\omega }_{k}}t}}},$
где ${{d}_{k}} = \tfrac{1}{T}\int_{ - \tfrac{T}{2}}^{\tfrac{T}{2}} {} \xi (t){{e}^{{ - i{{\omega }_{k}}t}}}dt$, ${{\omega }_{k}} = \Omega k$, $\Omega = 2\pi {\text{/}}T$.
(1.3)
$\begin{gathered} \xi (t)\, = \,2 \cdot \frac{a}{2}\sum\limits_{k = 0}^{ + \infty } {{{{( - 1)}}^{k}}} \cdot \frac{1}{{(2k\, + \,1)\pi }} \cdot cos(2k\, + \,1) \cdot \frac{{2\pi }}{T}t\, = \\ = \;\frac{a}{\pi }\left( {cos\Omega t - \frac{1}{3}cos3\Omega t + \; \ldots } \right), \\ \end{gathered} $
где $\Omega = 2\pi {\text{/}}T$.

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

Рис. 2.

Пример осцилляций Гиббса: (a) – исходная функция $\xi (t)$, (б) – функция, соответствующая обрезанному Фурье-спектру функции (a).

Важно заметить, что максимальная амплитуда осцилляций Гиббса постоянна и не зависит от выбранной частоты обрезки ряда [6].

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

Одним из наиболее эффективных подходов, основанных на применении необучаемых моделей, является работа [7], в которой ложное оконтуривание подавляется с помощью поиска оптимальных субпиксельных сдвигов. Этот метод был выбран нами в качестве классического алгоритма в предлагаемом гибридном подходе. Для краткости будем назвать его – алгоритм Кельнера. Предложенный метод авторы сравнивали с другими классическими подходами, например, с медианной фильтрацией, фильтрацией Ланцоша – и визуально было доказано превосходство описанного в статье [7] алгоритма Кельнера.

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

(1.4)
$J(u) = \frac{1}{2}\mathop {\left\| {u - {{u}_{0}}} \right\|}\nolimits^2 + \lambda \int\limits_\Omega {{\text{|}}\nabla } u(x){\text{|}}dx \to \mathop {min}\limits_{u \in U} ,$
где ${{u}_{0}}$ – входное изображение с артефактами ложного оконтуривания; $u$ – искомое оптимальное изображение среди изображений класса U; $\Omega $ – область задания изображения; $\lambda $ – параметр регуляризации (может зависеть от ближайшей контрастной границы, как например в [8]).

Также, известны методы одновременного подавления осцилляций Гиббса и оценки степени подавляемого ложного оконтуривания c использованием разреженных представлений [9].

В настоящее время все большую популярность приобретают методы глубокого обучения в области математических методов обработки изображений. В первую очередь – свёрточные нейронные сети. Свёрточные сети представляют возможность строить отображения изображений в обучаемые признаковые пространства большой размерности, фильтровать полученные признаки и преобразовывать обработанные признаковые представления в целевое цветовое пространство. Свёрточные нейронные сети успешно решают массу прикладных задач обработки биомедицинских изображений и компьютерного зрения, таких как сегментация биомедицинских изображений [10, 11], уменьшение уровня шума на медицинских изображениях [12], и, в том числе, подавление осцилляций Гиббса на изображениях магнитно-резонансной томографии [1315].

Одной из актуальных свёрточных моделей подавления артефактов ложного оконтуривания на изображениях МРТ, является GAS-CNN [14], которая была выбрана в качестве базовой. С результатами, полученными этой моделью, сравнивается предлагаемый гибридный подход.

GAS-CNN представляет собой очень глубокий свёрточный ансамбль (32 свёрточных блока), являющийся развитием модели увеличения разрешения EDSR [16]. Основными отличительными особенностями GAS-CNN стали:

• добавление внешних пространственных связей (англ. skip connectios), подобных используемым в архитектуре U-Net [10] между кодировщиком и декодировщиком;

• уменьшение размерности признакового пространства, в котором происходит фильтрация изображения (в GAS-CNN оно равно 64);

• отказ от использования слоев уменьшения пространственной размерности, например, свёртки с шагом 2 или слоя субдискретизации с функцией максимума (англ. max pooling), в силу локальности артефактов ложного оконтуривания.

GAS-CNN отображает входное изображение с осцилляциями Гиббса в признаковое пространство глубины 64, в котором затем осуществляется фильтрация с помощью 32-х остаточных свёрточных блоков (англ. ResBlock) [17]. Сеть преобразует ставшее результатом свёрточной фильтрации признаковое представление в изображение целевого цветового пространства – оттенки серого (англ. greyscale), с помощью заключительного реконструирующего свёрточного слоя.

Авторы GAS-CNN провели сравнение своей модели с билатеральной фильтрацией [18], нелокальной фильтрацией (англ. non-local means filtering) [19] и свёрточной нейронной сетью GARCNN [13].

Тем не менее, авторы GAS-CNN не исследовали возможную избыточность обучения 32-х свёрточных блоков для решения задачи уменьшения осцилляций Гиббса на изображениях МРТ. В данной статье предлагается новый, гибридный метод подавления артефактов ложного оконтуривания на изображениях магнитно-резонансной томографии, заключающийся в совместном использовании обучаемой глубокой свёрточной нейронной сети и необучаемого классического метода подавления осцилляций Гиббса. В качестве обучаемой части используется упрощение модели GAS-CNN, а в качестве необучаемой классической части выступает алгоритм Кельнера. Обучаемая часть предлагаемой модели состоит из 9 свёрточных блоков, вместо 32-ух у GAS-CNN, однако, благодаря использованию вспомогательных признаков от предварительно обработанных алгоритмом Кельнера входных изображений, наша гибридная модель сохраняет обобщающую способность на тестовом наборе данных. Такой гибридный подход, основывающийся на совместном использовании обучаемого и классического алгоритмов, ведет к улучшению итогового результата как визуально, так и по метрике оценки качества изображений PSNR [20]. Уменьшение числа свёрточных блоков также уменьшает риск переобучения модели.

2. АЛГОРИТМ КЕЛЬНЕРА

Алгоритм Кельнера [7] – это неитеративный необучаемый метод подавления осцилляций Гиббса на изображениях, основанный на поиске оптимальных субпиксельных сдвигов с целью минимизации функционала, оценивающего величину осцилляций. В качестве такого функционала используется полная вариация:

(2.1)
$TV(f,\Omega ) = \int\limits_\Omega {{\text{|}}\nabla } f(x){\text{|}}dx,$
где f – функция, полная вариация (TV) которой вычисляется в области $\Omega $.

Метод имеет два параметра ${{k}_{1}}$, ${{k}_{2}}$, определяющих “радиус” $[{{k}_{1}},{{k}_{2}}]$ пиксельного окна, в котором минимизируется полная вариация при различных субпиксельных сдвигах. Алгоритм Кельнера позволяет эффективно уменьшать эффект Гиббса с минимальным итоговым размытием, и его вычислительная сложность невелика.

Метод Кельнера осуществляет сначала поиск оптимальных субпиксельных сдвигов по двум одномерным направлениям: по горизонтали и по вертикали, а затем совмещает полученные результаты в Фурье области по следующей формуле:

(2.2)
$K = {\text{iFFT}}\left\{ {{\text{FFT}}({{K}_{x}}) \cdot {{G}_{x}} + {\text{FFT}}({{K}_{y}}) \cdot {{G}_{y}}} \right\},$
где K – итоговый результат алгоритма Кельнера; ${{K}_{x}}$ – результат одномерного алгоритма Кельнера вдоль оси OX; ${{K}_{y}}$ – результат одномерного алгоритма Кельнера вдоль оси OY; FFT – прямое дискретное преобразование Фурье; iFFT – обратное дискретное преобразование Фурье; ${{G}_{x}}$, ${{G}_{y}}$ – весовые функции, имеющие вид:
(2.3)
${{G}_{x}}({{t}_{x}},{{t}_{y}}) = \frac{{1 + cos{{t}_{y}}}}{{(1 + cos{{t}_{x}}) + (1 + cos{{t}_{y}})}},$
(2.4)
${{G}_{y}}({{t}_{x}},{{t}_{y}}) = \frac{{1 + cos{{t}_{x}}}}{{(1 + cos{{t}_{x}}) + (1 + cos{{t}_{y}})}},$
где ${{t}_{x}} \in [ - \pi ,\pi ]$, ${{t}_{y}} \in [ - \pi ,\pi ]$.

Рассмотрим подробнее одномерный алгоритм Кельнера с математической точки зрения. Пусть I – исходный сигнал и ${\text{\{ }}{{C}_{k}}{\text{\} }}$k-я комплексная амплитуда. Тогда, исходный сигнал можно приблизить следующим разложением в дискретный ряд Фурье:

(2.5)
$I(n) = \frac{1}{N}\sum\limits_{k = 0}^{N - 1} {{{C}_{k}}} {{e}^{{\frac{{2\pi ik}}{N}n}}}.$

Для каждого пикселя n с помощью теоремы о сдвиге преобразования Фурье получим 2M его значений в окрестности для дальнейшего выбора оптимального субпиксельного сдвига:

(2.6)
${{I}_{s}}(n) = \frac{1}{N}\sum\limits_{k = 0}^{N - 1} {\left( {{{C}_{k}} \cdot {{e}^{{ - \frac{{2\pi ik}}{N} \cdot \frac{s}{{2M}}}}}} \right)} \cdot {{e}^{{\frac{{2\pi ik}}{N}n}}},$
где $s \in \left[ { - M,M - 1} \right]$, $s \in \mathbb{Z}$. Авторы метода Кельнера предлагают произвести сначала поиск оптимального сдвига вдоль положительного и отрицательного направлений независимо:
(2.7)
$TV_{s}^{ + }(n) = \sum\limits_{k = {{k}_{1}}}^{{{k}_{2}}} {{\text{|}}{{I}_{s}}} (n + k) - {{I}_{s}}(n + (k - 1)){\text{|}},$
(2.8)
$TV_{s}^{ - }(n) = \sum\limits_{k = {{k}_{1}}}^{{{k}_{2}}} {{\text{|}}{{I}_{s}}} (n - k) - {{I}_{s}}(n - (k - 1)){\text{|}},$
(2.9)
${{r}^{ + }}(n) = {\text{arg}}\mathop {min}\limits_s TV_{s}^{ + }(n),$
(2.10)
${{r}^{ - }}(n) = {\text{arg}}\mathop {min}\limits_s TV_{s}^{ - }(n),$
где ${{k}_{1}}$, ${{k}_{2}}$ – параметры алгоритма, определяющие “радиус” пиксельной окрестности, по которой минимизируется полная вариация при различных субпиксельных сдвигах. Авторы пришли к выводу, что для большинства рассматриваемых ими прикладных случаев оптимальными входными параметрами являются: ${{k}_{1}} = 1$, ${{k}_{2}} = 3$.

А затем из двух полученных оптимальных сдвигов ${{r}^{ + }}(n)$, ${{r}^{ - }}(n)$ выбрать тот, что обеспечивает минимальную полную вариацию:

(2.11)
$r(n) = {\text{arg}}\mathop {min}\limits_{{{r}^{ + }},{{r}^{ - }}} \{ TV_{{{{r}^{ + }}(n)}}^{ + }(n),TV_{{{{r}^{ - }}(n)}}^{ - }(n)\} $

Окончательное значение в полученной точке $n - r(n){\text{/}}(2M)$ вычисляется авторами с помощью линейной интерполяции.

3. ПРЕДЛАГАЕМЫЙ МЕТОД

Предлагаемая архитектура представлена на рис. 3.

Рис. 3.

Схема предлагаемой гибридной архитектуры подавления осцилляций Гиббса – DGAS9-CNN.

Будем называть ее DGAS9-CNN (англ. Dual Gibbs-ringing artifact suppression CNN).

Предлагаемая гибридная модель: $C \to {{I}_{{out}}}$, имеет следующую структуру (номера соответствуют номерам на рис. 3):

1. Вход: обрезанный спектр Фурье (C).

2. Получение изображения с артефактами Гиббса, являющегося входным в свёрточную нейронную сеть: $C \to {{I}_{{{{G}_{l}}}}}$:

• дополнение нулями по периметру (англ. zero-padding) C, для соответствия размерам референсного изображения без осцилляций Гиббса ${{I}_{{gt}}}$;

• применение обратного преобразования Фурье.

3. Получение изображения с осцилляциями Гиббса, являющегося входным в алгоритм Кельнера: $C \to {{I}_{{{{G}_{c}}}}}$:

• применение обратного преобразования Фурье к входному спектру C без дополнения его нулями по периметру в силу требований классического метода.

4. Применение алгоритма Кельнера к ${{I}_{{{{G}_{c}}}}}$: ${{I}_{{{{G}_{c}}}}} \to {{I}_{K}}$.

5. Повышение разрешения IK  до соответствия с размерами референсного изображения: ${{I}_{K}}\, \to \,{{I}_{{K \uparrow }}}$.

6. Модуль преобразования изображения ${{I}_{{{{G}_{l}}}}}$ в признаковое представление глубины 64.

7. Модуль преобразования изображения ${{I}_{{K \uparrow }}}$ в признаковое представление глубины 64.

8. Модуль глубокого кодирования признаков изображения ${{I}_{{{{G}_{l}}}}}$ в признаковом пространстве глубины 64.

9. Модуль глубокого кодирования признаков изображения ${{I}_{{K \uparrow }}}$ в признаковом пространстве глубины 64.

10. Модуль восстановления признаков целевого изображения ${{I}_{{out}}}$, в котором агрерируются признаки изображений ${{I}_{{{{G}_{l}}}}}$ и ${{I}_{{K \uparrow }}}$.

11. Модуль реконструкции целевого изображения без эффекта Гиббса ${{I}_{{out}}}$.

Предлагаемая гибридная модель DGAS9-CNN имеет 9 структурных остаточных блоков: 3 блока в глубоком кодировщике признаков изображения ${{I}_{{{{G}_{l}}}}}$, 3 блока в глубоком кодировщике признаков изображения ${{I}_{{K \uparrow }}}$ и 3 блока в модуле восстановления признаков целевого изображения ${{I}_{{out}}}$. Схема остаточного блока предлагаемой модели приведена на рис. 4 $(a)$.

Рис. 4.

(а) – остаточный блок; (б) – блок агрегации признаков. Conv – свёрточный слой; BN – пакетная нормализация; ReLU – функция активации ReLU; Concat – слой конкатенации вдоль последней размерности (форма тензора описывается по стандарту: BHWC).

В [17] показано, что использование остаточных блоков в качестве структурных архитектурных составляющих позволяет преодолеть одну из основных проблем обучения глубоких свёрточных нейронных сетей – затухание градиента [21]. Стоит отметить, что выбранная базовая модель GAS-CNN также состоит из остаточных блоков, однако в отличие от GAS-CNN используемый нами остаточный блок включает в себя слои пакетной нормализации [22] (англ. batch normalization), позволяющие еще более эффективно вести обучение глубоких архитектур [23].

Особенностью модуля восстановления признаков целевого изображения ${{I}_{{out}}}$ является блок агрегации признаков изображений ${{I}_{{{{G}_{l}}}}}$ и ${{I}_{{K \uparrow }}}$. Его схема приведена на рис. 4 $(b)$. Он представляет свёртку с ядром 3 × 3 предварительно конкатенированных и нормализованных представлений изображений ${{I}_{{{{G}_{l}}}}}$ и ${{I}_{{K \uparrow }}}$.

Извлечение признаковых представлений изображений ${{I}_{{{{G}_{l}}}}}$ и ${{I}_{{K \uparrow }}}$ предлагаемая архитектура реализует с помощью одного свёрточного слоя (3 × 3) с последующей пакетной нормализацией и активацией ReLU.

Повышение разрешения результата алгоритма Кельнера IK реализуется интерполяцией методом ближайшего соседа с последующими двумя обучаемыми свёрточными слоями с ядрами (3 × 3).

Заключительный модуль реконструкции реализован через два последовательных свёрточных слоя с ядрами (3 × 3), осуществляющих:

• промежуточное проецирование обработанных признаков глубины 64 на подпространство меньшей размерности – 16;

• отображение полученного тензора в целевое изображение без эффекта Гиббса ${{I}_{{out}}}$.

4. ЭКСПЕРИМЕНТЫ

Предлагаемая гибридная модель DGAS9-CNN была обучена и протестирована на синтезированных изображениях с фиксированной степенью осцилляций Гиббса, полученных прямым и обратным преобразованиями Фурье из набора данных IXI11 референсных изображений ${{I}_{{gt}}}$ по следующей схеме [14]:

1. Применить преобразование Фурье к очередному референсному изображению ${{I}_{{gt}}}$ (256 × 256) из набора данных IXI.

2. Обрезать спектр Фурье таким образом, чтобы сохранилась лишь центральная часть частот, составляющая $\tfrac{1}{9}$ площади (см. изображение C на рис. 3).

3. Получить изображение (${{I}_{{{{G}_{l}}}}}$) с осцилляциями Гиббса, являющееся входным в свёрточную нейронную сеть:

• дополнить нулями по периметру спектр C до соответствия с размерами референсного изображения ${{I}_{{gt}}}$;

• примененить обратное преобразование Фурье.

4. Получить изображение (${{I}_{{{{G}_{c}}}}}$), являющееся входным для алгоритма Кельнера:

• применить обратное преобразование Фурье к обрезанному спектру $C$.

В результате тренировочный, валидационный и тестовый наборы данных представляют собой множество триплетов: (${{I}_{{gt}}}$, ${{I}_{{{{G}_{l}}}}}$, ${{I}_{K}}$). Тренировочный сет состоит из 10 427 триплетов, валидационный – из 2016 триплетов и тестовый – из 2617 триплетов [14].

Следует отметить, что набор данных IXI содержит МРТ-данные трех различных модальностей: T1, T2, PD, поэтому для приведения значений пикселей изображений к единому диапазону [0, 1] была использована стандартная Min-Max нормализация каждого триплета (${{I}_{{gt}}}$, ${{I}_{{{{G}_{l}}}}}$, ${{I}_{K}}$).

Авторы базовой статьи GAS-CNN [14] не предоставили кода и весов, так что для корректного сравнения мы обучили все сравниваемые здесь модели самостоятельно. Программная реализация сравниваемых здесь архитектур написана на языке Python 3 с использованием фреймворка глубокого обучения Tensorflow 1.1422.

Модели обучались с помощью метода оптимизации Adam [24] с параметрами: ${{\beta }_{1}} = 0.9$, ${{\beta }_{2}}$ = = 0.999, $\varepsilon = 1{\text{e - }}08$ на графическом ускорителе GPU NVIDIA GeForce RTX 2080 Ti. Сравниваемые здесь модели оптимизировались с точки зрения классической регрессионной, устойчивой к выбросам, функции потерь ${{l}_{1}}$ с ${{l}_{2}}$ регуляризацией ($\gamma = {{10}^{{ - 4}}}$). Скорость обучения менялась полиномиально:

(4.1)
$lr(x) = (l{{r}_{0}} - l{{r}_{1}}) \cdot {{\left( {1 - \frac{x}{M}} \right)}^{p}} + l{{r}_{1}},$
где $x$ – очередная итерация обучения; $l{{r}_{0}}$ – исходное значение скорости обучения; $l{{r}_{1}}$ – конечное значение скорости обучения; $M = {{N}_{e}} \cdot {{N}_{d}}{\text{/}}bs$ – общее количество итераций обучения; ${{N}_{e}}$ – число эпох обучения; ${{N}_{d}}$ – число триплетов в тренировочном наборе данных; bs – пакет триплетов, подаваемых модели на очередной итерации обучения (англ. batch size); $p$ – степень полинома, с которой ведется полиномиальное уменьшения скорости обучения.

Непосредственно, были использованы следующие значения представленных выше параметров: $l{{r}_{0}} = {{10}^{{ - 4}}}$, $l{{r}_{1}} = 0$, ${{N}_{e}} = 1000$, $bs = 20$, $p = 0.3$.

Важно отметить, что для борьбы с потенциальной проблемой переобучения применялась не только ${{l}_{2}}$ регуляризация, но и аугментация входных данных: случайными вращениями на углы, кратные 90, и отражениями по горизонтали и по вертикали. Модели обучались на случайно вырезаемых в процессе обучения патчах размера 48 × 48.

Отдельное внимание было уделено возможной избыточности обучения 32-х свёрточных блоков для решения задачи подавления осцилляций Гиббса на изображениях магнитно-резонансной томографии.

Для этого, сначала была обучена сильно упрощенная версия GAS-CNN с 6-ю свёрточными блоками (помимо модулей выделения признаковых представлений входного изображения и модуля реконструкции целевого изображения – эти модули сохраняют свою архитектуру для всех исследуемых здесь моделей) – GAS6-CNN, повторяющая остаточную архитектуру GAS-CNN, однако имеющая 6 свёрточных остаточных блоков вместо исходных 32-х. Полученные средние значения метрики оценки качества изображений PSNR на тестовом наборе данных для GAS-CNN и GAS6-CNN отражены в табл. 1, а динамика средних значений PSNR на валидационном наборе данных во время обучения показана на рис. 5. Видно, что такое значительное уменьшение числа слоев модели при равных условиях обучения привело к потере обобщающей способности на тестовом наборе данных.

Таблица 1.

Средние значения PSNR на тестовом наборе данных из 2617 изображений и средние времена обработки одного изображения на ЦПУ

Модель PSNR Время, с
Исходные изображения 29.79
GAS-CNN 32.04 1.03
GAS6-CNN 31.67 0.25
DGAS9-CNN 32.27 0.68
DGAS9-CNN 31.61 0.45
Рис. 5.

Средние значения PSNR на валидационном наборе из 2016 изображений во время обучения.

Затем, добавив к GAS6-CNN дополнительные признаки, кодирующие обработанные классическим алгоритмом Кельнера входные изображения, уровень обобщающей способности модели на тестовом наборе данных повысился и сравнился с результатом работы GAS-CNN, имеющей почти в 3.5 раза больше свёрточных остаточных блоков. (см. табл. 1 и рис. 5).

Для оценки вклада алгоритма Кельнера в гибридный метод, мы дополнительно обучили модель $\overline {{\text{DGAS}}9{\text{ - CNN}}} $, полностью повторяющую структуру предлагаемой DGAS9-CNN, однако без использования алгоритма Кельнера. В ней изображение IK совпадает с ${{I}_{{{{G}_{c}}}}}$. Полученный результат, как можно видеть из табл. 1 и рис. 5, уступает DGAS9-CNN и на тестовом наборе данных, и на валидационной выборке.

На центральном процессоре Intel(R) Core(TM) i7-8700 базовый подход GAS-CNN обрабатывает одно изображение 256 × 256 за 1.03 сек, в то время как предлагаемая нами гибридная модель DGAS9-CNN – в полтора раза быстрее, за 0.68 с (см. табл. 1), обеспечивая сравнимую с GAS-CNN обобщающую способность на тестовом наборе данных (см. табл. 1).

Визуально оценить результаты применения рассматриваемых алгоритмов можно на рис. 6. Результат применения алгоритма Кельнера имеет наименьшую полную вариацию среди рассматриваемых алгоритмов, однако, минимизируя полную вариацию, метод также подавляет и детали на изображениях, что приводит к приросту только в 1 дБ метрики PSNR (см. табл. 2). Результаты применения предлагаемого гибридного алгоритма DGAS9-CNN и базового GAS-CNN имеют также меньшую полную вариацию в сравнении с исходными входными изображениями, однако они значительно лучше восстанавливают детали на изображениях, увеличивая PSNR: на 8.87 и 8.76 дБ, соответственно (см. табл. 2).

Рис. 6.

Визуальное сравнение алгоритмов.

Таблица 2.

Средние значения PSNR и полной вариации на тестовом наборе данных из 2617 изображений. Изображения приведены к единому размеру ${{I}_{K}}$

Модель PSNR TV
Исходные изображения 20.70 659.68
Алгоритм Кельнера 21.68 540.18
DGAS9-CNN 29.57 611.25
GAS-CNN 29.46 620.09

5. ЗАКЛЮЧЕНИЕ

В статье представлена гибридная модель DGAS9-CNN подавления осцилляций Гиббса на изображениях магнитно-резонансной томографии. Предлагаемый метод состоит из агрегации классического неитеративного подхода подавления осцилляций Гиббса – алгоритма Кельнера, и легкой свёрточной нейронной сети, упрощающей архитектуру сети GAS-CNN. Предлагаемая гибридная модель имеет в 3.5 раза меньше свёрточных остаточных блоков, в сравнении с GAS-CNN, однако позволяет достичь лучшие, чем GAS-CNN результаты, обрабатывая изображения в 1.5 раза быстрее.

Список литературы

  1. Senyukova O., Zubov A. Full anatomical labeling of magnetic resonance images of human brain by registration with multiple atlases // Programming and Computer Software. 2016. V. 46. № 6. P. 356–360.

  2. Gottlieb D., Orszag S. Numerical Analysis of Spectral Methods: Theory and Application. SIAM, 1977. 176 p.

  3. Gray A., Pinsky M. Gibbs phenomenon for Fourier-Bessel series // Expositiones Mathematicae. 1993. V. 11. 123 p.

  4. Pinsky M.A. Fourier inversion for piecewise smooth functions in several variables // Proceedings of the American Mathematical Society. 1993. V. 118. № 3. P. 903–910.

  5. Pinsky M.A. Pointwise Fourier inversion in several variables // Notices of the American Mathematical Society. 1995. V. 42. № 3. P. 330–334.

  6. Малла С. Вэйвлеты в обработке сигналов. Пер. с англ. М.: Мир, 2005. 671 с.

  7. Kellner E., Dhital B., Kiselev V.G., Reisert M. Gibbs?ringing artifact removal based on local subvoxel?shifts // Magnetic Resonance in Medicine. 2016. V. 76. № 5. P. 1574–1581.

  8. Sitdikov I.T., Krylov A.S. Variational Image Deringing Using Varying Regularization Parameter // Pattern Recognition and Image Analysis: Advances in Mathematical Theory and Applications. 2015. V. 25. № 1. P. 96–100.

  9. Umnov A.V., Krylov A.S. Sparse Approach to Image Ringing Detection and Suppression // Pattern Recognition and Image Analysis: Advances in Mathematical Theory and Applications. 2017. V. 27. № 4. P. 754–762.

  10. Ronneberger O., Fischer P., Brox T. U-net: Convolutional networks for biomedical image segmentation // International Conference MICCAI 2015. 2015. P. 234–241.

  11. Sinha A., Dolz J. Multi-scale self-guided attention for medical image segmentation // IEEE Journal of Biomedical and Health Informatics. 2021. V. 25. № 1. P. 121–130.

  12. Krylov A., Karnaukhov V., Mamaev N., Khvostikov A. Hybrid Method for Biomedical Image Denoising // Proceedings of the 2019 4th International Conference on Biomedical Imaging, Signal Processing. 2019. P. 60–64.

  13. Wang Y., Song Y., Xie H. et al. Reduction of Gibbs artifacts in magnetic resonance imaging based on Convolutional Neural Network // 2017 10th International Congress on Image and Signal Processing, Biomedical Engineering and Informatics (CISP-BMEI). 2017. P. 1–5.

  14. Zhao X., Zhang H., Zhou Y. et al. Gibbs-ringing artifact suppression with knowledge transfer from natural images to MR images // Multimedia Tools and Applications. 2019. P. 1–23.

  15. Penkin M., Krylov A., Khvostikov A. Attention-based Convolutional Neural Network for MRI Gibbs-ringing Artifact Suppression // CEUR Workshop Proceedings. 2020. V. 2744. P. 1–12.

  16. Lim B., Son S., Kim H. et al. Enhanced deep residual networks for single image super-resolution // Proceedings of the CVPR IEEE Conference. 2017. P. 136–144.

  17. He K., Zhang X., Ren S., Sun J. Deep residual learning for image recognition // Proceedings of the CVPR IEEE Conference. 2016. P. 770–778.

  18. Zhang M., Gunturk B.K. Multiresolution bilateral filtering for image denoising // IEEE Transactions on Image Processing. 2008. V. 17. № 12. P. 2324–2333.

  19. Manj’on J.V., Coupé P., Buades A. et al. Non-local MRI upsampling // Medical Image Analysis. 2010. V. 14. № 6. P. 784–792.

  20. Hore A., Ziou D. Image quality metrics: PSNR vs. SSIM // 2010 20th International Conference on Pattern Recognition. 2010. P. 2366–2369.

  21. Hanin B. Which neural net architectures give rise to exploding and vanishing gradients? // Advances in Neural Information Processing Systems. 2018. P. 582–591.

  22. Ioffe S., Szegedy C. Batch normalization: Accelerating deep network training by reducing internal covariate shift // arXiv preprint arXiv:1502.03167, 2015.

  23. Santurkar S., Tsipras D., Ilyas A., Madry A. How does batch normalization help optimization? // Advances in Neural Information Processing Systems. 2018. P. 2483–2493.

  24. Kingma D.P., Ba J. Adam: A method for stochastic optimization // arXiv preprint arXiv:1412.6980, 2014.

Дополнительные материалы отсутствуют.