ПРИКЛАДНАЯ МАТЕМАТИКА
УДК 517.9
Н. Б. Ампилова, Е. И. Петренко
ОБ ОЦЕНКЕ ЭНТРОПИИ СИМВОЛИЧЕСКОГО ОБРАЗА ДИНАМИЧЕСКОЙ СИСТЕМЫ
В статье рассматривается метод оценки энтропии символического образа динамической системы, который представляет собой ориентированный граф, построенный по системе. На графе строится инвариантная мера на простых циклах, определяющая стационарный процесс. Оценка энтропии такого стационарного процесса позволяет оценить энтропию символического образа. Реализован алгоритм построения такой меры, приведены результаты вычислений.
1. Символический образ динамической системы. Символический образ [1] динамической системы представляет собой конечный ориентированный граф, построенный по выбранному покрытию фазового пространства. Вершины графа соответствуют элементам покрытия. Для динамической системы /, покрытия {М*}, i = 1,...,п, на символическом образе строится дуга (г, Д, если /(М*) П А'^ ф 0. Символический образ характеризуется следующими параметрами:
- максимальный диаметр ячейки й\
- нижняя грань г - наименьшее из расстояний между образами ячеек и теми ячейками, с которыми образы не пересекаются;
- максимальный из диаметров образов ячеек д.
Справедлива следующая
Теорема 1 [2].
1. Если последовательность {гД есть путь на символическом образе С и х*. Е М(Д), тогда последовательность {хД есть е-траектория гомеоморфизма / для любого е > д + с! В частности, если последовательность {гД есть периодический путь на символическом образе, то последовательность {хД есть е-периодическая траектория.
2. Если последовательность {гД есть путь на символическом образе С, тогда существует последовательность {хД, х*. Е М(гД, которая есть е-траектория гомеоморфизма / для любого е > с!
3. Если последовательность {хД есть е-траектория гомеоморфизма /, е < г и х*. Е М(гД, то последовательность {гД есть допустимый путь на символическом образе С. В частности, если последовательность {хД есть е-периодическая траектория, то последовательность {гД есть периодический путь на символическом образе С.
© Н. Б. Ампилова, Е. И. Петренко, 2008
Для получения более точного приближения к системе строится последовательность символических образов, каждый из которых соответствует последовательному подразбиению начального покрытия заданной области фазового пространства. Символический образ является конечным приближением системы и отражает ее глобальную структуру.
Энтропия оценивает сложность поведения динамической системы, измеряя рост числа траекторий. Используя символический образ динамической системы, мы оцениваем рост числа допустимых путей и получаем энтропию символического образа (аналог топологической энтропии), которая, в свою очередь, позволяет оценить энтропию исходной системы [2].
На символическом образе можно построить инвариантную меру, применяя технику стационарных процессов. Энтропия, вычисленная по такой мере (аналог метрической энтропии), дает характеристику асимптотического поведения путей на символическом образе и тем самым траекторий исходной системы.
2. Марковские цепи на графе. Рассмотрим множество 5 состояний некоторого процесса и взвешенный мультиграф С = (V, Е) с множеством вершин, соответствующим множеству состояний. Каждому состоянию приписывается некоторая вероятность. Множество вершин графа С соответствует множеству состояний с положительными вероятностями, т. е. V = {I Е > 0}. Веса дуг - вероятности переходов из
состояния в состояние.
Пусть г(е), 4(е) обозначают соответственно начало и конец дуги е. Обозначим /х(г(е)) вероятность находиться в начале дуги е, ц(е|г(е)) - вероятность перехода в вершину 4(е) при условии, что мы находимся в вершине г(е).
Определение 1. Пусть задан граф О = (У,Е) и Е[ обозначает множество дуг, выходящих из вершины I. Марковская цепь на графе определяется начальным распределением вероятностей ц(1) > 0,1 Е V и условными вероятностями ц(е|г(е)) > 0 для ее Е так, что У,, ( //(/) = 1; Ее61?1 м(е1-0 = IV! Е V.
Вероятность находиться на дуге е
ц(е) = ц(г(е))ц(е\г(е)). (1)
Вероятность пути ехе-2 ■ ■ ■ еп определяется следующим образом:
Ме1е2 ■ ■ -е„) = М*(е1))Ме1|*(е1))Ме2|*(е2)) ■ ■ -Меп|*(е„))- (2)
Из соотношения (1) следует, что
^2 м(е) = м(*(е)) ^2 м(е|г(е)) =м(*(е)),
ееЕце) ееЕце)
т. е. вероятность находиться в вершине равна сумме вероятностей (мер) выходящих из нее дуг.
Начальные распределения вероятностей образуют вектор-строку р, так что р/ = 1л(1). Величины ц(е\1) образуют стохастическую матрицу переходов Р с элементами
Ри = £ ц(е\1),
ее£/
где и/ обозначает множество дуг из I в Л.
3. Стационарные процессы на графах и их энтропия. Пусть для марковской цепи на графе С = (V, Е) заданы вектор начального распределения вероятностей р и матрица переходов Р.
Определение 2. Марковская цепь называется стационарной, если выполняется равенство рР = р.
Нетрудно заметить, что условие стационарности означает, что для каждой вершины сумма вероятностей (мер) по всем входящим в нее дугам равна сумме вероятностей по всем выходящим из нее дугам. Вероятность пути длины п на С найдем, следуя (2). Для каждого п составим вероятностный вектор из вероятностей путей длины п и вычислим его энтропию к{р^). Энтропия стационарного процесса /л определяется [3]
Известно [4], что величину, определяемую (3), можно вычислить непосредственно, в терминах начального распределения вероятностей и условных вероятностей:
4. Построение меры ц.
Определение 3. Ориентированный граф С называется сильно связным, если для любой упорядоченной пары вершин (I, 7) существует путь из I в ,1.
Ориентированный граф можно разбить на непересекающиеся компоненты сильной связности, связи между которыми изображаются с помощью диаграммы порядка. Меру можно задать либо на каждой компоненте, а на всем графе определить как выпуклую комбинацию мер на компонентах, либо задавать ее на всем графе. Для простоты иллюстраций опишем построение на компоненте связности.
Пусть С?т - компонента связности графа С. Найдем все простые циклы в Ст. Обозначим их щ,..., ии- Пусть 1{ обозначает число дуг в цикле щ. Припишем каждой дуге в цикле щ величину, равную 1/^. Для каждого из циклов выберем некоторый коэффициент (вес) г>£, так, чтобы Хл=1 = Тогда каждой дуге в цикле щ можно приписать величину Для дуги, принадлежащей в циклам, ее мера определяется
как Х^=1 'иь /1г, ■ Дугам, не входящим ни в какие циклы, приписывается нулевая вероятность.
Мера вершины определяется как сумма мер всех выходящих из нее дуг. Условная вероятность дуги рассчитывается из соотношения /х(е) = /х(г(е))/х(е|г(е)). Заметим, что если из вершины выходит только одна дуга, то ее условная вероятность равна 1.
Нетрудно проверить, что такая мера задает стационарный процесс и является инвариантной относительно оператора циклической подстановки в циклах.
5. Возможные конфигурации циклов. По построению символический образ не является мультиграфом. На нем возможны следующие расположения циклов.
1. Циклы не имеют пересечений (рис. 1).
Нетрудно заметить, что в этом случае энтропия равна нулю, поскольку все условные меры дуг равны 1 (из каждой вершины выходит одна дуга).
2. Циклы имеют общие точки (дуги или вершины) попарно (рис. 2).
Эти стационарные процессы имеют одинаковую энтропию, которая вычисляется по формуле
как
(3)
Мм) = - м(*(е))м(е|*(е))1о§(^(е|г(е))).
еЄЕ(С)
(4)
(УьЬ) (V2.ll) (Уэ,1з) (У4,Ц)
ЛЫ = -Е( |ю6
гРп ^г+1 ^г+1 Рп+1 ^
щ/и + ^+1/^+1 и+1 щ/и + ^+1/^+1 /
где 8 - ЧИСЛО ЦИКЛОВ.
3. Вложенные циклы (рис. 3).
Рассмотрим конструкцию, состоящую из двух наборов вложенных циклов (С и Г>), имеющих общую дугу АВ (или, что равносильно, общую вершину). Предположим, что количество ребер в циклах набора С меняется от то до 2, а в наборе Б от п до 2 с шагом
— 1, т. е. присутствуют все возможные циклы. Обозначим веса циклов в части С за
</,„..... 1/2- в части I) за сп.так что Хл=2 "I" Е]=2 из = I • Меры дуг в этих
циклах определяются как ит/т,ит-1/(т — 1),... ,1^2/2 и г>п/п, г>п_1/(п — 1),..., г?2/2. Введем обозначения Бг = ип/п + г>п_1/(п — 1) + ... + г>г/2, £2 = ит/т + ит-1/(т — 1) + ... + иг/2. Тогда ц(АВ) = Бг + 5г, и энтропия может быть вычислена по формуле
С
В
Рис. 3. Вложенная конструкция.
Рис. 4- Два пересекающихся цикла.
6. Пример. Рассмотрим граф, состоящий из двух циклов - с и с1 (рис. 4). Множество состояний 5 = {1,2,3,5,6,7,8}.
Циклы имеют общее ребро - дугу 2 > 3. Положим г>с = I /2. г./ = 1/2. Для цикла с величина 1/1,1 равна 1/4, для цикла й - 1/5. Меры (вероятности) дуг в первом цикле равны 1/8, кроме дуги 2 —>• 3, меры дуг второго цикла - 1/10. Мера дуги 2 —>■ 3 равна сумме мер, которые эта дуга получает в каждом содержащем ее цикле, т. е. 1/8+1/10 = 9/40. Поскольку /л(е) = /г(г(е))/г(е|г(е)), то, зная меры дуг и меры вершин, вычисляем условные меры (вероятности) дуг. Нетрудно заметить, что почти все условные вероятности равны 1, кроме вершины 3. Здесь /1(3) = 9/40, /1(3 —> 5) = 1/8 : 9/40 = 5/9, ц(3 -► 8) = 1/10 : 9/40 = 4/9.
Вектор р = (1/8,9/40,9/40,1/8,1/10,1/10,1/10),
Р.
и условие стационарности выполняется. Отсюда
Н(р) = -9/40(4/91оё(4/9) +5/91оё(5/9)) = — (1/81оё(5/9) + 1/101оё(4/9)) и 0.22299211.
7. Зависимость энтропии от выбора г>*. Вопрос о выборе весов г>£, при котором энтропия принимает наибольшее значение, приводит к решению задачи на нахождение условного экстремума. Для приведенного в п. 6 примера нужно найти условный экстремум функции
Ф(и1,и2, А) = Ыу1,у2) + Хф(у 1,и2), где ^>(г>1,г>2) = г>1 + г>2 — 1. Это приводит к соотношению
/ 0 1 0 0 0 0 0 \
0 0 1 0 0 0 0
0 0 0 5/9 0 0 4/9
1 0 0 0 0 0 0
0 1 0 0 0 0 0
0 0 0 0 1 0 0
V 0 0 0 0 0 1 0 /
из которого получаем уравнение для v\
!L \ 2 / 1-^1
i______I = _____________ь______
1 — V\ I \ Vl j 1—V\
(0 1 0 0 0 0 0 \
0 0 1 0 0 0 0
0 0 0 1 0 0 1
1 0 0 0 0 0 0
0 1 0 0 0 0 0
0 0 0 0 1 0 0
10 0 0 0 0 1 0 )
При h = 4,1,2 = 5 уравнение (5) имеет вид 2869wf + 3840«i — 1024 — 5120wf+ 2560wf = 0, откуда v% ps 0.4828949891, г>2 » 0.517105011, h{p) ps 0.22318030. Таким образом, в этом примере при назначении циклу с большим периодом большего веса энтропия возрастает. К сожалению, для числа циклов более 3 нахождение условного экстремума становится затруднительным. Этот вопрос будет рассмотрен при обсуждении оценки энтропии для отображения Хенона.
Известно [3], что энтропия любого стационарного процесса (3) на сильно связном графе G может быть оценена сверху числом log Л, где Л - максимальное собственное число его матрицы смежности Аа■ В рассмотренном примере
Ag
Решение характеристического уравнения Л7 — Л3 — Л2 — 1 = 0 дает вещественный корень Атах к, 1.236505703, что, в свою очередь, приводит к оценке h(p) 0.3062688937.
8. Описание алгоритма нахождения циклов. Задача построения инвариантной меры на символическом образе сводится к задаче перебора простых циклов на графе. Мы предлагаем приближенный алгоритм, который имеет сложность 0(т + п), где пит- число вершин и дуг графа соответственно.
Пусть G = (V,E) - ориентированный граф. Обозначим через е = (г —> j) дугу между вершинами i.j.
Определение 4. Будем говорить, что вершина с является потомком вершины р (или с является предком вершины р), если р - родительская вершина с, либо если родитель с есть потомок р.
Алгоритм использует поиск в ширину и обходит все компоненты связности. Одновременно строится дерево поиска, где существование ребра р —> с означает, что среди вершин, достижимых из р в графе G, была рассмотрена вершина с. Заметим, что для работы алгоритма достаточно хранить только ссылку из вершины-потомка (с) на вершину-родителя (р).
Узел дерева поиска представляется парой (вершина, указатель на родительский узел дерева поиска), ref обозначает операцию вычисления указателя на пару, null -пустой указатель, маркер вершины дерева поиска. Очередь таких пар обозначается q.
Схема алгоритма
1. Выберем произвольную вершину г.
2. Положим в q пару (i,null).
3. Пока очередь q не пуста, выберем вершину (п,р) из q.
4. Если вершина п уже была рассмотрена, пропустить ее. Перейти на шаг 3.
5. Иначе для всех ребер п —>■ то Е Е выполнить
а) если вершина то уже рассмотрена и является потомком вершины п, то найден цикл - путь в дереве поиска от вершины п к вершине т. Вершины найденного цикла можно перебрать при помощи ссылок на родительскую вершину пар представления дерева поиска;
б) иначе если то не является потомком п, положить (m,ref(n,p)) в очередь.
Для проверки того, была ли вершина обработана алгоритмом, используется хэш-множество. Каждая вершина рассматривается не более одного раза.
Замечание. Отметим, что некоторые простые циклы могут быть не найдены. Например, для полного графа из трех вершин алгоритм отработает следующим образом:
1. Начинаем алгоритм с вершины 1. Узел дерева поиска п\ = (1, null). Из вершины 1 находим ребра в вершины 3 - добавляем в очередь вершину (3, re/ni); 1 - найден цикл; 2 - добавляем в очередь вершину (2, re/ni).
2. Рассматриваем следующую вершину в очереди - п-2 = (3, re/ni). Из вершины 3 находим ребра в вершины 3 - найден цикл; 1 - найден цикл; 2 - добавляем в очередь вершину (2,re/ri2).
3. Рассматриваем следующую вершину в очереди - пз = (2, re/ni). Из вершины 2 находим ребра в вершины 3 - вершина уже рассмотрена, но цикл не найден; добавляем в очередь (3, re/пз); 1 - найден цикл; 2 - найден цикл.
4. Рассматриваем следующую вершину в очереди - гц = (2,re/ri2). Вершина 2 была рассмотрена алгоритмом ранее.
5. Рассматриваем следующую вершину в очереди - щ = (3,re/пз). Вершина 3 была рассмотрена алгоритмом ранее.
6. Очередь пуста. Алгоритм завершен.
В результате работы алгоритма найдены циклы (1), (2), (3), (1 —> 2 —> 1), (1 —> 3 —> 1). Не найден цикл (2 —> 3 —> 2). Эксперименты показывают, что могут не найтись некоторые циклы, соответствующие вложенным конструкциям. Таким образом, некоторым дугам не будут приписаны меры.
9. Алгоритм построения меры на графе. Заметим, что формула вычисления энтропии стационарного процесса (4) может быть переписана в более удобном для вычислении виде
Нр) = Ем№)) log(/i(i(e))) - Е м(е) logMe)). (6)
iev евЕ
Определение 5. Итератором будем называть объект, который предоставляет последовательный доступ к элементам некоторого набора объектов [5].
Понятие итератора дает возможность абстрагироваться как от способа представления, так и от метода получения перебираемых объектов. В данном случае можно рассмотреть описанный выше алгоритм как итератор, который позволяет перебрать некоторые простые циклы графа. Точнее говоря, предлагаемый алгоритм работает с итератором дуг. Циклы обрабатываются по мере их нахождения, т. е. мера назначается очередному найденному циклу.
Схема алгоритма
В начале работы всем дугам приписывается значение 0. В процессе работы алгоритм приписывает им значения описанным ранее способом с выбором весов V{ с помощью функции g : No —> И.
1. Для каждого найденного итератором цикла I выполним следующее:
а) вычислим г = 9<,щ^;
б) для всех дуг е цикла I выполним //(< ) := //(< ) + г.
2. Перебор циклов завершен. Проведем нормировку /і, т. е. ш := ^ м(е): затем
еЄЕ
для каждой дуги е /і(е) := Вычислим энтропию по построенной мере, следуя (6). Были рассмотрены такие способы выбора функции д:
1) 5(0 = 1- При этом выборе каждой дуге в цикле назначается мера 1 После просмотра циклов її,..., 1п величина ш равна п и после нормировки мера дуги становится равной 1/іі * (1 /п). Эта ситуация соответствует выбору одинаковых
2) д(1) = I. Мера цикла іі до нормировки равна 1^, после нормировки она становится
і к)- Таким образом, большему циклу назначается больший вес;
3) д(1) = I2. Мера цикла выбирается как Іг2/ЕГ= і Ь2)-
10. Численный эксперимент. Алгоритм был применен к получению оценки энтропии отображения Хенона
при а = 1.4, Ь = 0.3.
Известно [6], что при заданных параметрах у системы (7) существует инвариантное множество в области [—1.5,1.5] х [—1.5,1.5]. Рассмотрим область D = [10,10] х [—10,10]. Полученный на 13-м шаге при помощи линейного метода [7] символический образ состоит из двух компонент сильной связности с количеством вершин 572 427 и 8 соответ-
РФВОХХИА
ъ 1 JtSCJtlJtlU«
Приведем результаты вычислений для различных функций д :
Постоянная...................0.460142
Линейная ....................0.267599
Квадратичная ................0.213950
Заметим, что применение описанных стратегий назначения большему циклу большего веса не дает увеличения энтропии. Для отображения Хенона известна оценка энтропии сверху - 0.4651 [2]. Таким образом, по экспериментальным данным оценка энтропии символического образа как энтропии стационарного процесса по инвариантной мере дает оценку снизу для исходного отображения.
Summary
Ampilova N. В., Petrenko Е. I. On estimation of the entropy of the symbolic image of a dynamical system.
The method of the construction of the invariant measure on an oriented graph is considered. The algorithm of the estimation of the entropy of a dynamical system using its finite approximation -symbolic image is developed and implemented. The numerical results are given.
Литература
1. Осипенко Г. С. О символическом образе динамической системы// Граничные задачи / Под ред. В. А. Алексеева. Пермь, 1983. С. 101-105.
2. Осипенко Г. С., Ампилова Н. Б. Введение в символический анализ динамических систем. СПб.: Изд-во С.-Петерб. ун-та, 2005. 236 с.
1 — ах + by
х
(7)
3. Lind D., Marcus B. An introduction to symbolic dynamics and coding. Cambridge: Cambridge Univ. Press, 1995. 478 p.
4. Petersen K. Ergodic theory. Cambridge: Cambridge Univ. Press, 1989. 370 p.
5. Gamma Erich, Helm Richard, Johnson Ralph, Vlissides John. Design patterns elements of reusable Object-Oriented Software. Addison-Wesley, 1995. ISBN 0-201-63361-2.
6. Henon M. A two-dimensional mapping with a strange attractor // Comm. Math. Phys. 1976. Vol. 50. P. 69-77.
7. Петренко E. И. Разработка и реализация алгоритмов построения символического образа // Электр, журн. Дифф. уравнения и процессы управления, http://www.neva.ru/journal 3/2007.
Статья рекомендована к печати проф. Л. А. Петросяном.
Статья принята к печати 21 февраля 2008 г.