УДК 004.9
ОЦЕНИВАНИЕ ЧИСЛА КЛЕТОК НА ИЗОБРАЖЕНИЯХ ЦИТОЛОГИЧЕСКИХ РАСТИТЕЛЬНЫХ ПРЕПАРАТОВ
М.С. Тарков, М.И. Осипов*
Институт физики полупроводников СО РАН им. А.В. Ржанова, г. Новосибирск E-mail: tarkov@isp.nsc.ru 'Новосибирский государственный университет E-mail: inq1201@mail.ru
Предложен алгоритм для подсчета вакуолей (растительных клеток) на изображении. Алгоритм позволяет выделить контуры клеток в массивах слипшихся клеток путем разбиения диапазона яркости изображения на слои с последующим применением в каждом слое полутоновой эрозии и дилатации, преобразования «верх шляпы» и пороговой сегментации. При подсчете клеток используется аппроксимация их контуров эллипсами. Алгоритм может быть использован в различных приложениях, например при определении скорости распада клеточных структур по серии изображений.
Ключевые слова:
Оценивание числа клеток, цитологические растительные препараты, анализ изображений.
Key words:
Estimating the number of cells, cytological herbal preparations, image analysis.
Введение
В настоящее время анализ изображений цитологических препаратов осуществляется преимущественно вручную, что сильно замедляет процесс, а плохое качество изображения или наличие на изображении большого числа клеток увеличивает вероятность ошибки при выяснении тех или иных характеристик препарата. В частности, при большом объеме (102—103) выборки подсчет клеток вручную может занять много времени, а также возможен некорректный подсчет.
Использование автоматических методов анализа изображений цитологических препаратов может значительно улучшить качество анализа и ускорить его, однако универсальный алгоритм разработать трудно в силу разнородности изображений. Клеточные структуры могут выглядеть совершенно по-разному, иметь прозрачную, полупрозрачную или непрозрачную структуру, обладать ярко выраженными границами или не обладать таковыми. Поэтому целесообразно разрабатывать алгоритмы под конкретные классы изображений, что повысит их эффективность, но, очевидно, применение их будет возможно только в пределах соответствующих классов.
Разработка алгоритмов для подсчета числа клеток ведется, и для некоторых классов изображений они уже есть [1, 2], но такие алгоритмы обычно основаны на сегментации, то есть целью является сопоставление клетки с конкретным сегментом изображения. Такой подход эффективен в том случае, когда клетки непрозрачны, то есть на изображении цвет такой клетки значительно отличается от цвета фона, что позволяет выяснить, является ли сегмент клеткой. К тому же при наличии шума на изображении, а при съемке цитологических препаратов это не редкость, возможно появление множества мелких сегментов, сходных по цветовым характеристикам с клетками, вследствие чего получится некачественный подсчет. Существуют классы изображений, для которых использование алгоритмов, основанных на сегмента-
ции, неприемлемо. Один из таких классов образуют изображения вакуолей (растительных клеток) [3].
Целью данной работы является разработка алгоритма для подсчета числа вакуолей на изображении. При разработке необходимо учитывать специфику изображений данного класса, поэтому необходимо решить ряд проблем, а именно (рис. 1):
1. Клетки могут слипаться или перекрываться, вследствие чего несколько клеток могут восприниматься как одна, поэтому необходимо решить проблему разделения слипшихся клеток.
2. Клетки полупрозрачны, поэтому цветовые характеристики клеток и фона практически совпадают, что препятствует использованию сегментации изображения классическими способами, так как затруднительно по цветам пикселей некоторого сегмента определить, является ли он клеткой или нет.
3. Необходимо решить проблему отделения клеток от фона, так как в некоторых случаях клетки имеют нечеткую границу.
Рис. 1. Фрагмент снимка растительных клеток
Используемый подход
Для решения задачи используется подход, предложенный в [1]. После некоторых модификаций, связанных со спецификой рассматриваемого нами класса изображений, этот подход применим к изображениям растительных клеток (вакуолей). Особенность подхода заключается в том, что после выделения контуров клеток исходное изображение больше не используется. Для дальнейшего анализа используется бинарное изображение контуров. Анализируются геометрические свойства полученных контуров, что позволяет выделять слипшиеся клетки даже в случае, когда граница между ними плохо различима.
В алгоритме можно выделить два этапа:
1. Выделение и предобработка контуров клеток. Выделение контуров осуществляется с помощью комбинированного использования морфологических преобразований эрозии и дилата-ции [4], а также порогового преобразования карты евклидовых расстояний. Контуры аппроксимируются многоугольниками: из каждого контура выделяются точки, являющиеся вершинами аппроксимирующего многоугольника. На контурах отыскиваются точки вогнутости. Этими точками контур сегментируется на части.
2. Выделение клеток с помощью эллипсов. Строятся эллипсы, аппроксимирующие контуры клеток. Выявляются эллипсы, которые не могут аппроксимировать клетку (слишком узкие или слишком длинные). Объединяются сегменты. Обрабатываются эллипсы, не отобранные на этапе селекции. Они объединяются с другими наиболее подходящими эллипсами. Производится коррекция выделения и подсчет эллипсов.
Алгоритм выделения контуров
Перед тем как приступить непосредственно к нахождению контуров, входное цветное изображе-
ние преобразуется в монохромное. Далее осуществляется следующая последовательность действий:
1. Применение морфологических операций эрозии и дилатации с ядром 7x7.
2. Выравнивание гистограммы.
3. Сглаживание изображения путем свертки его с гауссовым ядром 3x3.
4. Вычитание фона.
5. Выделение десяти слоев (интервалов яркости) размером 20 ([0,19], [20,39], ..., [180, 199]).
6. Каждый слой подвергается последовательному применению морфологических операций дила-тации и эрозии для удаления шумовых выбросов. Также производится пороговое преобразование евклидовой карты расстояний для выделения контуров. Выделение контуров для каждого слоя производится независимо.
Вычитание фона является необходимой процедурой. Её задача - не допустить, чтобы в рассматриваемые слои попали фрагменты фона (рис. 2). Рис. 2 показывает, что после вычитания фона значительная его часть не попадает в рассматриваемые далее слои яркости. Чтобы выделить фон на изображениях данного класса необходимо [3]:
1. Преобразовать исходное изображение в монохромное.
2. Осуществить медианную фильтрацию.
3. Выровнять гистограмму.
4. Взять дополнение.
5. Построить размыкание изображения. В качестве структурообразующего элемента используется круг с радиусом 25 пикселей.
6. Вычесть размыкание из изображения, полученного в ходе шагов 1-4 (преобразование «верх шляпы» [4]).
7. Осуществить пороговое преобразование (порог подобран эмпирически и равен 25) (рис. 3). Пороговое преобразование карты евклидовых
расстояний [5] в данном случае применяется к бинарным изображениям (каждый пиксель является либо черным, либо белым). Карта евклидовых рас-
стояний изображения составляется таким образом: для каждого пикселя вычисляется расстояние до ближайшего черного пикселя (если пиксель черный, оно равно нулю). Пороговое преобразование состоит в изменении цвета пикселя с белого на черный, если его расстояние до ближайшего черного пикселя не попадает в заданный промежуток значений.
Аппроксимация контуров многоугольниками
Аппроксимация контуров многоугольниками основана на двух соображениях. Первое - сгладить контуры, чтобы при поиске точек вогнутости не возникало много ложных точек. Второе - уменьшить число точек контуров. После аппроксимации, как правило, количество точек контура уменьшается. Как следствие, уменьшается количество информации о геометрии контура и сокращается объем последующих вычислений.
Для аппроксимации используется алгоритм Дугласа-Пекера [6]. Алгоритм рекурсивно делит линию. Входом алгоритма служат координаты всех точек контура кроме первой и последней, которые сохраняются неизменными. Далее алгоритм находит точку, наиболее удалённую от прямой, проведённой через первую и последнюю точки контура. Если точка находится на расстоянии меньше, чем s, то все точки, которые ещё не были отмечены к сохранению, могут быть выброшены из набора, и получившаяся прямая сглаживает кривую с точностью не ниже s. Если же расстояние больше s, то алгоритм рекурсивно вызывает себя на наборах точек от начальной до данной и от данной до конечной. Это означает, что данная точка будет отмечена к сохранению. По окончании всех рекурсивных вызовов результирующая ломаная строится только из тех точек, что были отмечены к сохранению.
Поиск точек вогнутости и сегментация контуров
Точки вогнутости - это точки, разделяющие части контура, которые принадлежат различным клеткам. Обозначим pc(xc,yc) текущую точку, а PPre(xFre,yFre) и p„Jxmxt,ymxt) предыдущую и следующую точки соответственно. Далее определим величину concavity (pc) в точке рс [1]:
a(Ppre , Pc) = tan- ((ypre - yc) /(x„re - % ));
a(Pnex,, Pc ) = tan-1((y„ex, - Ус )/(X,cx, - xc ));
I I
\a( Ppre , Pc ) - a( Pnc, , Pc )|. concavity (Pc) = j если \a(Ppre, pc ) - a(Pmxt, pc ) | < n,
n-\a(PPre , Pc ) - a(Pnex, , Pc ^ ина4е-
Условия, при выполнении которых точка рс считается точкой вогнутости, выглядят так [1]:
1. a1<concavity (pc)<a2;
2. Отрезок pFpnext не должен проходить внутри контура.
Значения a1 и a2 задаются пользователем. Сегментация контуров предусматривает разделение контуров на части найденными точками вогнутости (рис. 4). Отдельные части контуров предположительно являются частями контуров отдельных клеток. По этим фрагментам осуществляется построение эллипсов [1]. Однако сегменты контура в ходе дальнейшей работы алгоритма могут объединяться, поэтому их первоначальное количество в общем случае не равно числу клеток, обнаруженных в области, контур которой рассматривается.
Рис. 4. Пример контура из двух сегментов, разделенных точками вогнутости (точки вогнутости обведены белыми окружностями)
Выделение клеток
Основной целью аппроксимации эллипсами является разделение выделенных областей на отдельные клетки. На данном этапе мы имеем набор контуров, каждый из которых разделен на сегменты точками вогнутости. Изначально для каждого сегмента строится свой эллипс. Далее выполняется анализ геометрических свойств эллипсов и их расположения относительно друг друга. Целью анализа является уточнение результатов, полученных при первоначальном построении.
Построение эллипса производится для каждого фрагмента каждого контура. Сам фрагмент представляет собой набор точек, которые по предположению принадлежат контуру клетки. Однако эти точки образуют только часть реального контура, а построенный эллипс должен дополнить картину (рис. 5).
Рис. 5. Пример эллипсов, построенных по точкам сегментов Рассмотрим уравнение кривой второго порядка: (А, X) = 0,
где Л=[а,Ь,е,й,е,/], Х=[х2,ху,у2,х,у,1]. Если для точки (х;,у;) сформировать вектор Х=[х,2,ху;,уДх;,у;,1], то можно определить алгебраическое расстояние (Л,Х) от точки X до кривой (Л, Х)=0.
Коэффициенты для искомого эллипса находятся методом наименьших квадратов, то есть решается задача условной оптимизации [7]. Минимизируется сумма квадратов алгебраических расстояний при условии, что получившаяся в результате кривая второго порядка - эллипс:
К (Л) = Х (л, X, )2
1 + с) ■
Ъ/2 а/2
Ъ/2 с а/2 е/2
е/2 /
< 0,
Ъ2 - 4ас < 0,
где N - количество точек фрагмента контура.
Минимизация проводится по вектору коэффициентов квадратичной формы Л=[а,Ь,с,С,е/|, а ограничения-неравенства гарантируют, что в результате минимизации получится эллипс.
Селекция
Процедура селекции заключается в проверке двух условий для каждого эллипса. Если эти условия выполнены, эллипс участвует в следующем шаге - объединении сегментов. Остальные эллипсы обрабатываются отдельно.
Выбираются эллипсы, которые достаточно точно аппроксимируют контур, и те, которые имеют приемлемую форму для аппроксимации контура клетки, то есть не являются слишком узкими или слишком вытянутыми. Чтобы эллипс участвовал в следующем шаге, необходимо выполнение следующих двух условий [1]:
1.
Обозначим сегмент контура Ц={р;1,...,рй}, р =(х1ру ) а соответствующий этому сегменту эллипс F(ЛІ,X). Определим алгебраическое расстояние йй (Л;,Ц) от точек сегмента Ц до соответствующего ему эллипса:
1.
<Ив( 4, ц) = - X |( Л, X )|'
5}=1
Первое условие: величина йй(Л;,Ц) должна быть меньше заранее выбранного порога ШяТН.
2. Обозначим МахЛ и МЫЛ. длины большей и малой полуосей эллипса соответственно. Второе условие: отношение еЯаИо=МшЛ/МахЛ должно быть не меньше заданного порога еТН, что означает близость эллипса к окружности. Порог еТН необходимо выбирать, исходя из формы рассматриваемых клеток.
Объединение сегментов
На данном шаге объединяются сегменты, эллипсы которых построены по фрагментам контуров одной клетки (рис. 6). Берется два фрагмента, строится новый эллипс по точкам из этих фрагментов, проверяется выполнение условий селекции, и в зависимости от ситуации принимается решение: оставить два изначально построенных эллипса или заменить их новым эллипсом, объединив рассматриваемые сегменты.
Возможны следующие ситуации [1]:
1. Явно соприкасающиеся клетки должны быть разделены, т. е. сегменты объединяться не должны. Для идентификации данной ситуации необходимо пользоваться следующими признаками:
а) Расстояние между центром нового эллипса С„т и центрами старых эллипсов С и С2 должно быть больше, чем величина йМтТН, то есть
IIС -С_||> йЫтТН, I = 1,2.
При соблюдении этого условия эллипсы не объединяются. В качестве йМтТН используется величина, близкая к длине меньшей полуоси наименьшей клетки.
б) Расстояние между центрами старых эллипсов слишком велико для того, чтобы они были частями одной клетки, а именно
(=1
00 О о
||с1 - С2|| > 3• йЫтТН.
При выполнении этого условия эллипсы также не объединяются.
2. Два несмежных сегмента должны быть объединены, то есть нужно заменить два старых эллипса новым, построенным по точкам обоих сегментов. Этот случай идентифицируется следующими признаками:
а) Длины меньших полуосей старых эллипсов меньше величины йЫтТН, то есть эллипсы слишком малы для аппроксимации контура какой-либо клетки, поэтому их нужно объединить.
б) Эллипсы, построенные по сегментам контура схожи, то есть длины их больших и малых полуосей соответственно близки:
\МахА1 -МахА2\ < йЫтТН,
\МтА1 -МтА2\ < йЫтТН /2.
Такие эллипсы объединяются в один, так как, вероятно, они аппроксимируют одну и ту же клетку.
3. В остальных случаях необходимо пользоваться следующим правилом: если йЪ (Спт,Ьпа)<йЪ(СхЦ^) или ййС^Д^йадЛ), где 1т=1^1ъ то старые эллипсы заменяются новым.
Обработка эллипсов, не отобранных в ходе селекции
На данном шаге рассматриваются эллипсы, для которых условия селекции не выполнены (необработанные сегменты). Поскольку они являются частями контуров, то соответствующие им сегменты объединяются с наиболее подходящими сегментами эллипсов, прошедших процедуру селекции (т. е. с обработанными сегментами).
Пусть для текущего контура мы имеем т необработанных сегментов Ьь...,Ьт и п обработанных ЕЪ..,ЕЛ. Берется сегмент Ц, /=1..т, строится последовательность Р, точек Ц и Е, у'=1..п. По последовательности точек Р, строится эллипс Т, и считается величина йй (Т,р). Находится сегмент Ек, ке{1,.., п}, для которого эта величина наименьшая. Объединение Ц с Ек происходит только в том
случае, если эллипс, построенный по точкам двух сегментов, удовлетворяет условиям селекции, иначе сегменты не объединяются, а сегмент Lt исключается из дальнейшего рассмотрения.
Коррекция выделения клеток и подсчет
Растительные клетки (вакуоли) имеют форму круга. Поэтому для лучшего визуального восприятия результирующее выделение нужно сделать при помощи окружностей, а не эллипсов. По полученному набору эллипсов достаточно просто строится набор окружностей следующим образом: по каждому эллипсу строится окружность с центром, совпадающим с центром эллипса, а радиус окружности равен полусумме полуосей эллипса.
Контуры клеток выделяются из различных слоев яркости. При этом возможны попадания разных частей одной и той же клетки в разные слои. В этом случае произойдет наложение окружностей друг на друга. Данный шаг нужен для отслеживания чрезмерного наложения. Рассматриваются все окружности. Если площадь пересечения некоторого круга с другими кругами составляет более 75 % собственной площади, то данный круг удаляется и в дальнейшем не рассматривается. Результирующее число клеток равно числу оставшихся кругов.
Эксперименты
При разработке алгоритма большое внимание уделялось подбору различных параметров, от которых существенно зависит качество работы алгоритма. Подбор осуществлялся, исходя из геометрии клеточных структур и специфики рассматриваемых изображений. Алгоритм программно реализован в среде Visual Studio C++ с использованием процедур обработки изображений библиотеки OpenCV [8]. В табл. 1 приведены значения параметров алгоритма, при которых получены удовлетворительные результаты. Здесь е - параметр аппроксимации в алгоритме Дугласа-Пекера, а1 и а2 - параметры алгоритма поиска точек вогнутости (в градусах), disTh и dMinTh - пороги расстояний, eTh -порог отношения полуосей эллипса.
Рис. 7. Тестовые образцы изображений растительных клеток
Таблица 1. Значения параметров алгоритма подсчета клеток
Параметр Е 31 32 (НэТЬ еТЬ с1М1пТН
Значение 2 20 160 0,1 0,1 8
В табл. 2 и 3 приведены результаты работы алгоритма для изображений рис. 7, а-г при различных параметрах, используемых в процессе извлечения контуров. В скобках указано отклонение в процен-
тах от результата, полученного подсчетом вручную. На основе экспериментальных данных выбраны оптимальные параметры алгоритма, дающие минимальную по последовательности кадров суммарную ошибку подсчета клеток. При выборе способа разбиения слоев учитывался тот факт, что основные части клеток находятся в диапазоне яркости от 0 до 200. Слой в диапазоне яркости от 201 до 255 содержит фон и мелкие фрагменты клеток.
Таблица 2. Результаты работы алгоритма при различных размерах ядра операторов эрозии идилатации, применяемых при извлечении контуров
Размер ядра 3x3 5x5 7x7 9x9 Результат подсчета вручную
Изображение (рис. 7)
а 879 (0,1 %) 817 (7,1 %) 822 (6,5 %) 594 (32,5 %) 880
б 480 (18,2 %) 451 (11,0 %) 391 (3,6 %) 354 (12,8 %) 406
в 370 (37,9 %) 351 (28,1 %) 303 (10,5 %) 291 (6,2 %) 274
г 342 (69,3 %) 324 (60,3 %) 278 (37,6 %) 264 (30,6 %) 202
Таблица 3. Результаты работы алгоритма при различном количестве и толщине слоев яркости, на которые разбивается изображение при извлечении контуров
Кол-во слоевхтолщина (пикс.) 2x100 5x40 10x20 12x15 Результат подсчета вручную
Изображение (рис. 7)
а 508 (42,2 %) 683 (22,3%) 822 (6,5%) 612 (30,4%) 880
б 284 (30,0%) 410 (0,9%) 391 (3,6%) 365 (10,1%) 406
в 367 (33,9%) 322 (17,5%) 303 (10,5%) 287 (4,7%) 274
г 249 (23,2%) 292 (44,5%) 278 (37,6%) 265 (31,1%) 202
В качестве тестовых образцов использовались и зображения размером 1600x1200 пикселей (рис. 7). Алгоритм тестировался на компьютере с процессором Intel Core i3-2310M, 2.1 GHz. Время обработки одного изображения равно 7,3 с. На рис. 8 представлен результат работы алгоритма. Количество обнаруженных клеток - 822.
Рис. 8. Пример результата работы алгоритма выделения клеток
Заключение
Разработан алгоритм подсчета растительных клеток (вакуолей) на изображениях клеточных препаратов. Алгоритм основан на идее выделения контуров клеток с последующей аппроксимацией контуров эллипсами (окружностями). Известные подходы, основанные на сегментации изображений, неприменимы к подсчету объектов на изображениях растительных клеток, что подтверждено экспериментами. Реализован метод выделения контуров слипшихся клеток. В алгоритме применен ряд дополнительных операций, таких как вычитание фона и коррекция выделения клеток, что позволило значительно повысить качество работы алгоритма.
Экспериментально подобраны параметры, обеспечивающие приемлемое качество работы алгоритма. Алгоритм реализован в среде Visual Studio C++ с использованием процедур обработки изображений библиотеки OpenCV. Алгоритм примечателен тем, что вся работа с исходным изображением заканчивается после выделения контуров, степень четкости границ клеток (если они слипшиеся) не играет роли, так как дальше анализируется только геометрия контуров.
Авторы выражают благодарность В.Н. Нурминскому (СИФИБР СО РАН) за предоставленные изображения клеточных препаратов.
СПИСОК ЛИТЕРАТУРЫ
1. Bai X., Sun C., Zhou F. Splitting touching cells based on concave points and ellipse fitting // Pattern Recognition. - 2009. - № 42. -P. 2434-2446.
2. Ковригин А.В. Применение принципов построения систем машинного зрения в задаче анализа изображений клеточных структур // Научный журнал КубГАУ. - 2007. - № 29 (5). URL: http://ej.kubagro.ru/2007/05/pdf/03.pdf (дата обращения: 12.03.2012).
3. Тарков М.С., Тихонов Н.В., Половинкин В.Г. Оценивание состава изображений клеточных препаратов для медико-биологических исследований // Известия Томского политехнического университета. - 2012. - Т. 320. - № 5. - С. 130-133.
4. Гонсалес Р., Вудс Р. Цифровая обработка изображений. - М.: Техносфера, 2005. - 1072 с.
5. Danielsson P-E. Euclidean distance mapping // Computer graphics and image processing. - 1980. - № 14. - P. 227-248.
6. Douglas D., Peucker T Algorithms for the reduction of the number of points required to represent a digitized line or its caricature // The Canadian Cartographer. - 1973. - V. 10. - № 2. - P. 112-122.
7. Fitzgibbon A., Pilu M., Fisher R.B. Direct Least Square Fitting of Ellipses // IEEE Transactions on Pattern Analysis and Machine Intelligence. - 1999. - V. 21. - № 5. - P 476-480.
8. Bradski G., Kaehler A. Learning OpenCV. - USA: O’Reilly Media, Inc., 2008. - 576 p.
Поступила 03.03.2013 г.