1.2.4. Модель «хищник – жертва»
Одним из классических объектов приложения математической экологии является система хищник-жертва. Цикличность поведения этой системы в стационарной среде была показана ещё Вольтеррой (1931) и Лотки (1925). Их модель описывает эволюцию двухкомпонентной системы, обладающей временной и пространственной однородностью и представляет систему двух зацепленных уравнений типа:
Здесь V(X) и Р(Y) — популяции соответственно жертв и хищников; Популяция хищников увеличивается лишь за счёт поедания ими жертв – сами по себе они вымирают. Жертвы, наоборот, обладают неограниченным источником питания и при слабой связи рост их популяции ограничен внутренними причинами, которые представлены соответствующим параболичным членом.
Квадратичные члены в модели Вольтерры играют роль нелинейного трения - именно они приводят к появления на фазовой плоскости положений равновесия типа "фокус", т.е. энергия системы не сохраняется и, что существеннее, к неэквивалентности обоих направлений времени. В принципе так и должно быть, ведь реальная биологическая система, очевидно, не обладает подобной симметрией.
Таким образом, модель Вольтерра есть простейшая модель со слабой нелинейной связью. Её можно легко обобщить, разрушив пространственную или временную однородность. Это, однако, не вносит ничего нового в плане поведения системы как таковой; качественно новое поведение может быть реализовано лишь введением дополнительных сильно нелинейных взаимодействий в исходные уравнения. Другой возможностью для модификации модели является введение эффекта насыщения, что приводит к добавлению ещё одного уравнения первого порядка. Следует отметить, что это, строго говоря, реального смысла не имеет, поскольку любая динамическая система должна описываться чётным числом уравнений первого порядка. Поэтому данное обобщение представляет лишь математический интерес.
Ниже рассмотрена модифицированная модель Вольтерры: введём сильное нелинейное взаимодействие между хищником и жертвой, что позволит (при некоторых значениях параметров модели) заменить устойчивые фокусы на предельные циклы, которые и означают устойчивые осцилляции численности популяций около их равновесной населённости.
В любой реальной системе есть (или должно быть) некое характерное время, задающее все остальные временные параметры. Например, в случае математического маятника такое характерное время суть период колебаний, для обращения планеты вокруг своего Солнца – это год и т.д.. Значение такого параметра чрезвычайно велико: он позволяет сформулировать уравнения движения (то есть соотношения, определяющие временную эволюцию системы) в терминах рекуррентных соотношений между координатами системы в моменты, отстоящие друг от друга на наш временной параметр. Так, для приведённого примера с маятником эти соотношения имеют вид:
где х и p - координата и импульс маятника, ε – его энергия , а период принят за единицу. Данное уравнение описывает систему в некоторые фиксированные моменты времени (отстоящие на период в данном примере). Это очень удобно, если нас не интересует динамика системы между этими выбранными моментами. Такие системы носят названия систем с дискретным временем или отображений (в случае с осциллятором отображение двумерно). С точки зрения математики сведение дифференциального уравнения к отображению является всего лишь грубой аппроксимацией производной по времени. С естественнонаучной позиции это всегда можно сделать, поскольку рассматриваются реальные системы. Таким шагом по времени может быть выбрано среднее время жизни животного. Сечения, ортогональные оси времени, называются сечениями Пуанкаре и являются одними из самых распространённых методов исследования систем сложных дифференциальных уравнений. Осуществляется переход от системы дифференциальных уравнений к отображениям и обратно. Ниже все выражения будут иметь вид отображений.
В качестве примера рассмотрим модель Вольтерры с неким конкретным набором параметров. Фазовый портрет этой системы приведён на рис.7 и рис.7а.
Рис.7: Фазовый портрет модели Вольтерры. (1)
В данном описании система чрезвычайно устойчива и попадает в фокус за отнoсительно малое время. Видно, что этот фокус - общий для всех траекторий и является единственным положением равновесия в рассматриваемой системе. В действительности, такую стабильность системе обеспечивает нелинейное ограничение для хищников. Здесь существенен тот факт, что этот показатель неограничен по абсолютной величине и эффективно действует при больших населённостях популяции хищников.
Отсюда следует : при перенаселённости "эффект группы" приводит к снижению рождаемости и уменьшению популяции. Важно отметить, что здесь описано собственное ограничение для хищников, совершенно не зависящее от численности популяции жертв. В данном контексте, однако, действие этого ограничителя представляется тривиальным – такие системы математиками изучаются уже давно. Кроме того, отсутствие в системе неустойчивых равновесий делает её неполной.
Рис.7а: Фазовый портрет модели Вольтерры при других начальных условиях.
Здесь так же фокус является единственным положением равновесия в данной системе, что неприемлемо в реальных экосистемах.
Взаимодействие видов всего достаточно слабо - численность популяции хищников оказывается мало чувствительной к колебаниям численности жертв. Такая ситуация не всегда имеет место в природе, скорее наоборот, в реальных экосистемах связь хищников и жертв более чувствительна, что является одним из источников устойчивости экосистем (при слабой связи возможны сценарии, когда система "проскакивает" положение равновесия и уходит в область нестабильности на своей фазовой плоскости; далее её движение становится инфинитным – система уходит в бесконечность или, по крайней мере, состояние, далекое от положения равновесия). Взаимодействие видов факторизовано, а это означает, что если ввести некое выражение, являющееся функцией Гамильтона-Рэлея – Н- для данной системы, ввести "температуру" – T, как меру влияния термостата – внешней среды, то такое взаимодействие в неком производящем функционале системы (проще говоря, в статсумме ) сравнительно легко расцепляется (в том или ином приближении – приводит к появлению так называемых "виртуальных" взаимодействий) с перенормировкой всех коэффициентов в системе (1). Данная система распадается на две невзаимодействующие подсистемы: подсистему "хищников" и подсистему "жертв". Поэтому, представляет интерес рассмотреть взаимодействие, которое в общем случае не расцепляется ни при каких допустимых приближениях. Существенно, что взаимодействие в системе (1) не подавляет фокусы, что также свидетельствует о его квазилинейности. Это соответствует давно известному в теории нелинейных процессов правилу: взаимодействие, не подавляющее фокусы, всегда перенормируемо. Рассмотрим неперенормируемое взаимодействие, в классе ограниченных сильно нелинейных функций, то есть функций, которые не могут быть линеаризованы ни в какой существенной области изменения заданных переменных.
В классе ограниченных функций известно всего несколько элементарных: косинус, арктангенс и всевозможные их модификации. Арктангенс является довольно плавной функцией: он может быть легко линеаризован (кроме больших значений аргумента). Поэтому, целесообразно выбрать косинус, и его степени. Проблема состоит в том, что все эти функции будут периодическими, что существенно повлияет на эволюцию системы, так как внесёт в неё "паразитные" периоды, которые в ряде случаев могут подавить или исказить естественные циклы системы. Этой проблемы можно избежать, рассматривая функции типа sin(sin(x)). Такие конструкции существенно нелинейны, содержат в своём спектре большое число гармоник, примерно, равной амплитуды, что очень важно для поставленной цели: существование такого спектра свидетельствует о том, что данная функция хаотична и это позволяет исследовать новые свойства системы, связанные с наличием в ней хаоса.
Рассмотрим систему (1) с конкретным взаимодействием указанного типа:
Здесь – собственные амплитуды для каждой популяции, – константы связи, определяющие взаимоотношения видов.
Взаимодействие здесь имеет вид V = sin cos х2. Для того, чтобы убедиться в отсутствии какой-либо периодичности, нужно построить его Фурье-образ, для чего разложить его в ряд типа:
- прямое преобразование и
- обратное преобразование.
Суммирование по п и х ведётся по всей области их определения. Если в спектре Vn будет доминировать одна или несколько гармоник, мы можем говорить о квазипериодичной функции которая будет отвергнута, если же спектр будет содержать много одинаковых гармоник, то такая функция будет приемлемой, так как не будет содержать какого-либо одного ярко выраженного периода.
Для данного взаимодействия результаты разложения представлены на рис.8:
Рис.8: Разложение Фурье – образ «взаимодействия» между хищником и жертвой в системе (2).
Расстояние между линиями на рисунке равно элементарной частоте. Симметрия спектра относительно вертикальной оси говорит о вещественности исходной функции (так как взаимодействие есть вещественная функция, его Фурье-разложение обязано быть зеркально симметричным относительно оси ординат). Симметричность в данном случае свидетельствует о правильности проведённого анализа. Несмотря на достаточно беспорядочный спектр, функция в действительности не является хаотичной, а будет вести себя как квазипериодическая функция. Это отразится на фазовом портрете системы. Здесь важно то, что гармоники можно с большой степенью точности считать распределёнными в спектре беспорядочно, поэтому никакого дополнительного сильно выраженного периода, связанного с конкретным выбором функции взаимодействия, не возникнет. Выбранное взаимодействие обладает ещё одним важным свойством: оно ограничено при любых значениях переменных. Его введение делает пищевую цепь очень "хрупкой", то есть, любое, сколь угодно малое изменение численности жертв сильно сказывается на численности популяции хищников, то есть численность хищников оказывается сильно зависящей от населённости жертв.
Такая связь между популяциями почти всегда имеет место в природе. При увеличении численности жертв популяция хищников будет сначала возрастать (синус во взаимодействии мал), а затем вступают в силу естественные ограничения роста её численности - болезни, снижение рождаемости из-за перенаселённости и другие, причём все эти взаимодействия сильно нелинейны, как и всякий реальный процесс. При этом, если особи характеризуются достаточно большой плодовитостью, то при малом уменьшении популяции жертв, численность хищников быстро вернётся к равновесному состоянию. Совершенно очевидно, что рассматриваемое взаимодействие не отражает всех особенностей упомянутого процесса. В частности, оно не содержит первоначального быстрого повышения численности хищников, которое ведёт к перенаселённости. Тем не менее, оно содержит многие из упомянутых факторов ограничения роста в перенормированном виде. Поэтому, использование такой связи представляет интерес с точки зрения применения к неким более точным моделям, которые, однако, достаточно сложны, чтобы быть решёнными точно.
Переходя от экологических аспектов к математическим, необходимо сказать, что данное взаимодействие, несмотря на свою ограниченность, замечательно тем, что эффективно подавляет фокусы, обусловленные собственным ограничителем для хищников.
После того, как модель сформулирована, можно перейти непосредственно к численному моделированию системы (2). Проследим эволюцию системы при изменении параметров задачи.
Для установления общего вида фазового портрета системы, рассмотрим некий конкретный набор параметров, например: α = 0.98; β = 0.03; γ = 0.07; δ = 0.03; ε = 1.30; θ = 0.20; μ = 0.01. Результаты расчётов проиллюстрированы рис.9:
Рис.9: Фазовый портрет системы (2) для конкретного набора параметров. Чётко виден предельный цикл (жирная линия в левой части рисунка), на который выходят все фазовые траектории, несмотря на то, что некоторые из них испытывают довольно большие отклонения от него.
Рассмотренное поведение системы является универсальным, - у системы нет других положений равновесия, кроме данного.
Далее изображена эволюция системы в пространстве начальных условий, то есть при проходе от 1 до 4 начальные условия подбирались таким образом, чтобы "прочесать" весь предельный цикл (рис.10). Траектории, находящиеся внутри цикла разматываются к нему (см. 1— 3), а траектории, начавшиеся вне его, наоборот, сматываются к циклу (см. 4), что свидетельствует о его абсолютной устойчивости.
Рис.10: «Внутренность» предельного цикла – разные траектории наматываются на него - цикл абсолютно устойчив. Значения параметров те же, что и для рис.9. Для 1 начальные условия есть (1.4;1.4). Далее обе координаты увеличиваются на 0.2 на шаг.
Данная структура появилась вместо фокуса, существовавшего в отсутствии нелинейного члена (см. рис.7). Кроме того, исчезновение фокуса, по-видимому, также является универсальным явлением: при данном наборе параметров он не появляется ни при каких начальных условиях. То, что предельный цикл возник на месте фокуса является причиной столь быстрого «наматывания» траекторий на цикл, как это видно из рис.10(3,4). Если начать траекторию далеко от цикла, будут происходить некие её блуждания.
Эволюция предельного цикла при изменении параметров задачи по отношению к b. проиллюстрирована на рис.11. Здесь нелинейное взаимодействие полностью подавляет квазилинейный факторизованный член.
Видно, что поведение системы качественно не меняется (цифры в скобках – начальные условия, а цифры сверху – значения b). Параметр β здесь равен 0.06, 0.09 и 0.11 (в предыдущих рассуждениях он был равен 0.03). Однако, вид фазовых траекторий качественно не меняется и, поэтому, β можно принять равным 0.03.
Рис.11: Поведение системы при различных значениях параметра b при всех остальных неизменных.
Эволюция предельного цикла при изменении параметров задачи по отношению к b. проиллюстрирована на рис.11. Здесь нелинейное взаимодействие полностью подавляет квазилинейный факторизованный член.
Видно, что поведение системы качественно не меняется (цифры в скобках – начальные условия, а цифры сверху – значения b). Параметр β здесь равен 0.06, 0.09 и 0.11 (в предыдущих рассуждениях он был равен 0.03). Однако, вид фазовых траекторий качественно не меняется и, поэтому, β можно принять равным 0.03.
Далее, цикл оказывается инвариантным по отношению к изменению параметра θ (выше он был равнее 0.2). На рис.12 параметр θ = 0.87 и очевидна инвариантность, поскольку система уравнений (2) является, в известной степени, симметричной (при замене хищников на жертв и наоборот). Поэтому, в дальнейших исследованиях примем θ равным 0.2.
Рис.12. Фазовый портрет при q = 0.87. Видно, что предельный цикл качественно ничем не отличается от предыдущих случаев. Начальные условия: (0.8;0.8) .
Исследуем изменение параметра δ (рис 13,14 ), который был принят равным 0.03. Все остальные параметры остаются постоянными, принимая значения , указанные выше, в частности, γ= 0.07.
Рис.13. Изменение вида цикла при изменении начальных условий (в скобках) и при d=0.01.
Совершенно очевидно, что при малых δ вид предельного цикла меняться не должен. Так на самом деле и происходит. Ранее на рис.7 рассмотрен случай с δ = 0.01 и приведены различные начальные условия.
Как и следовало, уменьшение собственной нелинейности хищников не ведёт к существенной перестройке цикла. При увеличении δ предельный цикл размывается и система стохастизуется. На рис.14 приведены фазовые траектории для δ = 0.2, 0.3 и 0.4.
Видно, что всё пространство внутри цикла заполнено, как и должно быть для сильно нелинейной диссипативной модели, находящейся в режиме хаоса.
При меньших значениях δ с циклом ничего особенно интересного не происходит. Здесь, однако, выделяются значения, близкие к 0.05. При этих значениях предельный цикл (при некоторых начальных условиях) переходит в фокус и процесс является периодическим по начальным условиям.
Рис.14. Фазовый портрет системы при больших d (цифры на рис.).
Начальные условия везде (1;1).
Поскольку всё происходит при конкретном значении δ, такое явление является исключительно свойством рассматриваемой модели. Результаты расчётов на эту тему приведены на рис.15:
Р
ис.15. Вид фазовой плоскости системы при d=0.05 при разных начальных условиях, указанных на рис.; видна периодическая зависимость вида плоскости от начальных условий.
Рассмотренное обобщение модели Вольтерры с введением дополнительного нелинейного взаимодействия между хищником и жертвой имеет на своей фазовой плоскости предельный цикл, который означает колебания численности обеих популяций около некоторого равновесного значения. Данный предельный цикл является очень устойчивым и свидетельствует об устойчивости полученной экосистемы по отношению к внешним воздействиям.
1.3. Примеры применения компьютера в моделировании
Компьютерный практикум посвящен изучению закономерностей изменения численности двух нескрещивающихся популяций и , конкурирующих за общий ресурс, путем построения математических моделей этого процесса. Особенность предлагаемого подхода состоит в том, что результаты анализа поведения моделей конкуренции, традиционно интерпретируемые в экологических терминах, здесь используются для получения дополнительных выводов эволюционного характера.
Задача формулируется следующим образом: можно ли рассматривать эволюцию фенотипа как процесс максимизации некоторого количественного критерия приспособленности, определяемого фенотипом? В некоторых модельных ситуациях процесс эволюционного отбора может интерпретироваться таким образом, хотя критерии эволюционной оптимальности могут быть различными для разных ситуаций. Решающим является механизм торможения роста численности популяции при увеличении ее плотности.
Рассмотрим несколько случаев:
- модель неограниченного экспоненциального роста;
- модель остановленного экспоненциального роста;
- фазовый анализ.
1.3.1. Модель неограниченного экспоненциального роста
В основе модели экспоненциального роста лежит простое и естественное предположение о том, что скорость роста популяции пропорциональна численности этой популяции
(8)
где - численность популяции в момент времени , а - так называемая, удельная скорость роста численности (биотический потенциал), которую можно представить как разность , удельной рождаемости b и удельной смертности d. Это уравнение описывает базовое состояние популяции при отсутствии внешних влияний, аналогичное, в некотором смысле, состоянию покоя или равномерного прямолинейного движения в первом законе Ньютона. Уравнение (8) можно легко проинтегрировать методом разделения переменных, что даёт
(9),
то есть действительно, уравнение (8) описывает динамику экспоненциально растущей популяции.
Совместная динамика численности двух экспоненциально растущих популяций (x- резидента и y- мутанта), характеризуемых удельными скоростями роста rx и ry описывается системой двух дифференциальных уравнений
(10)
Система (10) - простейшая модель конкуренции двух популяций и , различающихся своими удельными скоростями роста. Обе популяции неограниченно увеличивают свою численность независимо друг от друга. При решении вопроса об изменении их относительных численностей во времени и , с учётом (9), можно получить следующие выражения для и :
и
.
Из этих выражений следует, что если rx > ry , то стремится к 1, а - к нулю, а если rx < ry , то наоборот, стремится к нулю, а к 1.
Для графического представления можно воспользоваться пакетом программ STATISTICA. Для этого нужно создать (в режиме File - New) новый пустой файл данных из 1001 строк и 3 столбцов с именами t, ksi, eta, после чего написать (в режиме Analysis - STATISTCA BASIC) и выполнить (кнопкой Execute) следующую программу:
Sequential;
rx:=0.04; ry:=0.05; x0:=99; y0:=1;
t:=caseno-1;
ksi:=1/(1+(y0/x0)*exp((ry-rx)*));
eta:=1/((x0/y0)*exp((rx-ry)* )+1);
После этого в режиме "Graphs - Stats 2D Graphs - Scatterplots - Graph Type Multiple - Fit Off" можно построить графики зависимости переменных ksi и eta (откладываются по оси ординат) от времени (откладывается по оси абсцисс).
На рис. 16 приведен соответствующий график. Очевидно, что для rx = 0.04 и ry = 0.05 (что соответствует примерно 4-х и 5-ти процентному приросту численностей популяций X и Y за единицу времени) доля численности популяции X в конечном итоге приближается к нулю, а доля популяции Y - к единице.
Полученные результаты представляют интерес не только с точки зрения сравнительного анализа динамики численности двух экспоненциально растущих популяций, они позволяют сделать выводы относительно возможного направления эволюции в рассматриваемой ситуации.
Рис. 16. Динамика долей численности () и () экспоненциального роста популяций резидента X и мутанта Y для случая rx = 0.04,
ry = 0.05 и x0 = 99 и y0= 1.
Удельные скорости роста rx и ry, которые можно рассматривать как фенотипические характеристики популяций X и Y, при постоянных внешних условиях определяются генотипами этих популяций. Любое изменение генотипа можно рассматривать, в конечном итоге, с точки зрения его влияния на удельную скорость роста новой популяции, соответствующей измененному генотипу. Пусть популяция X является исходной (резидентной), а популяция Y возникла в результате мутации (в этом случае естественно, для наглядности, задать начальную численность x0 гораздо больше y0 , как это, в частности, сделано в примере, показанном на рис. 16). Из проведенного выше анализа следует, что если ry > rx , то доля популяции, порожденной мутантом, в суммарной численности двух популяций будет стремиться к 1, то есть мутации, увеличивающие удельную скорость роста имеют селективное преимущество при естественном отборе. Таким образом, в ситуации неограниченного экспоненциального роста удельная скорость роста популяции r (параметр Мальтуса) является естественной мерой эволюционной (дарвиновской) приспособленности. В результате естественного отбора субпопуляции (кланы) с более высокой удельной скоростью роста неизбежно увеличивают свою относительную численность, что дает основания интерпретировать эволюцию в данной постановке задачи как процесс максимизации r. Конечно, в более реалистичных постановках задачи следует вводить дополнительные условия, предотвращающие неограниченное возрастание критерия приспособленности в процессе его эволюционной оптимизации, однако на данном уровне анализа это обстоятельство не имеет принципиального значения.
Достарыңызбен бөлісу: |