УДК 517.929
Запаздывание в регуляции популяционной
1
динамики — модель клеточного автомата
А. Ю. Переварюха
Санкт-Петербургский институт информатики и автоматизации РАН Санкт-Петербург 199178. E-mail: [email protected]
Аннотация. Обсуждается метод математической формализации запаздывающей регуляции биологических процессов, в том числе предложенная нами модификация дифференциального уравнения для критического варианта развития популяционных флуктуаций. Использование популяционных моделей с отклоняющимся аргументом x = rv(x)f (x(t — т)) в некоторых случаях противоречит данным экспериментальных наблюдений. Рассмотрена проблема сущностной интерпретации величины т в контексте причин её появления. Для демонстрации эволюции ситуации с комплексом реалистичных факторов временного последействия предложен алгоритм клеточного автомата с троичным состоянием клеток. Условия трансформации состояния клеток показывают, что в гораздо большей степени формализуемое в уравнениях Хатчинсона и «blowfly's equation» запаздывание относится к динамике взаимодействия вида и поддерживающей условия жизни среды. Взаимодействие не может отражаться константным параметром. Действие запаздывания весьма поверхностно можно отождествлять с характеристиками непосредственно биологического вида. Диффузионную составляющую в нелокальных модификациях уравнений логичнее отражать зависимой от времени реакции величиной.
Ключевые слова: уравнения с запаздыванием, критические сценарии популяционной динамики, механизмы регуляции, клеточный автомат.
Delay in the regulation of population dynamics — cellular automaton model
A. Yu. Perevaryukha
St. Petersburg Institute for Informatics and Automation of the Russian Academy of Sciences, St. Petersburg 199178 .
Abstract. The article discusses the technique of mathematical formalization delayed regulation of biological processes. The modification of the differential equation proposed by us earlier for a specific variant of development of the population process is compared with analogues. The use of population models with deviating argument X = rv(x)f (x(t — т)) in some cases, contradictory to the experimental results. The problem of the intrinsic interpretation of the value т in the context of the causes of its appearance is considered. To demonstrate the evolution of the situation with a complex of realistic factors of temporary aftereffect, an algorithm for a cellular automaton with a ternary cell state is proposed. The cell state transformation algorithm shows that, to a much greater degree, the delay formalized in the Hutchinson equations and "blowfly's equation" refers to the interaction dynamics of the species and the supporting environment of the medium, which can not be expressed by a constant parameter. The action of delay in a very abstract sense can be identified with a directly characteristic
1 Работа выполнена в рамках проекта РФФИ: № 17-07—00125 (СПИИРАН).
© А. Ю. ПЕРЕВАРЮХА
of a biological species. The diffusion component in the modifications of these equations is more logically to represent by the time-dependent response of the systems.
Keywords: delay equations, critical scenarios of population dynamics, regulation mechanism, cellular automaton.
MSC 2010: 34A09, 65L07
1. Введение
В исторически первой модели экспоненциального роста численности популяции N(t) = N(0)ert для описания изменений состояния малой популяций не предполагалась отдельная формализация действующих механизмов регуляции в функциональной форме, только на уровне корректировки базового значения репродуктивного параметра r. Параметр r — репродуктивный потенциал (он же «мальтузианский параметр») объяснялся как разница между мгновенной рождаемость и смертностью: r = а — fi. Полагалось априори r > 0 независимо от t. Уравнение N = rNt экспоненциального роста давно не рассматривается всерьёз, упоминается как курьёз из истории математической биологии [3], хотя очевидно из многих примеров инвазий чужеродных видов, что свойства репродуктивной активности просто не могут сохраняться Vt. В дальнейшем развивались двухвидовые модели взаимодействия. Одной из проблем стало подтверждение на практике наблюдений или в экспериментах существования замкнутых циклических траекторий для численности противоборствующих популяций, которые предсказывает известная модель Вольтерра. Оказалось, что урожайность поколений канадских зайцев из хрестоматийного примера колеблется сама по себе и даже вне ареала рысей. Колебания двух видов получил в экспериментах С. Утида для другого типа биологического взаимодействия «паразит-хозяин», но эти флуктуации жука и осы даже приблизительно не походили на орбитально устойчивый цикл. Для объяснения противоречий возникла гипотеза, что эффект саморегуляции действует не от текущего состояния, но от существовавшего в прошлом и проявляется через некоторое время. Такая гипотеза подходила для имевшегося математического аппарата из области моделирования упругих деформаций--дифференциальных
уравнений с отклоняющимся аргументом.
2. Формализация запаздывающей регуляции
В ряде наблюдений и лабораторных экспериментов подтвердилось [16], что флуктуации численности могут возникать у изолированных популяций (не испытывающей межвидового трофического взаимодействия). Для формализации появления популяционных колебаний одновидовой саморегулируемой системы Хатчинсоном [14] предложено логистическое уравнение с запаздыванием:
где K — поддерживающая емкость среды обитания. При увеличении величины запаздывания т (или репродуктивной r—характеристики) реализуется бифуркация Андронова-Хопфа [1]. Чаще всего появление запаздывания связывали с продолжительностью и этапностью формирования взрослых половозрелых организмов. При запаздывании немного большем критического значения возникший цикл быстро приобретает релаксационную форму с очень низкими и продолжительными минимумами. Для улучшения характеристик цикла модели колебаний предложено значительное количество модификаций уравнения (1) (наиболее известно «food-limited equation» [10], [11]), в том числе с периодически возмущаемыми коэффициентами в [12]. Экстраординарное усложнение говорит о том, что методика согласования моделей и данных наблюдений зашла в тупик. Можно включить непосредственно в показатель степени N(t) = N(0)erf(t) функцию f (x) : f (x + Z) = f (x), max f (x) > 0, min f (x) < 0 получив колебания без всякого N(t — т), но возникнет неразрешимый вопрос, откуда в экологических факторах (оперирующих рациональными числами) могут возникнуть тригонометрической формы зависимости...?
В экспериментах австралийского энтомолога Никольсона с лабораторной популяцией мухи возникали сложные флуктуаций поколений и происходили смены режимов поведения при изменении количества и регулярности поступления в ящики с мухами корма, по современному анализу у колебаний выделяются две частотные составляющие [5]. Однако, количество корма не является бифуркационным параметром для уравнения (1). Доступности питательных веществ в (1) скорее соответствует масштабирующая амплитуду колебаний величина K, отождествляемая чаще с пределом «экологической ниши». В среднем плодовитость и время прохождения стадий онтогенеза особей можно считать постоянными на период наблюдений. Очевидно, в эксперименте можно было имитировать перемены скорости восстановления ресурсов. Недостаток питания влияет на смертность личинок, но данная величина не обособлена в модификациях уравнения Хатчинсона.
Для лучшего согласования экспериментов требовался непосредственный учет мгновенной убыли от текущей численности в (дополнительный бифуркационный параметр) одновременно с подавлением конкуренцией репродуктивного потенциала. Было предложено принципиально другое уравнение с запаздыванием, известное как « Nicholson's blowflies equation» [4] с отвечающей за саморегуляцию эффективности воспроизводства экспоненциальной нелинейностью (Y,y, в = const > 0):
dN = YN (t — т )e-YN (t-T > — /3N (t), (2)
где т полагалось временем, необходимым для выхода их экзувия взрослой мухи, а величина 1/y характеризует состояние популяции с максимально эффективным воспроизводством [13]. С поведением модели мясных мух связано несколько открытых проблем [4]. При больших т решение показывает сложные колебания (график 2 на рис. 1, т = 15, в = 0.3). При меньшем запаздывании и малых в ^ 1 демонстрирует затухающих осцилляции (кривая 2 на рис. 1: т = 8, [3 = 0.015) —
сценарии мягкого завершения инвазионного процесса чужеродного вида с единственным пиком, не относящегося к экстремальному развитию вспышки. К такому
Рис. 1. Релаксационный цикл и затухающие осцилляции в (2).
сценарию относится, например, инвазия гребневика Миетгорзгз 1е1с1у1 в Черное море (рис. 2 из [15]), но динамика его численности после вселения в Каспийское море в известные модели не укладывается (подробнее в следующей статье). Значимого
Рис. 2. Инвазия гребневика Мпетгорзгз 1елёу1 в Черное море: 1 — динамика численности гребневика; 2 — объемы вылова анчоуса.
запаздывания в резком падении уловов анчоуса и их восстановлении при нараста-ISSN 0203-3755 Динамические системы, 2017, том 7(35), №2
ющей инвазии гребневика на графике не прослеживается. Характеристика y ^ 1 по-видимому была заимствована из известной дискретной модели расчета пополнения запасов популяций рыб Рикера после сезона нереста:
Nn+i = TNne-YNn - QNn, (3)
где 0 < Q < 1 — доля промыслового изъятия. Тут параметр Y > 1 несколько отличается по трактовке от r из (1), и выражает максимально возможную эффективность пополнения при e-lNn — 1. Для (3) реализуется сценарий хаотизации M. Фейгенбаума [7] через каскад удвоений периода цикла p = 2г при возрастании Y (или сокращении Q). В итерационной модели (3) тоже можно учитывать запаздывание, указывая в вычислениях Nn-p:
Nn+i = YNne-YNn-1 - QNn,
но такое включение Nn-1 приведет к качественно другим метаморфозам фазового портрета: появлению альтернативных циклов, которые сложно интерпретировать биологически. В современной теории формирования пополнения промысловых запасов рыб применяются непрерывно-дискретные модели выживаемости на интервале t £ [0, T] [17]. Существует мнение о том, что дискретные модели лучше непрерывных для описания флуктуаций. В итерациях xn+1 = Ф(хп; а) функций, удовлетворяющих критериям теоремы Сингера, естественным образом возникают периодические точки ф'а(хг) = Фп+р(хг) после потери устойчивости стационарной точки ф(х*) = х* при ф'(х*) = -1,а = а, и в дальнейших бифуркациях удвоения периода при а > a,p = 2г ,i = 1... то. Однако, помимо периода, циклические траектории {Фп(х0)} различаются порядком обхода составляющих их точек.
Необходимо упомянуть, методика моделирования расширялась с рассмотрением динамики пространственного распределения особей — включения диффузии. Обычно предполагалась так называемая «Fickian diffusion», где поток особей пропорционален их концентрации, что приводило к модификации уравнений с добавлением диффузионного слагаемого, для (1) например так:
dN (ХЛ)= rN (x,t) (1 - N (xj- T)) + вИ (4)
81 \ К ) дх2
где О — диффузионный коэффициент [9]. Однако у подобного подхода имеется принципиальное возражение — запаздывание и пространственная диффузия не являются независимыми, то есть особи не находились изначально вместе в некоторый момент в прошлом.
В предыдущей работе [2] мы модифицировали (1) следующим образом:
f = nN (: - - N(t - т)), (5)
где дополнили концепцию предельной поддерживающей емкости среды К действующим на выраженность регуляции предпороговым уровнем численности Н.
В новой модели релаксационный цикл после бифуркации Андронова-Хопфа оказывается переходным режимом существования. Последующее образование псевдопериодической неограниченной траектории трактуется как катастрофическая популяционная динамика, сценарий «проблемы Rapa Nui» загадочной гибели цивилизации о. Пасхи. Сценарий с разрушением среды характерен для островных популяций, как в случае исчезновения быстро размножавшихся северных оленей выпущенных на острове Беринга в конце XIX в., при очень низкой скорости восстановления кормовой базы — лишайника (в экологии говорят о P/B—коэффициенте, отношении продуктивности к биомассе). Подобные модификации модели представляют значимость для исследования редких сценариев популяционной динамики, которые мы называем экстремальными, как например масштабные вспышки численности вредителей, заканчивающиеся дефолиацией лесных массивов в субарктических регионах.
Таким образом, существует проблема экологического истолкования запаздывания т, определения связи его величины с какой-то непосредственной популяци-онной характеристикой или свойствами условий обитания. Можно назвать много случаев, когда временные масштабы обсуждаемых критических явлений даже приблизительно не совпадают в интервалом t £ [0, T] индивидуального развития. При повторяющихся с более чем четвертьвековыми промежутками пилообразных вспышках численности никакие сдвиги характеристик онтогенеза короткоцикло-вой бабочки еловой листовертки не могут быть причиной бифуркации, приводящей к появлению серии отстоящих друг от друга крайне многочисленных поколений [6]. Долговременные флуктуации характерны для тихоокеанской сельди, хотя этот вид способен достигать зрелости на второй год. В упомянутых экспериментах проявляется тот факт, что скорость восполнения ресурсов явным образом влияет на колебания, и вероятно должна отражаться в величине запаздывающей регуляции потребителей. Для понимания проблем ситуации, которую призваны моделировать уравнения, переведем непрерывную формулировку в задачу дискретной динамики, но не в итерационную, а в алгоритмическую.
3. Клеточный автомат запаздывающей регуляции
Непрерывность единиц измерения совсем не обязательна для популяционной динамики, потому возникают альтернативные методы описания процессов, автоматы и когнитивные графы. Клеточный автомат, знаменитый как «Conway's Game of Life» [8], был популяризован Джоном Конвеем для объяснения процессов самоорганизации в различных естественных науках и просто как наглядная модель при обучении программированию. Популярность объясняется всего двумя правилами расчета следующего поколения группы клеток, которые могут быть живыми либо мёртвыми. Если рядом с мёртвой клеткой три живые, то мёртвая становится живой. Живая остаётся в своём состоянии, если рядом есть две или три живые.
Игровое поле может быть ограниченным, замкнутым--в виде компьютерной
эмуляции поверхности тора (наиболее часто встречающийся вариант) или бесконечным, как изначально полагал Конвей. При этом простота правил с бинарным
разделением клеток обеспечивает огромное разнообразие форм и интересных вариантов расстановки живых клеток в первом поколении, которые приводят устойчивым или периодическим вариантам итоговой расстановки, называемых фигурами. Существует модификация игры с непрерывным пространством без клеток с очень впечатляющей визуальной реализацией.
Было бы интересно создать принципиально отличный вариант клеточного автомата, где можно было бы рассматривать включение и выключение эффектов запаздывания. Для иллюстрации проблемы параметрического выражения действия запаздывания в динамических моделях и как новый объект для поиска форм самоорганизации мы предложим математическую игру с более сложными правилами, где клетка автомата будет иметь не два, а три допустимых состояния. Возможно, параметр автомата, относящийся к сущности возникновения запаздывающей регуляции, будет интересно изменять по ходу игры, имитируя эволюционную адаптацию сообщества к появлению нового вида.
Определим следующие правила алгоритма новой игры:
I. Задано аналогичное игре «Жизнь» поле клеток (в замкнутом варианте). Каждая клетка имеет 8 смежных. Изначально в каждой клетке произрастает дуб. В клетке содержащей дуб может поселиться грызун.
II. Взрослый грызун каждый сезон может плодить давать г(= 2) потомка, если хотя бы в одной из смежных с его клеткой тоже живет грызун.
III. Если вокруг грызуна оказывается занято более 5 клеток, то грызун погибает от перенаселенности, а дуб становится свободным.
IV. Грызуны подтачивают корни дерева. Дуб падает и грызун погибает. До падания дуба грызун может размножиться не более к(= 3) раз.
V. Потомки грызуна станут взрослым через два сезона, если к этому моменту займут пригодный дуб, иначе он гибнет. Новый грызун занимает ближайший свободный дуб с наименьшим числом занятых соседних клеток, осматривая смежные по часовой стрелке.
VI*. Для занятия дуба грызун может мигрировать за сезон на расстояние не более 1(=2) соседних клеток (аналог скорости света в игре Конвея). [Опционное миграционное усложнение.]
VII. На месте упавшего дуба вырастает новый за 1(= 6) сезонов.
В правилах преобразования состояния клеток заложено сразу несколько аспектов, относящихся к явлению запаздывающей регуляции. Более того, эти факторы являются противоборствующими в игре. Параметр 1 мы предлагаем сделать «управляющим» для оценки воздействия в сценариях с замедлением темпов восстановления ресурсов, что не было в алгоритме «Жизни», но согласуется с экологической реальностью.
Заключение
В анализе алгоритма Конвея оказался интересен поиск нетривиальных начальных расстановок, но в динамике экологических процессов, стремящих к некоторому установившемуся стационарному или колебательному режиму, начальные условия не столь значимая характеристика. Новый алгоритм преобразования клеток включает три динамически взаимодействующих фактора: онтогенетическую задержку, необходимость восстановления ресурсов для дальнейшего развития популяции и диффузионную составляющую, как мы видим зависящую от темпа развития взрослого организма, т.е. оказывается D(T). Будем считать, что четвертый возможный фактор накопления отравляющих продуктов метаболизма, связанный с зависимостью от предшествующих состояний, преодолевается за счет диффузии. Предложенный автомат базируется на гипотезе, что запаздывающая регуляция не относится к свойствам вида или только среды, но представляет собой аддитивную характеристику процесса взаимодействия или прямого противоборства вида в конкретном биотическом окружении. Для развития методов моделирования перспективным выглядит пересмотр формализации предела экологической ниши, от чего отказались в «blowflies equation», но что можно вернуть, используя альтернативную Yxe-YX унимодальную функцию. В условиях нарастающей инвазии эта величина не может выражаться единственным параметром, не зависящим от времени, так как видом-вселенцем часто генерируются последовательно уменьшающиеся пики численности. В современных работах проблему описания экстремальной динамики инвазий не удается разрешить за счет метода множественности включенных величин с запаздывания [18] «multiple state-dependent delays». Альтернативный подход к моделированию заключается в описании непосредственного взаимодействия в системе ресурс потребитель, где можно отразить естественную экологическую характеристику среды P/B—коэффициент, но этот сценарий не годится для перехода к колебаниям у лабораторной популяции.
Список цитируемых источников
1. Борздыко В. И. Об исследовании популяционной модели Хатчинсона // Дифференциальные уравнения. — 1985. — T. 21. — C 316—318.
BorzdykoV. I. An investigation of Hutchinson's population model (in Russian). Differ. Uravn. 21, No.2, 316-318 (1985).
2. Переварюха А. Ю. Модель сценария популяционного кризиса в результате бифуркации Андронова-Хопфа // Динамические системы. — 2016. — Т. 6(34), №2. — C. 149159.
Perevaryukha A. Y. The scenario of the population crisis as a result of the Andronov-Hopf bifurcation (in Russian). Dinamicheskie Sistemy 6(34), No.2, 121-130 (2016).
3. BacaerN. A Short History of Mathematical Population Dynamics. — London: SpringerVerlag, 2011. — 160 p.
4. Berezansky L, Braverman E, Idels L. Nicholson's blowflies differential equations revisited: Main results and open problems // Applied Mathematical Modelling. — 2010. — Vol. 34. — P. 1405—1417.
5. BrillingerD. The Nicholson blowfly experiments: some history and EDA // Journal of Time Series Analysis. — 2012. — Vol. 33, Iss. 5. — P. 718-723.
6. Cooke B, NealiSV., RegniereJ. Insect Defoliators as Periodic Disturbances in Northern Forest Ecosystems // Plant disturbance ecology: the process and the response. — Burlington.: Elsevier, 2007. — P. 487-525.
7. FeigenbaumM. Universal behavior in nonlinear systems // Physica D. — 1983. — Vol.7. — P. 16-39.
8. GardnerM. Mathematical Games —- The fantastic combinations of John Conway's new solitaire game «life» // Scientific American. — 1970. — Vol. 223. — P. 120-123.
9. Gourley S. A., RuanS. Dynamics of the diffusive Nicholson's blowflies equation with distributed delay // Proceedings of the Royal Society of Edinburgh Section A Mathematics. — 2000. — Vol. 130. — P. 1275-1291.
10. Gopalsamy K, KulenovicM., Ladas G. Time lags in a «food-limited» population model // Applicable Analysis. — 1988. — Vol. 31. — P. 225-237.
11. Gopalsamy K. Global stability in the Delay-logistic Equation with discrete delays // Houston J. Math. — 1990. — Vol. 16. — P. 347-356.
12. GopalsamyK, KulenovicM., Ladas G. Enrivonmental periodicity and time delays in a «food-limited» population model // Journal of Mathematical Analysis and Applications. — 1990. — Vol. 147. — P. 545-555.
13. Gurney W., BlytheS. P., NisbetR. M. Nicholson's blowflies revisited // Nature. — 1980. — Vol. 287. — P. 17-21.
14. Hutchinson G. An Introduction to Population Ecology. — New Haven.: Yale University Press. — 1978, 260 p.
15. Kideys A. E. The invasive ctenophore Mnemiopsis problem in the Black and Caspian Seas // Biomare Newsletter. —. 2002. — Vol. 3. — P. 5-6.
16. Nicholson A. An outline of the dynamics of animal populations // Australian Journal of Zoology. — 1954. — Vol. 2, Iss. 1. — P. 9-65.
17. PerevaryukhaA. Y. Cyclic and unstable chaotic dynamics in models of two populations of sturgeon fish // Numerical Analysis and Applications. — 2012. — Vol.5. — P. 254-264.
18. RuanS. Delay Differential Equations in Single Species Dynamics // Delay Differential Equations and Applications. — Berlin.: Springer, 2006. — P. 477-517.
Получена 09.03.2017