Стационарные процессы на графах и анализ изображений Текст научной статьи по специальности «Математика»

Ампилова Наталья Борисовна

УДК 51-76

СТАЦИОНАРНЫЕ ПРОЦЕССЫ НА ГРАФАХ И АНАЛИЗ ИЗОБРАЖЕНИЙ

Аннотация

В работе описан метод классификации изображений, иллюстрирующих процесс распространения вещества в среде, в частности изменения, происходящие при действии сверхмалых доз лекарственных препаратов. Изображение рассматривается как решетка пикселей заданной интенсивности. По изображению строится ориентированный граф, так что каждый узел соединен с N соседями. Всем выходящим из узла дугам приписывается значение интенсивности узла, деленное на N (для точек границы - на N - 1). Построенный поток нормируется. Для полученной таким образом марковской цепи методом Шелейховского-Брэгмана строится стационарное распределение, которое максимизирует взвешенную энтропию. Именно значение взвешенной энтропии выбирается как классификационный признак при анализе изображений, соответствующих различным дозам вещества. Приведены результаты численных экспериментов.

Ключевые слова: Марковские цепи, стационарный процесс на графе, аппроксимация инвариантной меры, динамические системы, максимизация взвешенной энтропии.

1. ВВЕДЕНИЕ

Современные методы исследования и классификации изображений часто используют текстурный, морфологический и фрактальный анализ, а также всевозможные сочетания этих методов, поскольку одной характеристики бывает недостаточно. Используются также методы нейронных сетей, вейвлет-анализ. Как правило, для каждого класса изображений в процессе анализа удается найти сочетание методов, которые наиболее разумно использовать в данном случае. Так, методы фрактального анализа успешно используются для классификации текстур наряду с уже ставшими традиционными статистическими признаками, введенными Р. Хараликом.

Анализ изображений биомедицинских препаратов широко использует текстурный

© Ампилова Н.Б., 2013

анализ [9] (использование цвета как определяющего признака, выделение тексто-нов - специальных визуальных элементов), морфологические методы (анализ изображений ультразвуковой диагностики), а также фрактальный и мультифрактальный анализы. Применение фрактальных методов к анализу различных классов биомедицинских препаратов было описано в [2]. В работе [8] мультифрактальный спектр как классификационный признак был использован для анализа изображений крови, гистологических препаратов и соединительной ткани. В данной работе предлагается использовать некоторую энтропийную характеристику для анализа и классификациии изображений диффузионных процессов. Такие процессы наблюдаются при действии лекарственных препаратов в различных дозах. Наш метод основан на представлении исследуемого изображения в виде ориентированного гра-

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

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

2. ОСНОВНЫЕ ОПРЕДЕЛЕНИЯ

Рассмотрим множество 5 состояний некоторого процесса и взвешенный граф О = (V,Е) с множеством вершин, соответствующим множеству состояний. Каждому состоянию приписывается некоторая вероятность. Множество вершин графа О соответствует множеству состояний с положительными вероятностями, то есть V = {Iе т(1) > 0}. Дугам тоже приписываются некоторые значения, которые соответствуют вероятностям перехода из состояния в состояние.

Пусть ¡(в), ?(в) обозначают соответственно начало и конец дуги в. Обозначим ¡л(1 (в)) вероятность находиться в начале дуги в, ¡л(в\1 (в)) - вероятность перехода в вершину I (в) при условии, что мы находимся в вершине ¡(в).

Определение 1. Пусть задан граф О = (V, Е) и Е1 обозначает множество дуг, выходящих из вершины I. Марковская цепь на графе определяется начальным распределением вероятностей т(1) > 0, I е V

и условными вероятностями m(e\i(e)) > 0 для e е E так, что

h^v m(D =1 ; le€El m(e\D =1 "1 ev.

Вероятность находиться на дуге e (мера дуги) определяется как m(e) = m(i(e))m(e\i(e)). Вероятность пути e^2 ...en определяется следующим образом:

m(eie2 •••en) = m(i (ei))m(ei\i (ei))

m(e2\i (e2)) ••• m(en\i (el)). (1)

Заметим, что из определения меры дуги следует, что

^m(e) = m(i(e)) \ i(e)) = m(i(e)),

ee Ei(e) ee Ei(e)

то есть мера вершины равна сумме мер выходящих из нее дуг.

Начальные распределения вероятностей образуют вектор-строку p, так что pj = m (J ) Величины m(e\J) образуют стохастическую матрицу переходов P с элементами

Pjj = ^eeEj m(e \ 1 ), где EJj обозначает множество дуг из J в J•

Определение 2. Марковская цепъ называется стационарной, если вытолняет-ся равенство pP = p^

Нетрудно заметить, что условие стационарности означает, что для каждой вершины сумма вероятностей (мер) по всем входящим в нее дугам равна сумме вероятностей по всем выходящим из нее дугам.

Полезной характеристикой стационарного процесса является энтропия, которая вычисляется как предельное значение энтро-

" 1 i ( ")\ (n)

пий h(pm ) вероятностных векторов pm ,

составленных из мер всевозможных путей длины n на G, определяемых по формуле (1).

Энтропия стационарного процесса m определяется [10] как

h(m) = lim-h (p(;y). (2)

n

Известно [12], что величину, определяемую по формуле (2), можно вычислить непосредственно в терминах начального распределения вероятностей и условных вероятностей:

h(m) = - | i(e))log(m(e | i(e))),

^ (3)

или в эквивалентной форме

h(m) = -Xm(i (e))iog(m(i (e))) -

ieV

-Xm(e)log(m(e)). (4)

eeE

3. ПОСТАНОВКА ЗАДАЧИ

Для изображений, показывающих процесс распространения вещества, можно предложить следующий метод анализа, основанный на его характеристиках яркости (интенсивности). Рассмотрим изображение как решетку пикселей {x (i, j)}, i = 1...M, j = 1...N. Каждой точке решетки сопоставляется значение интенсивности пикселя -I(i, j), где I(i, j) е [0, L]. Как правило, L = 255. (Заметим, что в дальнейшем нулевые значения яркостей мы будем заменять некоторым малым значением). Для каждого узла рассмотрим его ближайших соседей (окрестность). В качестве окрестности можно выбирать 4 или 8 соседей. Отношение соседства изображается так: x(i, j) и x (k, l) соседи, если между ними есть ребро. При таком построении получаем неориентированный граф. Преобразуем его в ориентированный, вводя дуги x(i, j) ® x(k, l) и x(k, l) ® x (i, j) вместо соответствующего ребра. Присвоим дугам значения следующим образом. Для всех точек решетки (кроме точек границы) рассматриваем узел x(i, j) и его соседей x(i, j + 1), x(i, j - 1), x (i - 1, j), x (i + 1, j) и всем выходящим дугам x(i, j) ® x(i, j + 1), x(i, j) ® x(i, j - 1), x (i, j) ® x (i - 1, j), x(i, j) ® x (i + 1, j) присвоим значение I (i, j)/4. Для точек границы присваиваем исходящим дугам I (i, j)/3. Таким образом, мера каждого узла равна сумме мер выходящих дуг. Если добавить условие нормировки V I(i, j) = 1, то мы полу-1,j

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

начальном распределении найденное стационарное распределение максимизирует взвешенную энтропию, определяемую как

I С, ])

V г, ju(i, j)ln

и О', ])

4. МЕТОД РЕШЕНИЯ

У задачи построения стационарного потока на графе очень интересная история и неожиданные области применения.

Рассмотрим следующую задачу линейного программирования.

Максимизировать функцию

V

Xj ln-

1, j у

при ограничениях

Xj ь 0;

X х] = а, X х] =ъ],

]=1 1=1 т т

Xа =Хъ];а,ъ} >0; Р]-0, XхУ- =1.

Эта задача возникала много лет назад в связи со специфическими расчетами градостроительного характера. Прогнозировалось устойчивое расселение жителей города по районам в зависимости от мест приложения труда. Иначе говоря, требовалось найти допустимое решение некоторой специальной транспортной задачи. Здесь хц обозначает число людей, живущих на территории г е I и работающих на предприятии; е /, а вышеприведенные соотношения

Хе1Х] = Ъ] , 1 е (*)

Х]^хи =а , г е 1 (**)

задают транспортные ограничения. Наборы {а,.} и {Ь;} заданы. Для вершины г эти ограничения означают, что сумма значений по выходящим из вершины дугам равна аг, а сумма по входящим в нее равна Ь{.

Ленинградский архитектор Г.В. Шелей-ховский [6] предложил искать решение, удовлетворяющее этим ограничениям, методом балансировки, который состоит в следующем. Выберем некоторый начальный набор х0 = р ] и, умножая строки матрицы X0 на подходящие множители, добьемся выполнения равенств (*). Обозначим полученную матрицу через X1. Теперь, умножая

X

столбцы этой матрицы на подходящие множители, получим матрицу X2, которая будет удовлетворять условиям (**), но не (*). Получаем последовательность матриц, которая довольно быстро сходится к некоторому пределу.

И.В. Романовский показал, что если этот метод сходится, то предельная матрица доставляет максимум взвешенной энтропии

X х V 1п р = X х V 1п р V - X х V 1п хг

геI, jеJ у геI, jеJ геI, ]е J

Такое распределение в некоторых случаях оказывается предпочтительным, поскольку значения присваиваются всем дугам графа.

Сходимость метода Шелейховского была доказана Л. Брэгманом в [4]. Вскоре Л. Брэгман разработал более общий метод для решения задач выпуклого программирования [5], который для задач с линейными ограничениями совпадает с методом Г. Шелейховского. В работе [3] было доказано, что этот общий метод:

1) сходится не медленнее, чем геометрическая прогрессия,

2) результат не зависит от порядка перебора строк и столбцов,

3) полученное распределение зависит от начальных условий.

Интересным частным случаем этой задачи является ситуация, когда значения а{ и bj совпадают. Тогда мы получаем так называемое стационарное распределение, когда для любой вершины сумма мер выходящих дуг равна сумме мер входящих.

Стационарный поток на графе оказался полезным при решении задач исследования динамических систем. В работе [11] Г. Осипенко описал применение алгоритма балансировки к задаче построения приближения инвариантной меры динамической системы с помощью графа символического образа. Практическая реализация алгоритма построения приближения к инвариантной мере была описана в [7]. Стационарное распределение (или сбалансированный граф) возникает и в задачах теории управления, например в задаче о достижении консенсуса группой взаимодействующих агентов [1].

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

Мы решаем задачу анализа изображения с помощью построения стационарного распределения на графе, построенном по изображению, и вычислении взвешенной энтропии. Кроме того, для полученного стационарного распределения мы вычисляем его энтропию по формуле (4).

5. АЛГОРИТМ ПОСТРОЕНИЯ

СТАЦИОНАРНОГО РАСПРЕДЕЛЕНИЯ

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

bвg(n) = {в е Е, в = (п, ]), ] е V;} (5)

впй(п) = {в е Е, в = (¡,п), I е V.} (6)

тош(п) = X („) (7)

тп(п) = X ее^ („) т(е) (8)

Схема алгоритма следующая.

• Назначить всем дугам произвольный ненулевой вес. Положить т(в) = 1 для всех в е Е.

• Построить приоритетную очередь Q из вершин графа, принимая в качестве приоритета #(п) вершины п ее дисбаланс

я(п) = \тош(п)- т>пШ

• В цикле: извлечь очередную вершину п из Q.

- Если #(п) < е, то завершить обработку п и выйти из цикла. (В силу приоритетности извлекаемых вершин, это неравенство выполнено для всех оставшихся элементов.)

- Иначе,

* вычислить тош(п)-min(n),

г=4тош («Vтт (п) •

* Меры всех дуг, исходящих из п, разделить на у.

* Меры всех дуг, входящих в п, умножить на у.

* Если одно из значений №оШ(п), тп(п), у либо достаточно велико, либо достаточно мало, выполнить общую нормировку потока.

• Выполнить общую нормировку потока. Алгоритм завершен.

6. ЧИСЛЕННЫЕ ЭКСПЕРИМЕНТЫ

Результаты работы алгоритма для тестовых черно-белых изображений приведены на рис. 1 и 2. Рассматривались изображения размером 20 х 20 пикселей, при иллюстрации они были увеличены для удобства восприятия, так что каждый квадратик соответствует пикселю изображения.

Алгоритм был применен также для классификации изображений фармакологических растворов серебра различной концентрации, полученных с помощью атомно-си-лового микроскопа. Цветные изображения сначала преобразовывались к полутоновым. Эксперименты проводились для изображений трех разных классов, соответствующих различным концентрациям. Для иллюстрации выбрано по одному представителю каждого класса. Примеры расположены по порядку возрастания доз серебра в растворе -нулевая, малая и большая концентрации. (рис. 3, 4 и 5 соответственно.) Для каждой пары изображений слева расположено исходное изображение, а справа изображение, соответствующее полученному стационарному распределению. Нужно заметить, что различия в картинках не очень заметны глазу, как для тестовых изображений, поэтому при вычислениях формировались так называемые «карты отличий», где отмечались несовпадающие точки.

Далее по тексту ка обозначает взвешенную энтропию, а к5 - энтропию стационарного распределения.

Вычислительный эксперимент проведен Е. Петренко.

XX

Рис. 1. Тестовый пример с диагоналями. к„ = 0,215955; к = 2,699519

Рис. 2. Тестовый пример с черно-белыми полосами. ка = 0,515709; к5 = 2,336137

Рис. 3. Нулевая доза вещества в растворе. й = 0.000002; к = 2,997057

а ' я '

Рис. 4. Малая доза вещества в растворе. й = 0.000044; к = 2,987622

а ' я '

Рис. 5. Большая доза вещества в растворе. й = 0.000043; к = 2,991464

а 'я '

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

Эксперименты показывают, что взвешенная энтропия может рассматриваться как классификационный признак, отличающий изображения по степени концентрации раствора, а именно растворы без серебра и

Литература

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

1. Амелина Н. О., Фрадков А.Л. Метод усредненных моделей в задаче достижения консенсуса / Стохастическая оптимизация в информатике (под ред. О.Н. Граничит), 2012. Т. 8, № 1. С. 3-39.

2. Ампилова Н.Б., СоловьевИ.П. Алгоритмы фрактального анализа изображений // Компьютерные инструменты в образовании, 2012. № 2. С. 19-24.

3. БрэгманЛ.М. Релаксационный метод нахождения общей точки выпуклых множеств и его применение / Дисс. на соиск. ст. к. ф.-м. н. Л., 1966.

4. Брэгман Л.М. Доказательство сходимости метода Г.В. Шелейховского для задачи с транспортными ограничениями // Журнал вычисл. математики и матем. физики, 1967. Т. 7, № 1. С. 147-156.

5. Брэгман Л.М. Релаксационный метод нахождения общей точки выпуклых множеств и его применение для решения задач выпуклого программирования // Журнал вычисл. математики и матем. физики, 1967. Т. 7, № 3. С. 620-631.

6. Шелейховский Г. В. Композиция городского плана как проблема транспорта / Рукопись, 1946.

7. Ampilova N., Petrenko E. On Application of Balance Method to Approximation of Invariant Measures of a Dynamical System // Эл. журнал «Дифференциальные уравнения и процессы управления», 2011. № 1. С. 72-84 / http://www.math.spbu.ru/diffiournal/RU /numbers/2011.1 (дата обращения 15.02.2013).

8. Ampilova N., SolovievI., Shupletzov Y. Multifractal spectrum as a classification sign for biomedical reparations images / Proc. 7 Int. Conf. CEMA12, 8-10 Nov. 2012. Athens, Greece. P. 98-101.

9. The handbook of pattern recognition and computer vision (Eds.C.H. Chen, L.F. Pay, P. Wang) // World Sci.Publ.Co., 1998. P. 207-248.

10. LindD., Marcus B. An introduction to symbolic dynamics and coding. New York, 1995.

11. Osipenko G. Symbolic images and invariant measures of dynamical systems // Ergod.Th.& Dynam.Sys. (2010)/ Vol. 30. P. 1217-1237.

12. Petersen K. Ergodic Theory. Cambridge Univ.Press. Cambridge, 1989.

Abstract

The method of a classification of images concerning to a substance propagation process is proposed. The image is considered as a lattice formed by pixels of given intensity. Then an oriented graph corresponding to the image is constructed in the following way: every vertex (pixel) is connected with N neighbours. For a given vertex all outcoming edges have a value (pixel intensity/N), for boundary vertex - pixel intensity/(N - 1). The constructed flow is normed. For obtained markov chain by the Sheleihovsky-Bregman method a stationary distribution is constructed, which maximizes weighted entropy. It is weighted entropy that is used as a classifying sign when images with different doses of a substance are analyzed. The results of numerical experiments are given.

Keywords: Markov chain, stationary process on graph, approximation of an invariant measure, dynamical systems, maximization of weighted entropy.

Ампилова Наталья Борисовна, кандидат физико-математических наук, доцент кафедры информатики математико-механи ческого факультета СПбГУ, атрИога@та^. зрЬи. ги

© Наши авторы, 2013. Our authors, 2013.