CC BY
153
УДК 004.932.2
АЛГОРИТМЫ СПЕКТРАЛЬНО-ТЕКСТУРНОЙ СЕГМЕНТАЦИИ СПУТНИКОВЫХ ИЗОБРАЖЕНИЙ ВЫСОКОГО ПРОСТРАНСТВЕННОГО РАЗРЕШЕНИЯ
И. А. Пестунов, С. А. Рылов
SPECTRAL-TEXTURAL SEGMENTATION ALGORITHMS FOR SATELLITE IMAGES WITH HIGH SPATIAL RESOLUTION
I. A. Pestunov, S. A. Rylov
Предлагаются алгоритмы спектрально-текстурной сегментации изображений. Приводятся результаты экспериментов на модельных и реальных изображениях, подтверждающие вычислительную эффективность предложенных алгоритмов и высокое качество получаемых картосхем.
Algorithms for multispectral image segmentation based on joint use of spectral and textural features are proposed. Results on the model and real images experiments confirming their efficiency are presented.
Ключевые слова: спутниковые изображения высокого разрешения, сегментация изображений, спектральные и текстурные признаки.
Keywords: satellite images with high spatial resolution, image segmentation, spectral and textural features.
Введение
Сегментация является одним из важнейших этапов анализа цифровых изображений [1, с. 346 - 362]. Она заключается в разбиении изображения на непересекающиеся области на основе однородности (похожести) их спектральных и (или) пространственных (текстура, форма, контекст и т. п.) характеристик. Методы сегментации изображений нашли широкое применение во многих прикладных областях, включая дистанционное зондирование Земли из космоса [2, с. 31 - 42].
С каждым годом растет число запускаемых спутников, обеспечивающих поставку многоспектральных изображений высокого пространственного разрешения (10 м и менее) [3, с. 11 - 17]. Их особенность состоит в том, что значительная часть информации о представленных на них объектах заключена в текстуре изображения. Применительно к таким изображениям, традиционные методы сегментации, учитывающие лишь спектральные характеристики, оказываются непригодными, поскольку приводят к чрезмерной фрагментации однородных по текстуре областей. Не являются исключением и методы, реализованные в известных пакетах программ обработки спутниковых изображений (ERDAS Imagine, ENVI и др.).
В связи с этим разработка методов и программно-алгоритмических средств так называемой спектрально-текстурной сегментации многоспектральных изображений в настоящее время является весьма актуальной задачей [4, с. 2479 - 2501].
В данной работе предлагается метод спектрально-текстурной сегментации спутниковых изображений, основанный на непараметрических алгоритмах кластеризации и свободный от проблемы введения единой метрики в пространстве спектральных и текстурных признаков.
Основная идея предлагаемого метода заключается в разделении процесса сегментации на два этапа. На первом этапе производится сегментация исходного изображения по спектральным признакам, а на втором
- сегментация картосхемы, полученной на первом этапе, по текстурным признакам. Такой подход позволяет избежать традиционных трудностей, связанных с введением единой метрики, которая необходима при кластеризации данных в пространстве разнородных признаков.
Ниже приводится описание алгоритмов, применяемых при реализации этого метода.
Алгоритм сегментация изображения по спектральным признакам
Для сегментации изображения по спектральным признакам применяется ансамблевый алгоритм кластеризации ECCA (Ensemble of Combined Clustering Algorithms) [5], краткое описание которого приводится ниже.
Пусть X = {xi = (x1,...,xd) е Rd, i = 1,N} - множество векторов спектральных яркостей, характеризующих d - спектральное изображение. Векторы xi лежат в прямоугольном гиперпараллелепипеде
W = [/1,r1] х ... х [/d,rd], где lj = minx еХ xj , rj = maxx еХ x j . Под сеточной структурой будем понимать разбиение пространства спектральных признаков гиперплоскостями:
xj = (rj — lj) ■ i / m + lj, i = 0,...,m
(m - число разграниченных участков по каждой размерности). Минимальным элементом этой структуры является клетка (замкнутый прямоугольный гиперпараллелепипед, ограниченный гиперплоскостями). Введём общую нумерацию клеток (последовательно от одного слоя клеток к другому).
Клетки Bi и Bj (i ^ j) являются смежными, если их пересечение не пусто. Множество смежных с B
клеток обозначим через AB . Плотностью DB клетки B назовём отношение DB = NB / VB , где NB - количество элементов множества X, попавших в клетку B; VB - объем клетки B. Клетку B будем считать непустой, если DB > т, где т - величина заданного порога. Все точки множества X, попавшие в клетки с плотностью меньше т, будем относить к «шуму». Обозначим множество всех непустых клеток через N . Непустая клетка Bt непосредственно связна с непустой клеткой Bj (Bt ^ B-), если Bj - максимальная
по номеру клетка, удовлетворяющая условиям:
B- = argmaxB е A DB и DB > DB .
j Bkе ABi Bk Bj Bi
Непустые клетки Bt и Bj непосредственно связны (Bt ^ Bj), если B. ^ Bj или Bj ^ Bt. Непустые клетки Bt и Bj связны (Bt ~ Bj), если существуют k1,...,kl такие, что k1 = i, kl = j и для всех
p = 1,...,1 — 1 выполнено B, ^ B, .
кр кр+1
Введение отношения связности порождает естественное разбиение множества непустых клеток на компоненты связности {Gi, i = 1,...,S}. Под компонентой связности будем понимать максимальное множество попарно связных клеток. Представителем компоненты связности G назовём клетку Y(G), удовлетворяющую условию:
Y(G) = argmaxBeG Db
(если несколько клеток удовлетворяют данному условию, то Y(G) выбирается из них случайным образом). Компоненты связности G' и G" смежные, если существуют смежные клетки B' и B" такие, что B' е G' и B" е G" . Смежные компоненты связности Gt и G- связны (Gi ~ G-), если существует набор клеток (путь)
P = Y = вч B„, •...•B,, = Yj}
такой, что:
1) для всех t = 1,...,l — 1 B, е G. U G и B, , B, - смежные клетки;
kt i j kt kt+1
2) minB, DBt / min(DB^ , DBy ) > T, T > 0 - порог объединения.
Определение. Кластером C назовем максимальное множество попарно связных компонент связности, то есть 1) для любых компонент связности Gi, Gj е C выполнено Gi ~ Gj; 2) для любых Gi е C,
G^ £ C верно Gt / G^.
С учетом изложенного, задача кластеризации заключается в разбиении множества N на совокупность кластеров {с., i = 1,...,M} таких, что N = U^C. и C. П Cj = 0 при i ^ j; число кластеров M заранее неизвестно.
Далее опишем эффективный метод решения этой задачи, основанный на ансамблевом подходе. Предлагаемый метод опирается на сеточный алгоритм CCA [6], где m - число разбиений, T - порог
объединения компонент связности, т - порог шума). В работе этого алгоритма можно выделить три основных этапа.
1. Формирование клеточной структуры. На этом этапе для каждой точки xi е X определяется содержащая её клетка, вычисляются плотности DB всех клеток и выявляются непустые клетки.
2. Выделение компонент связности G1,...,GS и поиск их представителей Y(G1),...,Y(GS).
3. Формирование кластеров C1,...,CM в соответствии с указанным выше определением на основе выделенных компонент связности.
Алгоритм CCA(m, T, т) является вычислительно эффективным в пространстве признаков небольшой размерности (< 6) [6, с. 149 - 152], его сложность составляет O(dN + dmd), где N - число классифицируемых объектов, d - размерность пространства признаков.
Однако CCA относится к классу алгоритмов с фиксированной сеткой и поэтому результаты его работы существенно зависят от параметра m, который определяет масштаб элементов сеточной структуры. На практике эта неустойчивость результатов значительно затрудняет настройку параметров алгоритма.
Известно [7, с. 583 - 617], что устойчивость решений в задачах кластеризации может быть повышена благодаря формированию ансамбля алгоритмов и построению на его основе коллективного решения. При этом используются результаты, полученные разными алгоритмами либо одним алгоритмом с различными значениями параметров. Кроме того, для формирования ансамбля могут быть применены разные подсистемы переменных. Ансамблевый подход является одним из наиболее перспективных направлений в кластерном анализе [8, с. 651 - 666].
В описываемом алгоритме ECCA для формирования ансамбля используются результаты выполнения алгоритма CCA(m, T, т) с различными значениями параметра m, а для формирования итогового коллективного решения - применить метод, основанный на нахождении согласованной матрицы подобия (или различия) объектов [9]. Этот метод может быть описан следующим образом.
Пусть с помощью некоторого алгоритма кластеризации (m = m(©)), зависящего от случайного вектора параметров © е 0 (здесь 0 - некоторое допустимое множество параметров), получен набор частных решений Q = {Q (1),...,Q(l \...,Q(l)] , где Ql) - l -й вариант кластеризации, содержащий M(l) кластеров.
Обозначим через H(©l) бинарную матрицу H(©l) = {Ht,(©l)} размерности N х N , которая для l -й группировки вводится как
Го, если объекты отнесены в один
Hi,3 (©l) = j ^ кластер,
j иначе
После построения L частных решений можно сформировать согласованную матрицу различий 1 L
H = {h. j} H , = L E H,, (©l).
L l=i
где i,, = 1 ,., N . Величина , равна частоте классификации xi и x, в разные группы в наборе группи-
ровок ^ . Близкое к нулю значение величины означает, что данные объекты имеют большой шанс попадания в одну и ту же группу. Близкое к единице значение этой величины свидетельствует о том, что шанс оказаться в одной группе у объектов незначителен.
В нашем случае m = CCA(m, T, т), где число разбиений m е {mmin,mmin + 1,....mmin + L} , а объектами
классификации будут являться представители компонент связности Y(G1),..., Y(GS).
После вычисления согласованной матрицы различий для нахождения коллективного решения применим стандартный агломеративный метод построения дендрограммы, который в качестве входной информации использует попарные расстояния между объектами [10]. При этом расстояния между группами будем определять по принципу «средней связи», т. е. как среднее арифметическое попарных расстояний между объектами, входящими в группы. Процесс объединения продолжается до тех пор, пока расстояние
Td
между ближайшими группами не превысит заданное пороговое значение d , принадлежащее отрезку [ 0,1 ]. Этот метод позволяет выделять иерархическую структуру кластеров, которая упрощает процесс интерпретации результатов.
В соответствии с описанным выше методом разработан и программно реализован на языке Java ансамблевый алгоритм ECCA(mmin, L, T, т, Td). Для работы алгоритма необходимо задать значения пяти параметров: mmin, L, T, т, Td . Многочисленные экспериментальные исследования, проведенные на модельных и реальных данных, показали, что при использовании 8 - 10 элементов ансамбля получаемые результаты являются устойчивыми к выбору параметра сетки mmin . Параметр T оказывает слабое влияние на результат кластеризации. При обработке изображений параметр T выбирался равным 0,8, а порог шума т е{ 0;1} . Алгоритм ECCA позволяет получать иерархическую структуру данных. Проведенные
исследования показали [5], что параметр Td, задающий глубину иерархии, достаточно выбирать из множества {0, 0.1,..., 0.9}.
С помощью алгоритма ECCA производится сегментация изображения по спектральным признакам, результатами которой является картосхема изображения.
Сегментация картосхемы по текстурным признакам
Алгоритм кластеризации ECCA позволяет разбить изображение на однородные по спектральным признакам области без какого-либо учета их текстурных характеристик. Поэтому картосхема, получаемая после первого этапа обработки спутникового изображения высокого пространственного разрешения, представляет собой, как правило, чрезмерно раздробленную картину. На ней присутствует множество спектральных классов, которым нет соответствия среди реальных (информационных) классов, интересных пользователю.
На втором этапе обработки производится сегментация картосхемы по текстурным признакам и снижение ее раздробленности.
Для описания соответствующего этому этапу алгоритма введем необходимые обозначения.
Пусть картосхема, полученная в результате применения алгоритма ECCA, состоит из Nr x Nc = N элементов. Разобьем ее на квадратные фрагменты размером к x к элементов. Для простоты положим, что числа Nr и Nc кратны к . Каждому из таких фрагментов fij, i = 1,..., Nr / к; j = 1,..., Nc / к , сопоставим вектор ztj = (z1,...,zM) є RM . Здесь zj = Nj/к2 , где Nj - число элементов фрагмента j, относящихся к кластеру с номером l, l = !,•••, M.
Вектор zij фактически содержит в себе информацию о текстурной композиции элементов картосхемы в пределах фрагмента fj. Поэтому кластеризация векторов z^ приведет к сегментации картосхемы по текстуре.
В качестве меры близости между двумя такими векторами x = (x1,...,xM) и y = (y1,...,yM) є RM
M
примем метрику m(x, y) = 1 - ^ min(xi, yi).
i=1
Процедура кластеризации векторов z^ выполняется в три этапа. На первом этапе осуществляется поиск центров кластеров с помощью простого алгоритма DINA [11, с. 45]. На втором этапе полученные центры анализируются на предмет близости друг к другу. Если какие-либо два центра находятся на расстоянии меньше r (заданный радиус сферы), то они объединяются в один центр, соответствующий центру тяжести их внутренних точек. На третьем этапе осуществляется классификация векторов путем отнесения их к ближайшему центру, лежащему на расстоянии не превышающем r .
Кластеры, число элементов которых меньше порога т (задаваемого в процентах от общего числа фрагментов картосхемы), относятся к «шуму». Попавшие в «шум» фрагменты обрабатываются поэлементно. Каждый элемент «шумового» фрагмента относится к тому из соседних расклассифицированных фрагментов, в котором чаще встречаются представители спектрального кластера, содержащего обрабатываемый элемент.
Границы однородных текстурных областей, выделенных на картосхеме в ходе второго этапа сегментации, являются очень грубыми. Для их уточнения последовательно рассматриваются пограничные элементы всех выделенных текстурных областей и содержащие их фрагменты. Рассматриваемый пограничный элемент может либо остаться в текущей области, либо переместиться в смежную (если спектральный кластер, к которому он отнесен, чаще встречается в текстурном кластере смежного фрагмента, чем в содержащем его). После этого рассматриваемый элемент считается частью границы той области, к которой он отнесен. Процесс продолжается до тех пор, пока все пограничные элементы всех текстурных областей не будут рассмотрены.
Предложенный алгоритм текстурной сегментации картосхемы реализован на языке Java в виде модуля STSA с входными параметрами к, r, т .
Результаты экспериментальных исследований
В данном разделе представлены результаты сегментации модельных и реальных изображений с использованием модулей ECCA и STSA, подтверждающие вычислительную эффективность предложенного метода. Обработка проводилась на ПЭВМ с тактовой частотой 3.3 ГГц.
Во всех экспериментах при работе алгоритма ECCA использовался ансамбль из восьми элементов, а для модуля STSA значения параметров r и т полагались равными 0.4 и 1 соответственно.
На рисунке 1 представлены результаты сегментации цветного изображения размером 600 х 450 пикселей. При работе модуля STSA параметр k = 15 . Общее время обработки 0.3 с.
а б в
Рис. 1. а - исходное изображение; б, в - результаты сегментации по спектральным и спектрально-текстурным признакам соответственно
На рисунках 2 и 3 представлены примеры сегментации двух модельных RGB-изображений размером 800х800 и 858х572 пикселей соответственно. Картосхемы, полученные в ходе выполнения алгоритма ECCA, разбивались на фрагменты размером 12х12 пикселей для первого и 22х22 пикселей для второго изображения.
а б в
Рис. 2. а - исходное изображение; б, в - результаты сегментации до и после коррекции границ соответственно
Рис. 3. Исходное изображение (слева) и результаты сегментации (справа)
На рисунке 4 представлены результаты сегментации фрагмента снимка новосибирского Академгородка размером 2701х2458пикселей, полученного со спутника QшckBird 22 октября 2008 г. Для обработки использовались три спектральных канала (2, 3 и 4). Картосхема, полученная в результате применения ЕССА, разбивалась на фрагменты размером 12х12 пикселей. Общее время обработки 2.6 с.
Рис. 4. Слева - фрагмент синтезированного изображения QuickBird, построенного по 4, 3 и 2 каналам; справа - результаты сегментации. Выделено 6 кластеров (1 - песок; 2 - здания, сооружения, дороги; 3 -травянистая растительность; 4 - лиственные насаждения; 5 - сосновые насаждения; 6 - вода).
Заключение
В работе предложены двухэтапный метод и алгоритмы сегментации изображений высокого пространственного разрешения на основе совместного использования спектральных и текстурных признаков. Приведены результаты экспериментов на модельных и реальных изображениях, подтверждающие вычислительную эффективность предложенных алгоритмов и высокое качество получаемых картосхем.
Литература
1. Zhang, D. A review on automatic image annotation techniques [Текст] / D. Zhang, D. Zhang, M. Islam, G. Lu // Pattern Recognition. - 2012. - Vol. 45.
2. Dey, V. A review on image segmentation techniques with remote sensing perspective [Текст] / V. Dey, Y. Zhang, M. Zhong // ISPRS TC VII Symposium - 100 Years ISPRS, Vienna, Austria, July 5 - 7, 2010. IAPRS. - Vol. XXXVIII. - Part 7A.
3. Болсуновский, М. Л. Развитие систем ДЗЗ и информационно-аналитического обеспечения данными космической съемки: ближайшие перспективы [Текст] / М. Л. Болсуновский, Б. А. Дворкин // Геоинформатика. - 2010. - № 4.
4. Ilea, D. E. Image segmentation based on the integration of colour-texture descriptors. - A review [Текст] / D. E. Ilea, P. F. Whelan // Pattern Recognition. - 2011. - Vol. 44.
5. Ансамблевый алгоритм кластеризации больших массивов данных [Текст] / И. А. Пестунов, В. Б. Бери-ков, Е. А. Куликова, С. А. Рылов // Автометрия. - 2011. - Т. 47. - № 3.
6. Куликова, Е. А. Nonparametric clustering algorithm for large datasets [Текст] / Е. А. Куликова, И. А. Пестунов, Ю. Н. Синявский // Математические методы и распознавание образов: труды 14 конф.
- Изд-во MAKS Press, 2009.
7. Strehl, A. Clustering ensembles - a knowledge reuse framework for combining multiple partitions [Текст] / A. Strehl, J. Ghosh // The Journal of Machine Learning Research. - 2002. - Vol. 38.
8. Jain, A. K. Data clustering: 50 years beyond K-means [Текст] / A. K. Jain // Pattern Recognition Letters. -2010. - Vol. 31. - Is. 8.
9. Berikov, V. B. Construction of the ensemble of logical models in cluster analysis [Текст] / В. Б. Berikov // Lecture Notes in Artificial Intelligence. - 2009. - LNAI 5755.
10. Дуда, Р. П. Распознавание образов и анализ сцен. [Текст] / Р. Дуда, М. Харт. - М: Мир, 1976. - 559 с.
11. Загоруйко, Н. Г. Прикладные методы анализа данных и знаний [Текст] / Н. Г. Загоруйко. - Новосибирск: Изд-во Ин-та математики, 1999.
Информация об авторах:
Пестунов Игорь Алексеевич - кандидат физико-математических наук, доцент, заведующий лабораторией обработки данных Института вычислительных технологий СО РАН (ИВТ СО РАН), т. 8(383)334-91-55, pestunov@ict.nsc.ru
Pestunov Igor Alexeevich - Candidate of Physics and Mathematics, Associate Professor, HAed of the Laboratory for Data Processing at the Institute of Computational Technologies of the Siberian Branch of the RAS.
Рылов Сергей Александрович - магистрант Новосибирского государственного университета, т. 8(383)334-91-55, rylovs@mail.ru
Rylov Sergey Alexandrovich - Master’s programme student at Novosibirsk State University.
