Нелинейная динамическая модель системы запас-пополнение Текст научной статьи по специальности «Математика»

^^МОДЕЛИРОВАНИЕ СИСТЕМ И ПРОЦЕССОВ

УДК 629.075

НЕЛИНЕЙНАЯ ДИНАМИЧЕСКАЯ МОДЕЛЬ СИСТЕМЫ ЗАПАС-ПОПОЛНЕНИЕ

А. Ю. Переварюха,

аспирант

Санкт-Петербургский институт информатики и автоматизации РАН

Исследуется динамическая система, основанная на модели запас-пополнение Рикера. Делается вывод об ограниченности применения известных дискретных моделей пополнения популяций. Предлагается модифицированная непрерывно-дискретная модель воспроизводства, качественно отличающаяся от модели Рикера.

Введение

Контринтуитивные процессы, протекающие в сложных саморегулируемых биологических системах, не имеют простых объяснений, а «серебряная пуля» в арсенале математической теории динамики популяции либо просто не существует, либо является хорошо охраняемой тайной. В сложных ситуациях автору приходится использовать нетрадиционные для математической биологии подходы, например непрерывно-дискретные динамические системы. Необходимость в разработке описываемой в этой статье модели обусловлена целью анализа имеющихся данных наблюдений, а не является попыткой дать умозрительную популяционную интерпретацию некоему математическому аппарату.

Проблема стабильности биологических сообществ и отдельных популяций — одна из наиболее актуальных и обсуждаемых в теоретической экологии [1]. Природные системы могут иметь более одного устойчивого режима существования; если значения некоторых показателей находятся в пределах некоторого диапазона, незначительные возмущения поглощаются. Количественные показатели будут колебаться, но качественное поведение динамических процессов останется постоянным, пока внешнее воздействие не выведет систему за границы области устойчивого равновесия и переведет в иное качественное состояние. Если возможны два таких состояния и одно из них нежелательно, то задача управления состоит в поддержании параметров природной системы в диапазоне значений, удаленных от критических. Выведенная из равновесия экосистема может стабилизироваться в новом и непригодном для хозяйственного использования состоянии.

Для анализа устойчивости были предложены различные методы, в частности, использование

в качестве меры устойчивости аналога информационной энтропии. МакАртур (MacArthur) в 1955 г. предложил характеризовать стабильность сообщества величиной

Б' = --^р(в1 )1п р(в1),

1=1

где р(э) — вероятность переноса энергии по определенному пути. Чем больше Б', тем устойчивее считается сообщество. Но при использовании такого подхода возникают противоречия. Достижение максимума Б' возможно, когда равновероятен выбор любого пути переноса энергии, т. е. такое сообщество не будет обладать какой-либо иерархичностью и структурой, что не подтверждается наблюдениями над реальными сообществами. Энтропийная мера применима на ранних стадиях сукцессии, когда конкурентные отношения еще не столь жесткие и сообщество может быть рассмотрено как система со слабыми взаимодействиями [2].

Другой подход, развивающийся в последнее десятилетие, заключается в применении идей нелинейной динамики в экологических процессах и изучении качественного поведения динамических моделей [3]. В настоящей статье речь пойдет о нелинейных динамических системах, имеющих интерпретацию в теории запас-пополнение популяций рыб.

Продолжительное и экономически эффективное управление популяцией возможно только тогда, когда воспроизводство обеспечивает восстановление промыслового запаса. Амплитуда и частота колебаний численности обуславливаются изменчивостью процессов формирования пополнения и темпом обновления нерестового стада за счет новых возрастных классов. Концепции о линейной зависимости между запасом и пополнением Ф. Баранова и К. Бэра [4] не нашли подтверждения в последующих исследованиях.

Определение зависимости между запасом и пополнением показывает, при каких значениях запаса воспроизводство становится восполняющим запас, и может служить важнейшим критерием при прогнозировании уловов. Результат действия различных факторов смертности молоди выражается в итоговой численности поколения, полученного от нерестового запаса данной численности через определенный «интервал уязвимости» [0, Г].

Уравнение для численности пополнения (recruits) R предложено известным канадским исследователем У. Е. Рикером (W. E. Ricker) в 1953 г. [5] (до Рикера вопрос о возможности уменьшения пополнения при возрастании нерестового стада никем не рассматривался):

R = aS exp(-bS), (1)

где S (stock) — величина нерестового запаса; b — коэффициент, отражающий величину, обратную количеству выметанной икры, при котором число выжившей молоди максимально, соответственно, имеет смысл только b << 1; a — параметр. График зависимости представляет собой куполообразную кривую с единственным нетривиальным пересечением с биссектрисой координатного угла R = S. Количество икры определяется, исходя из средней плодовитости особей: E = XS. График зависимости числа рекрутов от численности производителей называется кривой пополнения (рис. 1).

Очевидно, что применение модели (1) вызывает ряд сложностей, наиболее значимая из которых состоит в том, что при увеличении количества отложенной икры выживаемость молоди будет стремиться к нулю, что не согласуется с наблюдениями биологов. Для ограничения возможности возникновения такого явления модель Рикера при имитационном моделировании популяций корректно использовать только до некоторого критического количества икры Ек, как делал А. Б. Казанский [6], при превышении Ек количество молоди зависит от соотношения постоянных коэффициентов модели или определяется константой. Фактически в таком подходе при имитационном моделировании формирования пополнения использу-

2 4 6 8 10 12

Число производителей • 102 ■ Рис. 1. Кривая запас-пополнение Рикера

ется только восходящая ветвь кривой запас-пополнение:

R=

aEexp(-bE), если E < 1/b

a,если E > 1/b

Ы

ГаЕехр(-ЪЕ), если Е < Ек или Б = \ .

[Ек,если Е > Ек

Не столь очевидные, но существенные проблемы для моделирования воспроизводства возникают при значительном уменьшении нерестового запаса. При крайне низких численностях производителей уравнение Рикера предсказывает увеличение эффективности воспроизводства популяции. Когда количество отложенной икры стремится к нулю, выживаемость пополнения N стремится к предельному значению:

dN

---= аехр(-ЪЕ)(1 - ЪЕ), Иш аехр(-ЪЕ)(1 - ЪЕ) = а.

dE е^о

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

Запас-пополнение Рикера как нелинейная динамическая система

Описания свойств моделей запас-пополнение как математических функций оказывается явно недостаточно. Один из пионеров отечественного гидробиологического моделирования В. В. Суханов писал следующее: «Оказалось, что собственно кривой Рикера результаты расчетов не дали. На плоскости Б х Б получилось сгущение точек, мало чем напоминающее вышеуказанную кривую. При тех параметрах, при которых популяция испытывала наиболее регулярные колебания, точки на графике образовывали нечто вроде неотчетливого эллипса. Применение метода наименьших квадратов для определения значений параметров а и Ъ привело к абсурдным результатам: коэффициент а оказался меньшим, чем единица, что невозмож-

но для популяции с ненулевой численностью. Полученные результаты говорят о том, что кривая Рикера не является хорошим приближением к эмпирическим данным» [7]. Для понимания того, что аспекты, отмечаемые многими критиками теории зависимостей запас-пополнение, на самом деле не содержат противоречий, необходимо провести дальнейшие исследования.

Представим процесс изменения состояния популяции, определяемый зависимостью запас-пополнение, в виде динамической системы — математического объекта, для которого можно указать набор динамических переменных, характеризующих состояние системы. Значения таких переменных в последующий момент времени рассчитываются из текущих значений по определенному правилу (называемому оператором эволюции). Динамическая система — это тройка (M, T, у), состоящая из фазового пространства M, времени T, оператора эволюции у. Причем для всех x е M и t, s е T выполняется условие

y(y(x, t), s) = y(x, s + t).

Множество {y(t)(x)}te T называют фазовой траекторией точки х. Графически эволюция динамической системы во времени представляется движением точек в фазовом пространстве. Важным понятием теории диссипативных динамических систем является аттрактор. Словарное значение слова «attractor» — то, что привлекает или притягивает. Под аттрактором мы будем понимать подмножество фазового пространства A с M, инвариантное относительно эволюции в системе: y(t)(A) = A для всех t е T, и такое, что существует окрестность U множестваA, в которой для всех у е U выполняется условие

lim у (t)( у) = A.

t^<^>

Простые аттракторы — устойчивое состояние

« « *

равновесия с неподвижной точкой x :

lim у (t)( у) = х*,

t^<^>

и устойчивый цикл, отвечающий режиму периодических автоколебаний: y(t)(x*) = y(t+p)(x*), p е T. Множество точек, приводящих к некоторому аттрактору, называется его областью или бассейном притяжения (basin of attraction). Рассмотрим динамическую систему как полугруппу итераций {y(i)}/ > о , и пусть Ro, R1, R2, ... — последовательность точек, описывающих эволюцию системы, определенных условием Rj+1 = y(Rj) при всех j > 0. Оператор эволюции в непрерывно-дискретной динамической системе запишем в следующей форме:

--(aN+1(0) + P)Ni+1(t); dNi±1

dt

dt

= Я

dNi

t-G

dt

Константы а, в, т соотносятся с константами a, Ъ формулы Рикера: a = Аехр(-вт), Ъ = ат. Каче-

ственное поведение системы зависит от параметра а, он является управляющим. Исследование динамической системы на основе (1) показало интересное и подчас весьма странное ее поведение. До определенного значения а, не превышающего бифуркационное, система стремится к точечному аттрактору Б* = 1п(а)/Ъ, даже если искусственно выводить ее из равновесия (рис. 2).

Первый метаморфоз поведения системы — отображенная на временной диаграмме бифуркация (рис. 3), происходит, когда производная в неподвижной точке перестает удовлетворять критерию „ 2 устойчивости при выполнении условия а > е :

y'(R ) - ae Ъ - Ъ—

у(R) - ae-ш - ЪRae-т;

a(1 - ln a)

In a .. In a

-^ , ln a -b-r-

ae

Ъ

In a

-1 - ln a,

1 - ln a--1, a-(

Теперь очевидно, что параметр Ь не влияет на топологические характеристики фазового портрета. Динамическая система стремится в устойчивое циклическое состояние с периодом 2 — глобальный аттрактор, состоящий из двух периодических точек — последовательности (ЯТ1, ЯТ2, ИТ1, ЯТ2, ...), областью притяжения которого является все фазовое пространство. Тот факт, что отображение Рикера имеет в этом диапазоне параметра цикл с периодом 2, говорит о том, что не стоит от построенного по эмпирическим данным графика ожи-

60 55 » 50 £ 45 со * 40 2 35 30 25 20 15

я

1 1 шш

НІШІ шш —1 1 IIIMJiLlLl ППпяяя

11“ и1 т

1

2

3

4

5

б 7 t • 102

Рис. 2. Притяжение к точечному аттрактору

80

70

і60 ®. 50

СО

2 40

0$ 30

9 10 11

■ Рис. 3. Бифуркация в системе на основе (1) ИHФOPMДЦИOHHO-УПPДBAЯЮШИE СИСТЕМЫ ^

t • 102

№ 2, 2008

ll

R, экз. 703 653 603 553 503 453 403 353 303 253 203 153 103 53 3

2 52 102 152 202 252 302 352 402 452 502 552 602 652 702 й,экз.

■ Рис. 4. Странный аттрактор в системе на основе (1)

2 52 102 152 202 252 302 352 402 452 502 552 602 652 702

дать характерной куполообразной кривой, как на рис. 1. Если далее увеличивать параметр а, будет происходить увеличение амплитуды колебаний, и по достижении следующего порогового значения а > 12,51 произойдет бифуркация удвоения периода и установится цикл с периодом 4. При дальнейшем увеличении параметра а будет происходить каскад бифуркаций удвоения периода. При а > 14,8 невозможно выделить устойчивых точек или замкнутого цикла — происходит детерминированный хаос, напоминающий стохастический процесс, очень чувствительный к начальным условиям. При i ^ го фазовая диаграмма будет выглядеть так (рис. 4). Траектория притягивается к подмножеству фазового пространства — «странному аттрактору», не имеющему внутри себя ни одной устойчивой траектории.

Как показали дальнейшие исследования, описанным М. Фейгенбаумом сценарием перехода к хаосу (1) не ограничивается, для некоторых диапазонов а > 14,8 наблюдаются «окна периодичности» — появляются устойчивые циклы нечетных периодов, т. е. для (1) возможны касательные бифуркации. В диапазоне значений параметра 18,474 < а < 18,564 появление ожидаемого странного аттрактора не наблюдается, возникает устойчивый цикл.

Дискретные отображения обладают свойствами, делающими их объектом исследования для теории динамического хаоса. Ли и Йорк в 1975 г. опубликовали статью «Period three implies chaos», в которой показали, что если одномерное отображение вида Rj+1 = y(Rj) при некотором значении одного из параметров имеет цикл периода три, то оно также имеет и бесконечное множество циклов других периодов. Даже самые простые дискретные модели, например: xn+1 = Asin(nx„), xn+1 = Axn(1 - xn) и т. п., могут приводить к хаотическим режимам динамического поведения. Появление хаоса в экологических моделях — отдельная проблема («Chaos

Paradox», «Paradox of Enrichment»), различные исследования по анализу эмпирических данных, например работа [8], не смогли подтвердить наличие хаоса в реальных популяциях. В работе [9] МакАлистер показал, как «хаотическое» поведение сменяется стабилизацией при наличии для популяции ограничений со стороны внешней среды.

Для построения кривых запас-пополнение были предложены довольно сложные преобразования исходных данных наблюдений. Наверное, никогда исследователи не бывают так настойчивы, как в случае, когда пытаются отыскать то, чего нельзя обнаружить. Рикер предложил прологарифмировать формулу (1) и строить кривую с использованием регрессии ln(R/S) на S для арифметической средней. Возникает естественный вопрос: имеют ли смысл методы построения кривой, если эмпирические данные о воспроизводстве не находящихся в стадии деградации популяций, для которых справедлива зависимость Рикера, будут представлены в виде сгущений точек на графике, мало напоминающих ожидаемую кривую. Именно эллипсы (как верно заметил В. В. Суханов) образуют точки на фазовом портрете, стремящиеся по спирали к устойчивому состоянию (фокусу) в его области притяжения.

Модификация модели запас-пополнение

Необходимость модификации модели запас-пополнение возникла при выявлении нелинейной зависимости численности «скатывающейся» молоди от численности нерестового стада в задаче моделирования процессов, приведших к деградации популяции осетровых Каспия. Численность пропущенных на нерест производителей севрюги за период наблюдений изменялась очень существенно: от 230 тыс. экз. в 1979-81 гг. до 15 тыс. в 2000 г. [10].

П. Летт на основе анализа воспроизводства трески показал, что двухпараметрические модели не

обеспечивают приемлемую интерпретацию данных наблюдений. Наличие ограниченных пищевых ресурсов в неявном виде учитывалось известными моделями, но общим недостатком существующих моделей является игнорирование декомпенсаци-онного фактора смертности и изменения пищевых потребностей по мере развития молоди, которое совпадает во времени с началом сезонного уменьшения кормовой базы. Декомпенсационные факторы, увеличивающие смертность при уменьшении плотности (наблюдаются, например, у насекомых, для которых характерно групповое поведение), снижают эффективность нереста, уменьшая количество икры, реально вступившей в репродуктивный процесс при Б^0. Введение такой зависимости становится очевидной необходимостью для популяций, подвергающихся «перелову». Целью модификации модели запас-пополнение стало создание гибкого математического аппарата для задачи согласования характера поведения имитационной модели динамики популяции со статистическими данными о популяциях, так как очевидно, что поведение модели зависит от математических особенностей выбранной функции воспроизводства. Куполообразная кривая подразумевает наличие диапазона оптимальной численности нерестового стада.

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

dN

= -(аю(Щ(О + 9(Б)Р^(О, N(0) = ХБ

dt dw gNn (^

, (2)

dt N (t) + (

-, 9(Б) =

1 - ехр(-сБ)

, к > п

где w(t) отражает изменение пищевых потребностей в зависимости от массы особей; g — ограничивающий фактор, учитывающий количество доступных кормовых организмов; убывающая функция 9(Б)^1 при Б^го. Графиком предложенной автором и исследованной в инструментальной среде разработки имитационных моделей AnyLogic5 новой модели является куполообразная кривая. Среда AnyLogic предлагает достаточный выбор численных методов с изменяющимся шагом интегрирования для работы с системами обыкновенных

дифференциальных уравнений в форме Коши и дифференциальных уравнений с отклоняющимся аргументом. При деградации нерестового стада кривая переходит в прямую линию с малым углом наклона. Кривая системы уравнений (2) имеет ниспадающую ветвь с уменьшающимся наклоном (рис. 5).

Полученная кривая адекватно воспроизводит динамику системы запас-пополнение и соответствует наблюдениям за воспроизводством популяций, для которых характерна значительная амплитуда колебаний численности производителей. В частности, для горбуши рек Аляски отмечено уменьшение относительного количества молоди как в периоды возрастания численности популяции, так и по мере ее уменьшения.

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

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

С—^

■ л

5

10

15

20 25

Б • 103, экз.

Рис. 5. Кривая предложенной автором модели (2)

ется. Хорошо документирован коллапс без последующего восстановления промысловых запасов форели и сельди Великих Озер [11] в 1950-е гг.,

особенности графика деградации которых подтверждают предположения автора о наличии дополнительного максимума на кривой пополнения.

Литература

1. Свирежев Ю. М., Логофет Д. О. Устойчивость биологических сообществ. М.: Наука, 1978. 353 с.

2. Свирежев Ю. М., Елизаров Е. Я. Математическое моделирование биологических систем // Проблемы космической биологии. 1972. Т. 20. С. 133-139.

3. Короновский А. А., Трубецков Д. И. Нелинейная динамика в действии. Саратов: Изд-во ГосУНЦ «Колледж», 2002. 324 с.

4. Шибаев С. В. Промысловая ихтиология. СПб.: Проспект науки, 2007. 400 с.

5. Ricker W. Stock and recruitment // Journal Fisheries research board of Canada. 1954. N 11. P. 559-623.

6. Казанский А. Б. Имитационная модель популяции рыбы с учетом промыслового воздействия // Применение методов имитационного моделирования в рыбохозяйственных исследованиях на внутренних водоемах. Л.: Госниорх, 1989. С. 33-47.

7. Суханов В. В. Исследование модели популяции нерки Oncorhynchus nerka в условиях изменчивой кормовой базы// Вопросы ихтиологии. 1973. Т. 13. Вып. 4. С.627-632.

8. Ellner S., Turchin P. Chaos in a noisy world: new methods and evidence from time-series analysis // Amer. Naturalist. 1995. N 145. P. 343-375.

9. McAllister, C. D., R. J. LeBrasseur Stability of enriched aquatic ecosystems // Science. 1971. N 175. P. 562565.

10. Довгопол Г. Ф., Вещев П. В. Оценка численности поколений севрюги Acipenser stellatus и основных факторов, влияющих на структуру ее популяции// Вопросы ихтиологии. 1993. Т. 33. С. 93-99.

11. Christie W. J. Changes in the fish species composition of the Great Lakes//Journal Fisheries research board of Canada. 1954. N 5. P. 827-853.

UlkHIIII1 ЩРІШІІ

Роман с авиацией

ИЗДАТЕЛЬСТВО «ПОЛИТЕХНИКА» ПРЕДСТАВЛЯЕТ: Андриевский А. В.

Роман с авиацией: Повесть; Технология авиакатастроф (Записки командира авиалайнера). СПб.: Политехника, 2008. 256 с. ^N-978-5-7325-0879-6

Автор этой повести, посвятивший авиации 50 лет, — инженер-пилот первого класса, налетавший на пассажирских самолетах Ил-12, Ил-14, Ил-18, Ту-154 в качестве пилота и командира корабля почти 14 000 часов. Последний выпускник спецшколы ВВС и первый выпускник Академии гражданской авиации представляет свою книгу: «Это повесть о тебе, мой юный романтик, пытающийся вырваться из тесного мирка своего двора на просторы воздушного океана, чтобы за штурвалом боевого самолета или реактивного пассажирского авиалайнера вступить в борьбу со стихией, открыть для себя новые земли и встретить верных, преданных друзей — настоящих мужчин, таких, как ты».

Книгу можно приобрести по адресу: 191023, Санкт-Петербург (ст. метро «Гостиный Двор»), Инженерная ул., дом 6, 3-й этаж, ОАО «Издательство «Политехника», с 10 до 18 час., кроме сб., вс. Тел./факс (812)312-44-95, тел. 571-61-44. Web: www.polytechnics.ru

14 S ИHФORMДЦИOHHO-УПRДВAЯЮШИE СИСТЕМЫ