CC BY
12
Известия АлтГУ. Математика и механика. 2019. № 4 (108)
УДК 517.938
0 дискретизации фазовых портретов д инамических систем*
Н.Е. Кириллова1, Л.С. Минушкина2
1 Институт математики им. С.Л. Соболева СО РАН (Новосибирск, Россия)
2 Новосибирский государственный университет (Новосибирск, Россия)
On Discretization of Phase Portraits of Dynamical Systems
N.E. Kirillova1, L.S. Minushkina2
1 Sobolev Institute of Mathematics Siberian Branch of the Russian Academy of Sciences (Novosibirsk, Russia)
2 Novosibirsk State University (Novosibirsk, Russia)
Рассматриваются динамические системы, моделирующие функционирование кольцевых генных сетей, в которых скорость изменения концентрации любого вещества зависит монотонно от концентрации вещества, предшествующего данному. Ранее исследовались вопросы существования, единственности, а также устойчивости периодических решений для динамических систем, уравнения которых содержат монотонно убывающие и монотонно возрастающие функции, и для описания фазовых портретов таких систем проводилось построение положительно инвариантной области. В случае больших размерностей возникает необходимость сужения построенной области. Для этого авторами был разработан алгоритм, основанный на задаче об «укладывании змеи в ящик» (snake-in-the-box problem) из теории графов и применимый для динамических систем различного вида и различных размерностей. Алгоритм позволяет построить диаграмму переходов, соответствующую сдвигам вдоль траекторий динамических систем. Приведены примеры использования алгоритма в случае шестимерной системы со ступенчатыми функциями и в случае десятимерной системы с гладкими функциями. Результаты построения были применены в исследовании существования циклов у этих систем. Ключевые слова: динамические системы, обратные связи, фазовые портреты, инвариантные области, кольцевые генные сети, существование циклов, стационарная точка.
In this paper, we consider the dynamical systems that model the functioning of circular gene networks, in which the rate of the concentration change of any substance monotonically depends on the concentration of the substance preceding this one. Previously, the problems of existence, uniqueness and stability of periodic solutions for dynamical systems, whose equations contain monotonically decreasing and monotonically increasing functions, were studied, and a positively invariant domain was constructed to describe the phase portraits of such systems. In the case of large dimensions, it becomes necessary to narrow the constructed domain. For this, the authors developed an algorithm based on the snake-in-the-box problem from graph theory and applied to dynamical systems of various types and different dimensions. The algorithm allows us to construct a transition diagram corresponding to the shifts along the trajectories of dynamical systems. The examples of using the algorithm in the case of a six-dimensional system with step functions and in the case of a ten-dimensional system with smooth functions have been also given. The results of this construction have been applied in the study of the existence of cycles for these systems.
Key words: dynamical systems, feedbacks, phase portraits, invariant domains, circular gene networks, existence of cycles, equilibrium point.
DOI 10.14258/izvasu(2019)4-12
Введение. Периодические решения динамических систем в математике соответствуют периодическим режимам функционирования кольцевых генных сетей в биологии. Наличие таких осциллирующих режи-
*Работа поддержана грантом 18-01-00057 РФФИ.
мов является важной характеристикой происходящих в генных сетях биохимических процессов. Ранее существование и устойчивость циклов изучались в [1-4]. Подробнее биологическая интерпретация описана в [5].
1. Математическая постановка. Рассматривается система, каждое уравнение которой имеет вид:
^jj = fj(x-i)" k
(1)
где 3 = 1,..., п; 3 - 1 := п при 3 = 1, а могут являться как ступенчатыми функциями следующего вида
Lj (y) =
Г (y) =
aj kj, 0 < y < 1; > 1,
4 0 <y< 1;
. bj lj, У > 1,
так и гладкими монотонно возрастающими или монотонно убывающими функциями. Эти функции характеризуют положительные или, соответственно, отрицательные обратные связи в моделях кольцевых генных сетей.
В случае гладких функций количество монотонно убывающих функций нечетно. Построим для этого случая в положительном октанте параллелепипед Qn следующим образом:
Qn := П[0, Aj].
j=i
max fj (xj-i)
Здесь Aj = ---. Следующие две леммы
kj
для различных систем вида (1) были установлены
в [6].
Лемма 1. Qn является инвариантной областью системы (1).
Лемма 2. Система (1) имеет единственную стационарную точку S°.
Для описания поведения траекторий системы (1) разобьем инвариантный параллелепипед Qn гиперплоскостями, проходящими через точку So = (xl,x°,... ,хП), на 2n более мелких параллелепипедов, которые мы будем называть блоками и нумеровать бинарными мультииндексами следующим образом:
{X € Qn \ xi x°; x2 x°; . . . ; xn xn
"} ,
(2)
где X = (ж1,..., хп), £1,..., £п € {0,1}, и отношения порядка задаются следующим образом: символ ^о означает <, а символ означает >.
В случае системы со ступенчатыми функциями инвариантной областью является описанный выше параллелепипед Qn, где Ау = ау или Ау = Ьу, который разбивается на блоки гиперплоскостями, проведенными через точку Е = (1,...,1).
В этой точке, а также при любом пересечении гиперплоскостей ху = 1 функции в правых частях системы имеют разрыв первого рода.
Описанное выше разбиение инвариантной области задает ориентированный граф, вершинами которого являются вершины п-мерного единичного куба, а ребра обозначают переходы из блока в блок. Нумерация каждого блока, т.е. вершина ориентированного графа, имеет вид Е = {£1£2...£п}, а из каждой вершины искомого цикла должно выходить только одно ребро.
«Змеей в ящике» называют (п, ¿)-цепь при d = 2, где п — это длина слов в цепном коде, а d — расстояние Хемминга, определенное в [7]. В графе единичного куба цепным кодом является путь по ребрам куба. Задача заключается в нахождении цикла в графе без самопересечений и хорд и длины этого цикла. Для максимальной длины таких цепей в [7] были получены верхние оценки.
Мы говорим, что валентность блока Е = {£1,..., £п} равна I, если число его соседних блоков Еу такое, что Е ^ ЕЕ равно I.
Лемма 3 (Казанцев М.В. [8]). Возможные переходы из блока в блок для систем, описываемых схемой (1), можно построить, пользуясь правилом: для убывающих функций любая пара чисел {00} переходит в {01}, а пара {11} переходит в {10}; для возрастающих функций наоборот — {01} ^ {00}, {10} ^{11}.
Замечание. При переходе из блока в блок их валентность либо сохраняется, либо понижается на 2.
Поэтому для решения поставленной задачи достаточно найти блок валентности 1, т.е. вершину в графе, из которой выходит только одно ребро.
2. Алгоритм построения блоков валентности один. 1 шаг. Пусть первый индекс начального блока равен нулю, т. е. набор его индексов выглядит следующим образом: Е1 = {0£2 .. .£п}. Чтобы соблюсти условие одновалентности, построим начальный блок так, чтобы изменилась, например, только первая пара индексов {0 £2 }, а остальные индексы остались неизменными.
2 шаг. Далее возможны две ситуации:
а) если второе уравнение системы (1) содержит убывающую функцию, то для того чтобы первая пара индексов поменялась, нужно, чтобы £2 = 0;
б) если второе уравнение системы (1) содержит возрастающую функцию, то для того чтобы первая пара индексов поменялась, нужно, чтобы £2 = 1.
3 шаг. Чтобы блок оставался одновалентным, т. е. чтобы из него больше не было никаких переходов, нужно заполнить оставшиеся индексы таким образом, чтобы ни один из них при переходе в следующий блок не изменился.
x
Известия АлтГУ. Математика и механика. 2019. №4(108)
Если на втором шаге срабатывает условие а) и мы получаем, что индекс е2 равен 0, тогда:
а) если третье уравнение системы (1) содержит убывающую функцию, то для того чтобы первая пара индексов не менялась, нужно, чтобы е3 = 1;
б) если третье уравнение системы (1) содержит возрастающую функцию, то для того чтобы первая пара индексов не менялась, нужно, чтобы £з = 0.
В случае если при выполнении второго шага выполнился пункт б) и е2 = 1, то:
а) если третье уравнение системы (1) содержит убывающую функцию, то для того чтобы первая пара индексов не менялась, нужно, чтобы е3 = 0;
б) если третье уравнение системы (1) содержит возрастающую функцию, то для того чтобы первая пара индексов не менялась, нужно, чтобы £з = 1.
Дальнейшие шаги выполняем аналогично шагу 3 до тех пор, пока не будет найден последний блок, из которого возможен переход только в блок £1, и система не замкнется.
Далее приведем примеры построения блоков валентности 1 для конкретных динамических систем.
3. Пример построения диаграммы для системы размерности 6. Опишем построение блока валентности 1 и диаграммы сдвигов вдоль траектории системы, правая часть которой имеет вид LГLГLГ. Симметричный случай подобной динамической системы был предложен в [9].
Точное максимальное значение длины «змеи» и цикла в булевом шестимерном кубе было получено и равно 26, согласно [7], тогда в инвариантном параллелепипеде содержится не более 26 блоков валентности 1.
Пусть для определенности в блоке, с которого начинается диаграмма, £1 = 0, т.е. х1 < 1. Построение можно начать и с £1 = 1, поскольку траектория рано или поздно пересечет гиперплоскость х1 = 1.
Далее требуем, чтобы переход осуществлялся через грань у1 = 1, т.е. значение е2 в следующем блоке должно поменяться. Определим знак производной во втором уравнении системы, которое имеет вид:
^ = Г1(х1) - Ьт.
В этом уравнении Г1(х1) = 0, поскольку для любой точки внутри искомого блока х1 < 1. Тогда при у1 > 1 траектория выходит из блока только через грань у1 = 1 и е2 = 1. Заметим, что при переходе значения остальных е^ для следующего блока не меняются.
Для переменной х2 имеем L2(y1) = 0, тогда производная отрицательна, при х2 < 0 траектории системы не выходят через грань х2 = 1. На
этом шаге нумерация искомого блока имеет вид {010е4е5еб }.
Продолжая аналогичные рассуждения для остальных индексов, получаем блок {010011}. Проверим последний индекс: при е6 = 1 Д(у3) = L1(y3), у3 > 1, а значит, через грань х1 = 1 траектория выходить не может.
Таким образом, можно построить диаграмму сдвигов вдоль траектории системы, начиная с найденного блока валентности 1:
... ^
{010011} ^ {000011} ^ {001011} ^ {001111} ^ {001101} ^ {001100} ^ {101100} ^ {111100} ^
{110100} ^ {110000} ^ {110010} ^ {110011}
^____
Следует отметить, что если зафиксировать значение индекса для первой переменной (е1 = 0) и определить направление перехода, выбрав грань (у1 = 1), можно построить блок валентности 1 единственным образом. Если же выбрать на первом шаге другую переменную и другую грань, по которым осуществляется сдвиг, и зафиксировать другое значение выбранного е^ (не 0, а 1), то получим блок, принадлежащий этой диаграмме. Следовательно, выше указаны все блоки валентности 1, поэтому для поиска цикла инвариантную область можно сузить до объединения этих 12 блоков.
4. Пример построения диаграммы для системы размерности 10. В работе [6] рассматривалась десятимерная модель кольцевой генной сети, у которой правые части уравнений задаются функциями: = L2j-1,f2j =
^ ,3 = 1710.
С помощью описанного выше алгоритма строится диаграмма
{1100110011} ^ {0100110011} ^ {0000110011} ^ {0010110011} ^ {0011110011} ^ {0011010011} ^ {0011000011} ^ {0011001011} ^ {0011001111} ^ {0011001101} ^ {0011001100} ^ {1011001100} ^ {1111001100} ^ {1101001100} ^ {1100001100} ^ {1100101100} ^ {1100111100} ^ {1100110100} ^ {1100110000} ^ {1100110010} ^ ...
переходов траекторий этой системы по блокам валентности 1. Пользуясь леммой 3, заметим, что переход из блока £1 = {1100110011} можно осуществить только один, по первой паре чисел: пара мультииндексов {11} соответствует убывающей функции, а значит, переходит в пару {01}. Таким образом, делаем вывод, что блок Е1 переходит в блок Е2 = {0100110011}. Откуда, согласно лемме 3, можно совершить только один переход, по второй паре чисел: {10} переходит в {00}, соответственно, третий блок Е3 будет состоять из следующего набора мультииндексов: {0000110011}.
И так далее, продолжая аналогичные рассуждения, мы замыкаем цепочку блоков последним блоком {1100110010}, откуда траектории переходят в блок Е\.
В работе [6] подобные построения были проведены и для восемнадцатимерной системы вида (1), для которой с помощью дискретизации (2) фазового портрета была установлена неединственность цикла.
Аналогичная диаграмма для несимметричной динамической системы приводится в [10]. Во всех случаях построение диаграмм из блоков валент-
ности 1 проводилось с помощью описанного выше алгоритма. Подобные динамические системы рассмотрены ив [11,12].
Данные рассуждения справедливы для любой размерности систем вида (1), в которых количество убывающих функций нечетно.
Заключение. Применение разработанного алгоритма позволило получить результаты, описанные в ряде других работ (см. [6,10]). Авторы искренне благодарны В.П. Голубятникову за постановку задач, приведенных в работе, а также за полезные обсуждения и замечания.
Библиографи чески йсписок
1. HastingsS., TysonJ.J.,Webster D.Existenceofperiodic solutions for negative feedbacks cellular control systems // Journal ofDifferectial Equations. 197B. V. 2C.
2. Колесов А.Ю., РсзовИ.Х., СадовничийВ.А. Периоди-чегкие вешениявипа бещщих волн в кольцзсыхг/ннылне-тяе //ИзвАстияААН.серАя мотематаческая. 2016. T.C0. № 3.
3.Глызин С.Д., Колесов А.Ю., Розов Н.Х. Периодические решения типа бегущих волн в кольцевых цепочках однонаправленно связанных уравнений // Теоретическая и математическая физика. 2013. Т. 175. № 1.
А. Gla2s L., fa/tercacf JH. SBable O/cületions inM/tlumatiral bCodeas tS Biological Csrtrol Systems // J. of Math. Biology. 1978. V. 6.
Гс Murray J.D. Mathematical Ыо1оgy:I. An introduttion, third edition — Springer. 2007.
6. Голубятников В.П., Кириллова Н.Е. О циклах в моделях функционирования кольцевых генных сетей // Сиб. жур. чистой.прпкл.матемашхт.Г:18. Т. 18. №i.
7. Евдокимов А.А. Цепныекоды и snakeMnThe-box problem // УченыезапискиКазансеого университета, фи-зико-маткммеическиевеуки. 2014.T.156, кн.3.
8. Казке/вв М.М. Онекстдбыи; свойствах rj^E^^^c^^ доме-нов динамичепкихдистем/3 СиН.жпмидлдедустрггасьной математики. 2015. Т. 18.№ 4.
9. Elowitz M.B., Leibler S. A Synthetic Oscillatory Network of Transcriptional Regulators // Nature. 2000. V. 403.
10. Bukharina T.A., Furman D.P., Golubyatnikov V.P., KirfflevaN.E. et al. Mathematical and numetica1 models oftwo as.mmetikgene netwepks // SEitlR. 2018. V.15.
1L.arksH.T., Mahaffy J.Mftability of cydic gene modnls foesystemrinvolvmg acpdens.on // the5r. BioL 197p.V.74.
12. Аюпова Н.Б., Голубятников В.П., Казанцев М.В. О сущеатвоопнии 1еиклавьдной нееимметричнпй пнпели мосскулпуногйуеепупосулптора12 Соб.журпплсышвсйе-овлпеей математики.2017. Т. 20. № 2.
