Содержание к диссертации
Введение
1 Актуальность задачи массовой обработки медицинских рентгенографических снимков. Обзор существующих средств обработки графической медицинской информации 10
1.1 Задача массовой обработки рентгенографических медицинских изображений и ее актуальность 10
1.2 Обзор существующих программных средств обработки медицинских изображений 12
1.3 Анализ известных методов обработки цифровых медицинских изображений 23
1.4 Формальная постановка задачи исследования 27
1.5 Выводы 28
2 Методы детектирования объектов и уточнения их границ на рентгенографических медицинских изображениях 30
2.1 Метод детектирования объектов на рентгенографических медицинских изображениях 33
2.2 Метод уточнения границ объектов на рентгенографических медицинских изображениях 59
2.3 Осуществление оценки расстояний между ключевыми точками объектов 68
2.4 Выводы 74
3 Программный комплекс для анализа и интерактивной диагностики медицинских рентгенографических снимков 76
3.1 Библиотека обработки изображений FIPL 76
3.2 Программный комплекс XRIDA 95
3.3 Выводы 123
4 Опыт практического применения программного комплекса XRIDA и библиотеки обработки изображений FIPL 125
4.1 Интерактивная диагностика по единичному рентгенологическому исследованию коленного сустава 127
4.2 Сравнительный анализ возрастных изменений коленного сустава у детей 4–6 лет 132
4.3 Массовая обработка рентгенографических изображений коленного сустава с применением механизма фильтрации по диагностическим критериям 137
4.4 Выводы 144
Заключение 146
Список рисунков 150
Список таблиц 151
Термины и обозначения 152
Литература 152
- Обзор существующих программных средств обработки медицинских изображений
- Метод уточнения границ объектов на рентгенографических медицинских изображениях
- Программный комплекс XRIDA
- Сравнительный анализ возрастных изменений коленного сустава у детей 4–6 лет
Обзор существующих программных средств обработки медицинских изображений
Вместе с тем, в силу относительной простоты получаемых изображений, а также благодаря тому, что основные методы рентгеновской диагностики были разработаны довольно давно, зачастую использование компьютерных средств ограничивается визуализацией изображения на экране без каких-либо автоматизированных средств их обработки или анализа. Обзор программных продуктов медицинской визуализации показывает, что способы работы с рентгенографическими изображениями ограничены в лучшем случае инструментами ручных количественных измерений и не предполагают применения современных средств обработки изображений. Кроме того, зачастую, во многих лечебных учреждениях в России используется устаревшее материально-техническое обеспечение, что приводит к тому, что диагностика осуществляется «обычными» методами, при помощи пленочных аппаратов.
Роль рентгенографической диагностики, однако, остаётся существенной, что, в частности, обусловлено повсеместным распространением рентгеновских аппаратов. Рентгенография также обладает некоторыми преимуществами по сравнению с более «тяжеловесными» средствами диагностики и зачастую является наиболее простым и надежным средством.
В настоящее время в лечебных учреждениях России на регулярной основе проводится ряд исследований с использованием проективной рентгенографии. При этом массовый поток изображений, подлежащих анализу, обрабатывается специалистами вручную, что занимает существенное время и создает высокую нагрузку на рентгеновские кабинеты.
С учетом сказанного, представляется весьма актуальной разработка методов и программных средств массовой интерактивной обработки цифровых рентгеновских снимков.
Объектом исследования являются математические и программные средства цифровой обработки медицинских изображений. Целью настоящей работы является повышение эффективности диагностического процесса при диагностике и профилактике различных заболеваний у населения.
Общая научная задача исследования состоит в разработке методов и средств автоматизированной диагностики по медицинским рентгенограммам, обеспечивающих уменьшение количества рентгенограмм, обрабатываемых вручную при заданных исходных данных о диагностических критериях и допустимой доле ложно-отрицательных результатов.
Научная новизна
Наиболее существенные положения, выдвигаемые для защиты
1. Существующие математические методы обработки цифровых изображений не учитывают особенности медицинских цифровых рентгенограмм и не дают приемлемого для диагностики качества детектирования и выделения границ. Дальнейшее развитие таких методов в целях повышения качества обработки имеет существенное практическое значение для повышения эффективности диагностического процесса.
2. Использование предлагаемых в работе методов обнаружения объектов, уточнения их границ и автоматизированной проверки диагностических критериев при их программной реализации позволяет существенно (в эксперименте c конкретным критерием — до 90%) уменьшить количество рентгенограмм, требующих внимания специалиста и ручной обработки.
Наиболее существенные новые научные результаты, выдвигаемые для защиты:
1. Метод детектирования объектов на медицинских рентгенограммах, отличающийся от известных учётом точного направления вектора градиента интенсивности изображения при минимизации расстояния от шаблона искомого объекта до границ, обнаруженных на изображении. 2. Новый метод фильтрации градиентного препарата изображения, использующий сведения о приблизительном местоположении целевого объекта, полученные на этапе детектирования.
3. Метод уточнения границ объекта на базе известного метода активных контуров, отличающийся использованием сведений о приблизительном местоположении целевого объекта для подавления границ посторонних объектов.
4. Модель диагностического критерия, отличающаяся формализованным описанием диагностического критерия как предиката над множеством относительных расстояний между ключевыми точками на рентгенограмме.
5. Программный комплекс медицинской визуализации, реализующий предложенные методы, отличающийся наличием процедур автоматического контроля качества вводимого шаблона и его корректировки.
Апробация работы. Основные положения работы прошли обсуждение на научных семинарах механико-математического факультета (июль 2013, март 2014) и НИИ многопроцессорных вычислительных систем (июль 2013) Южного федерального университета. Результаты работы представлялись на следующих конференциях: Научно-методическая конференция «Современные информационные технологии в образовании: Южный федеральный округ» (Южный федеральный университет, г. Ростов-на-Дону, 15–16 апреля 2010 г.), X Всероссийская научная конференция молодых ученых, аспирантов и студентов «Информационные технологии, системный анализ и управление» (Южный федеральный университет, г. Таганрог, январь 2013 г.), XX научная конференция «Современные информационные технологии: тенденции и перспективы развития» (Южный федеральный университет, г. Ростов-на-Дону, 24–26 апреля 2013 г.), 9th International Conference on Computer Science and Information Technologies (CSIT-2013), (НИИ информатики и проблем автоматизации НАН Армении, г. Ереван, Армения, 23–27 сентября 2013 г.),
Метод уточнения границ объектов на рентгенографических медицинских изображениях
Важной информативной составляющей рентгенографического изображения являются границы объектов. Под границей понимается линия, отделяющая данный объект от фона или от других объектов. В основе всех методов, предлагаемых в настоящем исследовании, лежит естественное предположение, что в точках изображения, принадлежащих границе объекта, имеет место перепад интенсивности. Обратное не всегда верно, и точки, в которых интенсивность меняет свое значение, не всегда принадлежат к границам целевого объекта.
Таким образом, первым этапом обработки изображения является выявление точек существенного перепада интенсивности, как наиболее вероятных кандидатов на принадлежность к границам объектов. Для предварительной подготовки изображения и обнаружения перепадов интенсивности применяется ряд известных стандартных методов обработки изображений.
.Здесь «о» означает операцию свертки, а 0IK — нормированная матрица гаус-совского фильтра со среднеквадратичным отклонением и размерами х . Для вычисления свертки на границах изображения исходное изображение непрерывно продолжается за границы, т. е.
Такой подход позволяет упростить формулы и избежать особой обработки при вычислениях свертки на границе изображения.
Применение гауссовского фильтра сглаживает перепады интенсивности и частично устраняет шумы, присутствующие на изображении.
Следующий этап составляет собственно детектирование перепадов. В обработке изображений широкое распространение получили дифференциальные методы детектирования перепадов, состоящие в численной аппроксимации значений градиента интенсивности.
В настоящей работе для детектирования перепадов использовались широко известные операторы Собеля [59, C. 112-114] и Шарра [60, C. 144-150]. Данные операторы применяются к размытому изображению : соответственно горизонтальная и вертикальная матрицы используемого дифференциального оператора. Можно воспринимать как матрицу, элементами которой являются векторы из R2, либо как две матрицы, со 37 держащие горизонтальную и вертикальную компоненты вычисленного вектора градиента в каждой точке. Полученную матрицу будем называть градиентным препаратом изображения
Применение дифференциальных операторов к изображению позволяет вычислить приближенное значение производной изображения. Однако в силу природы дифференциальных операторов, значения модуля градиента изображения будут возрастать не только в местах максимальных перепадов интенсивности, но и в некоторой окрестности таких перепадов, что будет порождать «толстую» линию границы (см. рис. 2.2, b).
Детектирование перепадов. Исходное изображение (a), модуль градиента изображения, вычисленного при помощи оператора Собеля (b), результат работы детектора Кенни (c). С точки зрения обработки изображений, наибольший интерес представляют как раз места максимальных перепадов интенсивности. Для выделения этих участков и получения «тонких» линий применяется детектор Кенни [57].
Детектор Кенни подавляет немаксимальные значения перепада интенсивности по данному направлению перепада. Результат применения фильтра Кенни в данной точке совпадает со значением исходной матрицы, если это значение является локальным максимумом по соответствующему направлению, и равно нулю в противном случае. В литературе детектор (или фильтр) Кенни иногда рассматривается как самостоятельный алгоритм выделения границ, входными данными для которого являются полутоновые изображения. В этом случае в понятие детектора Кенни включается предварительная фильтрация и вычисление градиента изображения при помощи какого-либо дифференциального оператора. В настоящей работе под детектором Кенни будем понимать только фильтр, применяющийся к градиенту изображения с целью подавления локально немаксимальных значений перепадов интенсивности.
Детектирование целевых объектов на изображении предполагает наличие некоторых априорных сведений об искомом объекте в том или ином виде. В настоящей работе в качестве носителя таких сведений предлагается шаблон, содержащий сведения о форме границы искомого объекта, а также некоторые дополнительные данные.
Сведения о границе объекта хранятся в виде контура, то есть последовательности точек плоскости. Любые точки, рассматриваемые на цифровом изображении, могут иметь лишь целочисленные координаты. Однако в целях избежания накопления ошибок округления при последовательном применении к контуру различных трансформаций будем рассматривать контур как последовательность точек из Ш2.
В процессе работы с контурами регулярно возникает необходимость вывода (отображения) контура на изображении, а также оценка каких-либо параметров изображений в точках Z2, ближайших к точкам контура. В связи с этим представляется целесообразным ввести функцию, действующую из М2 и ставящую каждому вектору в соответствие ближайший целочисленный индекс из Z2. Рассмотрим такую функцию [] : М2 — Z2:
Контур С, образующий шаблон, представляет собой приближение границы целевого объекта ломаной, то есть его точки должны принадлежать этой границе, причем так, чтобы ребра контура также повторяли ход границы и покрывали всю интересующую часть границы объекта. Контур шаблона считается замкнутым и не может иметь самопересечений.
Последовательность чисел F будем называть полем направлений шаблона. Каждый элемент F отвечает одной точке шаблона и содержит приближенную ориентацию вектора градиента интенсивности в соответствующей точке изображения. Ориентация измеряется как величина угла между положительным направлением горизонтальной оси координат и направлением данного вектора, выраженная в радианах и приведенная к полуинтервалу [0,2-7г). Подробнее о получении значений этих направлений сказано ниже.
Предполагается, что один шаблон служит для обнаружения конкретного класса объектов на изображении. Из анатомии известно, что форма элементов скелета человека в состоянии, соответствующем медицинской и биологической норме, определяется тремя параметрами: возрастом человека, полом и типом конституции (астенический, нормостенический, гиперстенический). Применительно к детектированию объектов, это означает, что шаблон, вообще говоря, должен строиться отдельно для каждой возрастной, половой и конституционной группы
Программный комплекс XRIDA
В настоящее время библиотека OpenCV является де-факто стандартом для построения приложений и проведения исследований, посвященных обработке изображений. Несмотря на это, по целому ряду причин было принято решение построить собственную библиотеку обработки изображений, целиком написанную для платформы .NET. Во-первых, предлагаемые в настоящей работе методы используют весьма ограниченный набор общеизвестных алгоритмов, большая часть которых не требует высокой производительности и не трудоёмка в реализации. Таким образом, в разрабатываемом приложении не требуется вся мощь библиотеки OpenCV и её использование оказалось бы избыточным. Во-вторых, применение OpenCV обязало бы к использованию языка C++ в качестве основного языка разработки. Это связано с тем, что программные средства, обеспечивающие возможность использования OpenCV в приложениях .NET, разрабатываются на добровольных началах независимыми разработчиками и не поддерживаются в актуальном состоянии. При этом разработка под C++ гораздо более трудоёмка по сравнению с .NET. Использование OpenCV также обязало бы распространять вместе с приложением различные версии этой библиотеки, скомпилированные для разных архитектур ПК. Это уменьшило бы переносимость программного комплекса, увеличило бы объём дистрибутива и накладные расходы на развертывание приложения. Кроме того, так как программный комплекс разрабатывался прежде всего для использования в медицинских учреждениях, особое внимание уделялось построению удобного пользовательского интерфейса, а в этом плане платформа .NET предоставляет обширнейшие возможности при относительной простоте разработки.
Библиотека обработки изображений FIPL написана на языке C# под платформу .NET. Библиотека была разработана с целью упростить разработку приложений, связанных с цифровой обработкой изображений, для платформы .NET. Библиотека написана в объектно-ориентированном стиле и содержит классы, соответствующие основным понятиям цифровой обработки изображений: изображения, области, индексы, преобразования, фильтры и др. FIPL полностью разработана автором диссертации [68] (работа опубликована в соавторстве, автором диссертации предложена идея библиотеки и ее архитектура, а также выполнена реализация библиотеки) и является проектом с открытым исходным кодом [69], распространяемым по лицензии GNU GPL v.3 [70].
Библиотека сконструирована таким образом, чтобы позволять создавать относительно легко читаемый код. При этом была произведена попытка создания такой структуры классов и их программных интерфейсов, которая бы максимально соответствовала терминологии, обозначениям и понятийному аппарату, принятым в обработке изображений. Благодаря этому при программировании становится возможным больше внимания уделять логике алгоритма, нежели тонкостям программной реализации.
Классы библиотеки организованы в несколько пространств имен: FIPL — корневое пространство имен. Содержит все основные классы для работы с изображениями, их представлениями, преобразованиями и др. FIPL.Calculators — содержит классы-калькуляторы, то есть классы, инкапсулирующие арифметические операции над элементами изображений различных типов. FIPL.Exceptions — содержит специфичные для FIPL исключения, генерируемые методами библиотеки. FIPL.Drawing — классы, инкапсулирующие методы рисования на изображениях. FIPL.Mathematics — классы, содержащие реализации математических функций, необходимых для работы алгоритмов библиотеки, но отсутствующие в стандартной поставке .NET. Также инкапсулирует матричную алгебру и модели вероятностных распределений. FIPL.Geometry — методы геометрии на плоскости. FIPL.Animation — классы, связанные с работой с анимацией. FIPL.Windows — классы, инкапсулирующие ввод и вывод изображений при помощи пространства имен System.Drawing, а также элементы пользовательского интерфейса для представления изображений в Windows Forms. FIPL.Visualizers — содержит визуализаторы изображений для отладчика Visual Studio. FIPL.GTK — классы, инкапсулирующие ввод и вывод изображений при помощи Gtk#, а также элементы пользовательского интерфейса для представления изображений в Gtk#. FIPL.Tools — вспомогательные классы для работы с объектами библиотеки. FIPL.Tests — unit-тесты классов библиотеки. Граф зависимостей между пространствами имен библиотеки представлен на рис. 3.2. Рисунок 3.2: Граф зависимостей между пространствами имен библиотеки FIPL Библиотека строится в соответствии с парадигмой SOLID [72]. В связи с этим везде, где только возможно, применяется принцип инверсии зависимостей [73], что предполагает наличие разветвленной структуры интерфейсов.
Принцип разделения интерфейса [74], входящий в парадигму SOLID, предписывает создание большого количества мелких специфичных интерфейсов вместо небольшого количества универсальных.
Использование интерфейсов также позволяет реализовывать шаблон проектирования «примесь» (mixin) [75], который с появлением в языке C# методов-расширений становится одним из наиболее используемых шаблонов в C# [76]. Наконец, использование интерфейсов вкупе с методами расширения [77, 460–470] позволяет в какой-то мере компенсировать отсутствие в языке C# множественного наследования классов (напомним, что язык C# допускает множественное наследование интерфейсов).
Рассмотрим основные интерфейсы, которые используются для формирования структуры классов. Интерфейс ISize описывает сущность, имеющую два линейных размера. Он экспортирует два свойства, возвращающие длину и ширину объекта в целочисленных значениях. С ним тесно связан маркерный [78, cc. 179–180] интерфейс IBox, наследующий ISize и служащий для указания на объект, расположенный в начале координат. Различие между этими двумя интерфейсами станет очевидным, когда мы рассмотрим методы-расширения, связанные с каждым из них.
От интерфейса IBox наследуется интерфейс IMatrix, который также является маркерным и, соответственно, не привносит новых членов. IMatrix определяет сущности, являющиеся «матрицами», т. е. имеющими целочисленные размеры и содержащие в себе элементы, которые могут быть прочитаны.
Маркерный интерфейс ICanvas также наследует интерфейс IBox и аналогичен IMatrix, за тем исключением позволяет реализациям осуществлять запись элементов. Ключевыми и наиболее часто используемыми интерфейсами являются обобщенные интерфейсы IMatrix T и ICanvas T . Именно на уровне абстракции этих интерфейсов работают большинство методов библиотеки.
IMatrix T наследует интерфейсы IMatrix, IGetIndexed Index,T , IGetIndexed Int32,Int32,T . Этот интерфейс абстрагирует объекты, имеющие размеры, и индексируемые (на чтение) двумя целочисленными координатами или одним индексом из Z2 (причем элементы имеют тип T).
ICanvas T наследует интерфейсы ICanvas, ISetIndexed Index,T , ISetIndexed Int32,Int32,T . Этот интерфейс абстрагирует объекты, имеющие размеры и индексируемые (на запись) двумя целочисленными координатами или одним индексом из Z2 (причем элементы имеют тип T).
Два предыдущих интерфейса объединены в интерфейсе IImage T , который, таким образом, абстрагирует объекты имеющие размеры, а также элементы типа T, которые можно читать и изменять.
Иерархия базовых интерфейсов приведена на рис 3.3. Также важным интерфейсом является IPaddingMatrix T , который наследует IMatrix и представляет собой продолжение изображения. Этот интерфейс определяет способ продолжения матрицы (непрерывное, дополнение константой и др.), а также константу типа T для продолжения.
Структура данных Index представляет собой индекс из Z2. Одним из ключевых классов в библиотеке является класс Image T , который реализует интерфейсы IImage T и ICloneable [77, С. 349]. Этот класс хранит элементы изображения в массиве в памяти. Класс Image T наиболее часто используется для представления изображений, загруженных с жесткого диска, а также для возврата результатов преобразований. Данный класс является сериализуемым [77, CC. 804–813], что позволяет сохранять его экземпляры на жестком диске или в бинарных полях базы данных.
Сравнительный анализ возрастных изменений коленного сустава у детей 4–6 лет
Осуществление собственно детектирования объектов производится классом Detector. Этот класс осуществляет вычисления расстояния между контуром, трансформированным в соответствии с элементом пространства гипотез, и картой границ изображения.
Метод СalculateDistance осуществляет вычисление расстояния для переданного в качестве элемента экземпляра TransformationHypothesis. Другими аргументами этого метода являются: карта направлений изображения, преобразование евклидовых расстояний карты границ изображения, диаграмма Вороного карты границ изображения, исходный контур и поле направлений исходного контура.
Метод Process осуществляет, собственно, детектирование по алгоритму 2.1.2. На вход методу Process передаются следующие аргументы: диаграмма Вороного карты границ изображения, преобразование евклидовых расстояний карты границ изображения, карта направлений изображения, исходный шаблон и сетка в пространстве гипотез. Для минимизации расстояния используется метод Minimize класса HypothesisGrid, который минимизирует результат метода СalculateDistance. Итеративный подход алгоритма 2.1.2 не реализован; вместо этого всегда производится ровно две итерации детектирования и построения решетки. Как показывает практика, двух итераций достаточно для получения результатов на рассматриваемых размерах изображений. Для сгущения сетки вокруг результата первой итерации используется метод TightenAround. Результатом работы метода является некоторый экземпляр TransformationHypothesis из заданной сетки.
Для осуществления фильтрации градиентного препарата необходимо вычислить пороговую функцию (2.50). Для вычисления этой функции необходимо знать преобразование евклидовых расстояний контура, карту расстояний преобразованного контура, а также диаграмму Вороного.
Если при вычислении преобразования евклидовых расстояний использовать только точки шаблона, при фильтрации на карте границ образуются разрывы, создающие «пунктирный» эффект. Это связано с тем, что точки, принадлежащие границе целевого объекта и лежащие недалеко от ребер контура, при таком способе вычисления могут дать большое расстояние, если они располагаются далеко от вершин контура (например точки между соседними вершинами). И, хотя контур задается своими вершинами, можно утверждать, что именно ребра контура формируют визуальную картину границ целевого объекта. Следовательно, при оценке вероятности принадлежности некоторой точки к реальной границе целевого объекта следует учитывать прежде всего ее близость к ребрам контура, а не к его вершинам.
Оценка близости требует не только оценки геометрического расстояния, но и сравнения направления, ассоциированного с данной точкой (направление вектора градиента интенсивности изображения), и направления, связанного с ближайшей точкой контура (с учетом преобразования). И, хотя оценка геометрического расстояния от данной точки до прямой, порожденной соседними вершинами контура, не представляет сложности, возникает вопрос о том, какое направление следует ассоциировать с точками, принадлежащими ребрам контура.
Для вычисления преобразования евклидовых расстояний с учетом приведенных соображений был разработан и релизован прием «дорисовки» контура с интерполяцией значений направлений между отдельными вершинами. Первый этап алгоритма 2.1.1 вычисления преобразования евклидовых расстояний и диаграммы Вороного состоит в занесении нулевых значений в те точки результирующего изображения, которые содержатся в целевой области (контуре). Метод Calculate класса TemplateEDT вычисляет преобразование евклидовых расстояний, а также диаграмму Вороного и карту направлений в ближайших точках.
Первый этап состоит в рисовании контура на изображении типа Image bool методом BresenhamLine класса FIPL.DrawingDrawing. В результате получается булевское изображение, к которому затем применяется метод
GetDistanceTransform, как и в случае вычисления преобразования евклидовых расстояний для целевого изображения. Карта направлений вычисляется следующим образом: контур рисуется на изображении типа Image double , причем каждой точке рисуемого ребра приписывается направление, отстоящее на 90 против часовой стрелки от направления от текущей вершины к следующей.
Результат этого метода используется для фильтрации градиентного препарата, который представлен в виде экземпляра IMatrix R2Vector . Фильтрация осуществляется классом EdgeSuppressFilter. Данный класс реализует интерфейс библиотеки FIPL ITransform R2Vector,R2Vector . В качестве аргументов этого конструктора ему передаются параметры 0,,,, участвующие в (2.50). Результат вычисления преобразования в виде IMatrix TemplateEDT.EDTElement заносится в свойство EDT этого класса. Фильтрация осуществляется по формуле (2.51) при помощи модуля параллельных вычислений FIPL.