Распознавание эмфиземы лёгких по данным компьютерной томографии

Хроническая обструктивная болезнь лёгких (ХОБЛ) – это хроническое заболевание, характеризующееся персистирующим ограничением скорости воздушного потока. Заболеваемость и смертность от ХОБЛ продолжают расти во всем мире [1]. Компьютерная томография является ведущим методом для оценки морфологических изменений легочной ткани и бронхиального дерева. Наиболее частым проявлением ХОБЛ является эмфизема, которая характеризуется наличием перманентно расширенных воздушных пространств дистальнее терминальных бронхиол с разрушением альвеолярной стенки, что препятствует нормальному обмену кислорода и углекислого газа в крови. На компьютерных томограммах эмфизема проявляется фокальными участками или областями низкой плотности, как правило, без видимых стенок [2] (рис. 1).

а)

б)

Рис. 1. Срез компьютерной томографии лёгких: без патологии (а), с эмфиземой (б)

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

Среди программ по анализу томограмм наибольшим набором функций обладает система распознавания и оценки эмфиземы лёгких от компании General Electric «Thoracic VCAR» (Visualization Computer Reading) [3], которая позволяет провести углубленный анализ лёгочной ткани на предмет различных патологий. Thoracic VCAR является пакетом программного обеспечения для анализа изображений для AW WorkStation, которая, в свою очередь, поставляется только вместе с дорогостоящим оборудованием компании General Electric. Преимуществами Thoracic VCAR являются подсчёт количественных показателей, визуализация зон эмфиземы, построение гистограммы. Из недостатков следует подчеркнуть низкую распространённость и частоту использования, высокую стоимость, невозможность варьирования настроек программы, а также невозможность анализа срезов толщиной более 1 мм.

Нашей целью является исследование метода распознавания эмфиземы лёгких, основанного на вычислении поля локальных средних интенсивности с применением операций математической морфологии, а также создание кроссплатформенного программного обеспечения для анализа лёгких предложенным в работе методом. Данный метод достаточно прост в реализации, а также имеет небольшую вычислительную сложность по сравнению с аналогами при достаточной информативности. Программное обеспечение разработано в рамках проекта Автоплан – системы полуавтоматического планирования вмешательства [4]. Проект «Автоплан» базируется на MITK – систе- ме обработки и визуализации диагностических изображений в формате DICOM [5].

Выделение области интереса – лёгких человека

Метод Оцу (Otsu’s method [6]), благодаря своей универсальности, широко применяется при сегментации изображений в распознавании образов. Алгоритм основан на разделении пикселов полутонового изображения на «полезные» и «фоновые» путём нахождения порога, уменьшающего внутриклассовую дисперсию, которая определяется как взвешенная сумма дисперсий двух классов:

о, (t) = w (t) 02 (t) + w2 (t) 02 (t), где веса Wi - это вероятности двух классов, разделенные порогом t, о, - дисперсия этих классов.

Полученное бинарное изображение имеет следующий вид:

. ,х 1 1, X ( i , j , k ) e D , У ( i , j , k ) = 1

1 0, x ( i , j , k ) £ D ,

где D - область лёгких и воздуха, x ( i , j , k ) - пиксель исходного полутонового изображения, y ( i , j , k ) - пиксель полученного бинарного изображения (1). Каждый индекс i , j , k соответствует одной из осей - x , y , z соответственно. Таким образом, мы пробегаем по всем пикселам трёхмерного изображения.

Метод Оцу разделяет исходное изображение на два основных класса по интенсивности - воздух вокруг тела вместе с лёгкими и само тело. Воздух вокруг тела и лёгкие образуют область D , которая включает в себя область интереса O - лёгкие. Разделим общий класс по интенсивности на две области, используя тот факт, что лёгкие находятся внутри тела и не соединены с окружающим воздухом, который является в данном случае шумом. Для этого ставится точка в области лёгких, тем самым определяется область интереса. После этого обрабатывается область, связанная с этой точкой [7]. Таким образом, полученное изображение имеет вид:

| 1, x ( i , j , k ) e D и x ( i , j , k ) e O , y ( i , j , k ) = 1

[ 0, x ( i , j , k ) g D о x ( i , j , k ) g O ,

где O - связанная с заданной начальной точкой область.

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

Выделение очагов эмфиземы

Данная патология характеризуется однородными потемнениями в лёгких. Эмфизема находится в порогах интенсивности ~ [-1200, -950] по шкале Хаунс-филда [8], в то время как лёгкие в пределах ~ [-1200, -650]. Необходимо отметить, что текстура лёгких достаточно неоднородна и простой пороговой обработки здесь недостаточно, так как необходимо сегментировать очаги заболевания, а не все участки с заданной интенсивностью, которые в данном случае являются шумом.

Вариант двумерной обработки

В этом случае мы проходим по каждому срезу трёхмерного изображения, обрабатывая его независимо от других. Исходным является изображение одного слоя с выделенной областью интереса D : { x ( i , j )} i,■ e d . Далее вычисляем поле локальных средних интенсивностей [9]:

x A ( i , j ) = Л 2 x ( i + m , j + n ), ( i , j ) e D ,     (3)

| A (m, n )e A где A - окно локального усреднения, например: A = {(m, n): |m| < M, |n|

Если среднее значение интенсивности по окну меньше определённого порога, то окно считается очагом заболевания.

Пороговая обработка локальных средних интенсивности:

, x_A (i, j )< h^;

°, ^x_A (i’^j)>^h’

(i, j)eD,

где h - настраиваемое пороговое значение интенсивности.

Далее производим дилатацию бинарного изображения маской B:

|1, 3(m, n) e B : yA (i + n, j + m) = 1, y (i, j) = 1

[0, V(m, n) e B : yA (i + n, j + m) = 0,     (5)

(i, j)e D, где B = {(m, n): |m|

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

Вариант трёхмерной обработки

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

Дано трёхмерное изображение:

{x (i, j, k)}i,j=i:n .                                                     (6)

k=1, L

Область интереса: {x(i,j, k)}i,*ed.

Поле локальных средних:

xA (i, j, k ) = A 2 x (i + m, j + n, k +1), | A| (m, n,k )e A

(i, j, k)e D,                                               (7)

A = {(m,n, 1): |m| <M/2, |n| <M/2, |1| <L/2}.

Пороговая обработка локальных средних интенсивности:

У a (i, j, ^k ) =

1, Xa (i, j, k )^ h; , A

_                  (i, j,k)e D.

0, Xa (i, j, k)< h,

В качестве маски B мы выбрали квадратное окно. Дилатация бинарного изображения маской B:

y (i, j, k) =

1, 3(m,n, k) e B: y^(i + n, j + m,1 +1) = 1, 0, V(m,n,k) e B: y^(i + n, j + m,k +1) = 0,

(i, j, k )e D,

B = {(m,n, 1): |m| <M/2, |n| <M/2, |1| <L/2}.

Данный алгоритм является эффективным на средних и больших размерностях изображений (более 200 слоёв). В задачах малой размерности результат не имеет никаких преимуществ, так как размерность по координатам x и y в несколько раз превышает размерность координаты z, следовательно, статистическая зависимость соседних срезов достаточно слабая, и мы можем ей пренебречь в целях экономии времени выполнения программы.

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

Настройка параметров алгоритма

Введём понятие эффективности алгоритма:

E = min M, L T ,

где T – время выполнения алгоритма, M, L – размер окна сглаживания и дилатации. Алгоритм является эффективным, если соотношение его параметров дает высокую скорость работы.

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

Имеется «эталонное» изображение:

где G ′ – область, классифицируемая как патологическая.

Находим количество неверно классифицированных пикселей:

Nek = Z(i,j)ed|^X3k⁽i,^j)-^X₈k⁽i,^j)|, k = i^n ,        ⁽¹³⁾

где D – область интереса, k – номер изображения.

Далее находим вероятность ошибочной классификации пикселя на изображении (9):

N

Рош = -D^ ,                                    (14)

где |D| – количество пикселей в области интереса, N;_k - количество ошибочно классифицированных пикселей в k-м изображении.

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

Рош (м, h) = Z Рош k ( m , _h),                          (15)

где K – количество пациентов в обучающей выборке.

Экспериментальные исследования

Параметры алгоритма распознавания эмфиземы были получены на основе обучающей выборки, затем алгоритм был протестирован на контрольной выборке. Обучающая выборка состоит из 15 пациентов с различной степенью тяжести заболевания. Контрольная выборка состоит из 10 пациентов с лёгкой степенью тяжести. В качестве оценки качества работы алгоритма использовалась экспертная оценка врачей-рентгенологов по следующей схеме. На контрольной выборке срез за срезом была отмечена область заболевания, затем рассчитывалась ошибка перекрытия (Volume overlap error, VOE) по следующей формуле:

V = 1 - e k

N^(хэ„ ⁽i, ^j) A ^x₈, ⁽i, ^j)) ₇ — ------^k-------------k--------, k = 1, n ,

N(x^ (i, j) v x8 (i, j))

N(Хэ_к (i, j) a x8 (i, j))

где ^{k        k} – коэффициент сходства

N^{( Х}э„⁽i, ^j) ^{v X}8_k⁽i, ^j ))

x_э (i, j)

1, x (i, j) e G, .0, x(i, j) г G,

где G – область очагов эмфиземы, отмеченная вра-

чом-экспертом.

Алгоритм зависит от двух параметров: размер сглаживающей маски и порог отсечения области патологии. Оптимальные размеры матрицы Mопт и пороговое значение hопт определяются на основе экспертной выборки из 6 образцов. Для этого меняем

значение порога и размеры матрицы и получаем изображение с некоторой погрешностью:

1, x (i, j) e G, x8 / i, j) = I               • k 0, x (i, j) г G ,

Жаккара. Ошибка перекрытия является стандартным средством оценки качества медицинских объемов.

Двумерный срез трёхмерного изображения представлен на рис. 2.

Рис. 2. Срез изображения контрольной выборки

Находим оптимальные параметры алгоритма перебором параметров. Результаты исследований представлены в табл. 1.

Таким образом, получаем оптимальное соотношение параметров алгоритма, при котором вероятность ошибочной классификации пикселя Pош минимальна и равна 0,011 при параметрах hопт = –1020, M_опт =5.

Табл. 1. Нахождение оптимальных параметров алгоритма

M×M	Порог, h
	-1028	-1024	-1020	-1016	-1012
1×1	0,104	0,157	0,35	0,200	0,240
3×3	0,019	0,014	0,019	0,16	0,125
5×5	0,079	0,046	0,015	0,023	0,063
7×7	0,112	0,081	0,067	0,027	0,024
9×9	0,134	0,105	0,070	0,065	0,017
11×11	0,152	0,121	0,092	0,059	0,1

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

Необходимо отметить, что из-за различных анатомических особенностей лёгких, а именно: текстуры, интенсивности и объёма, нельзя однозначно выделить универсальный оптимальный порог распознавания. Для устранения данной неточности в плагине введён выбор порога оператором (по умолчанию ставится оптимальное значение).

После выполнения плагина программа визуализирует 3 проекции исходного изображения: аксиальное, сагиттальное и корональное. Также показана трёхмерная визуализация лёгких и эмфиземы. Результат работы программы – на рис. 3.

Рис. 3. Вид системы после работы плагина

Тёмным цветом визуализируются очаги заболевания, наложенные на исходное полутоновое изображение.

В приведённом случае при верхнем пороге –920 плагин наиболее точно распознал эмфизему и показал следующие характеристики: среднее значение – 809,525; процент эмфиземы 9,67774 при объёме лёгких 5,80445. VOE в данном случае 0,019, т.е. чуть более 2 %. Средняя VOE по контрольной выборке составила 0,076.

Среднее значение интенсивности и гистограмма (числовые характеристики) используются для оценки заболевания как более простой, но менее надежный способ. Эти факторы коррелируют с тяжестью заболевания, и по ним также можно сделать прогноз о состоянии больного и эффективности лечения, предположить развитие различных осложнений, например, хронического легочного сердца, которое зачастую приводит к летальному исходу у больных ХОБЛ [2].

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

Заключение

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

Работа выполнена в рамках проекта «Автоплан» (система автоматического планирования и контроля за ходом операции) по государственному контракту на выполнение научно-исследовательской и опытноконструкторской работы «Разработка технологии и организация производства систем автоматизированного планирования, управления и контроля результатов хирургического лечения» («4.3-Автоплан-2014»).