Моделирование порогового эффекта в эффективности воспроизводства Волжской популяции осетра Текст научной статьи по специальности «Математика»

Водные экосистемы

УДК 629.075

МОДЕЛИРОВАНИЕ ПОРОГОВОГО ЭФФЕКТА В ЭФФЕКТИВНОСТИ ВОСПРОИЗВОДСТВА ВОЛЖСКОЙ ПОПУЛЯЦИИ ОСЕТРА

О 2014 А.Ю. Переварюха

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

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

Ключевые слова', пороговые эффекты в динамике популяций, деградация волжских осетровых, дискретно-непрерывные динамические системы

ВВЕДЕНИЕ

Методы математического моделирования популяционных процессов развиваются с 20-х годов прошлого века и при существенном интересе исследователей к данной области накоплен набор общедоступных подходов. Можно выделить в классификации большие группы «вольтерров-ских» моделей, формализующих взаимодействие конкурирующих популяций на основе систем ОДУ, методы трофо-динамического подхода предназначены для описания переноса энергии по пищевым цепям. В ихтиологии и энтомологии развивались модели, ориентированные на определении механизмов саморегуляции численности при ограничивающем влиянии выделенных факторов среды. При наличии арсенала моделей практическая организация рационально использования биоресурсов на основе прогностической математического аппарата оценивается как сложная. Поставивший под угрозу исчезновения перелов осетровых рыб Каспийского моря можно объяснить ошибками и развалом системы регулирования эксплуатации, но коллапс рыбных запасов отмечался в странах экологически ответвленным обществом и действующей научно-обоснованной системой регулирования промысла, как в Норвегии [1]. Вынужденными решениями об остановке уже и так оскудевшего промысла удается лишь более или менее своевременно предотвращать дальнейшие необратимые последствия деградации биоресурсов.

К причинам отсутствия универсальной отработанной методики можно отнести то, что построение математического аппарата моделей изначально было слабо связано с данными наблю-

Переварюха Андрей Юрьевич, кандидат технических наук, старший научный сотрудник лаборатории Прикладной информатики и проблем информатизации общества СПИИРАН, madelf@pisem.net

дений. Теоретическая экология формулируется на уровне качественных принципов, а не строгих законов. Заимствованные из механики методы описания процессов в отсутствии строго формализованных законов оказались не столь эффективны в новой области применения. Начиная с уравнения П. Ферхюльста предполагалось, что популяция характеризуется некоторым репродуктивным потенциалом, r-параметром, что выразилось в том числе в представлении о К- и г-стратегиях выживания видов. Однако, наблюдения свидетельствуют, что r-параметр постоянен только в некоторых усредненных условиях [2].

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

Динамическая сложность популяционных моделей.

Службами учета Канады неоднократно отмечался нестабильный характер воспроизводства облавливаемых популяций лососевых рыб, мигрирующих для размножения в реки. Фиксировалось отсутствие пропорциональной зависимости R » kS между нерестовым запасом и сезонным пополнением. Возникла теория о роли различных факторов смертности в зависимости от плотности рыб на нерестилищах [3]. Личиночные стадии крупных осетровых и лососевых протекают при высокой плотности, даже если взрослых рыб не так много. Канадский биолог У.Е. Рикер, рассматривая причины уменьшения пополнения при возрастании нерестового стада, предложил теоретическую концепцию о существовании зависимо-

548

Известия Самарского научного центра Российской академии наук. 2014. Т. 16, № 5(1)

сти запаса и пополнения, и роли её свойств при долговременной эксплуатации [4]. Рикер предложил простую, как казалось, формулу расчета величины пополнения R = AS):

f(S) = aSexp(-bS\ (1)

где S - от английского «stock» величина нерестового запаса, Ъ - коэффициент, отражающий величину, обратную количеству производителей, при котором число выжившей молоди максимально, имеет смысл b = 1; а - параметр характеристика репродуктивного потенциала. График зависимости пополнения от численности производителей называется кривой пополнения. (1) задает куполообразную кривую с асимптотой lims f (S) = 0 и единственным нетривиальным

пересечением с биссектрисой координатного угла R=S.

Очевидно, что при применении модели (1) возникает ряд сложностей. При увеличении количества отложенной икры выживаемость молоди должно уменьшаться до неразличимых значений, что противоречит наблюдениям над аквариумными популяциями. Пополнение уменьшается до некоторого Q>0 = const. (1) при моделировании на ЭВМ сообщества гидробионтов часто использовалась как кусочно-гладкая, до некоторого критического количества Sc, при S > Sc количество молоди предикативно устанавливалось Q = const.

При низких численностях производителей (1) предсказывает увеличение скорости воспроизводства. При стремящемся к нулю количестве отложенной икры, эффективность пополнения стремится к предельному значению:

= a&cp(-bS)(\-bS), lima&cp(-bS)(\-bS) = a. (1a)

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

Модель (1) обладает интересными свойствами, о которых не имел представления ее разработчик и потому не предложил для них адекватной биологической интерпретации.

Представим процесс изменения состояния популяции, определяемый зависимостью запаспополнение в виде динамической системы, математического объекта, для которого можно указать набор величин: динамических переменных, характеризующих состояние системы. Значения таких переменных в последующий момент времени рассчитываются из текущих значений по закону, называемому оператором эволюции. Динамическая система представляется тройкой

(м,г,у, состоящей из фазового пространства М, времени Т, оператора эволюции v|/. Для VX еМ и^^еГ выполняется условие:

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

Множество (yw(x)}teJ, - фазовая траектория

точки х. Графически эволюция диссипативной динамической системы во времени представляется движением точек в фазовом пространстве к аттрактору, подмножеству фазового пространства А ^ М инвариантное относительно эволюции:

у w (Д) = А дня всех t еТ и такое, что существует окрестность U множества А. в которой для всех у е U выполняется lim y(i) (_у) = А. Регулярными аттракторами для отображения интервала у 1®1 считаются состояние равновесия с

неподвижной точкой х*: ПтУ\у) = х и устой-

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

Стремящиеся к одному аттрактору множество точек есть его область притяжения, которая может быть непрерывным подинтвалом I, так и образовывать фрактальное множество [5]. Рассмотрим динамическую систему как полугруппу

итераций {у(;'1};>о, и пусть i?o, i?j, i^,... после_

довательность точек эволюции системы, определённых условием Rj+i = \|/(i?; ) при всех j > 0.

Качественное поведение зависит от параметра а, он является управляющим. До определённого значения а, не превышающего бифуркационное, система стремится к точечному аттрактору R* = In ajb . Первый метаморфоз поведения

системы происходит, когда производная, вычисленная в неподвижной точке, перестаёт удовлетворять критерию устойчивости |у'(^ )| < 1. Для (1) это происходит при выполнении условия

а > е2:

\|/ '(R) = аеГьк -bRaebR,

, -Ь— In а -Ъ—

Чj\R ) = ae 6 ь

1-1па = -1, а=е2.

а( 1 - In а)

1-1па.

Параметр Ъ не влияет на топологические характеристики фазового портрета. Динамическая система стремится в устойчивое циклическое состояние с периодом 2: глобальный аттрактор, состоящий из двух периодических точек. Факт, что

—Ъх

хп+1 = ахпе " имеет в большом диапазоне пара-

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

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

549

Водные экосистемы

вого значения а> 12.51, произойдёт бифуркация удвоения периода и установится цикл р — 4. При дальнейшем увеличении а все быстрее реализуется каскад бифуркаций удвоения периода. Около а> 14,8 невозможно выделить устойчивых точек или замкнутого цикла. Происходит детерминированный хаос, напоминающий стохастический процесс очень чувствительный к начальным условиям. Траектория притягивается к подмножеству фазового пространства, получившему из-за своих фрактальных свойств название «странный аттрактор». Важно знать, что сценарием хаотиза-ции через бесконечный каскад удвоений, иссле-дованым на примере х;?+1 = цхп (1 ~Хп)

М. Фэйгенбаумом в [6], спектр нелинейных эффектов не ограничивается. В ограниченных диапазонах при а> 14,8 наблюдаются «окна периодичности» устойчивых циклов р Ф 2\ далее вновь испытывающих удвоения.

_Ъх

В хп+х = ахпе 11 наблюдаем последовательность касательных бифуркаций с появлением устойчивых циклов периодов р Ф 2\ Ли и Йорк в 1975 г. показали, что если одномерное отображение Rj+l = V|/(i? .) при некотором значении параметра а = b имеет цикл периода р — 3, то оно при b имеет и бесконечное множество циклов всех других периодов. Ранее А. Шарковский в [7] доказал, если одномерное отображение f: I ® I имеет цикл р = п, то оно также имеет циклы со

всевозможными периодами предшествующими числу /7, среди целых чисел, выписанных в порядке Шарковского, который завершает число три.

Рис. 1. Перемежаемость с циклом р = 3 в

Xn+1 = аХп*ЬХП ■

Появлению окон предшествует явление «перемежаемости» (отмечен стрелкой на рис. 1). Перемежаемостью считается режим чередования во времени почти регулярных колебаний (ламинарная фаза) с интервалами хаотического поведения (турбулентная фаза). Перед появлением циклических точек ветви графика третьей итерации функции /3(х) образуют так называемые каналы, в

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

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

Проблема интерпретации нелинейной ни-намики.

Появление хаоса интерпретировалось в биологическом контексте как последствие повышения репродуктивного потенциала популяции. В (1) это а, но для динамики подобных унимодальных

отображений при мы увидим набор резких

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

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

Не для всех рыб актуальна теория пополнения, в гораздо большей степени ее методы обоснованы для анадромных видов, заходящих из морей для нереста на речном грунте.

Эффект резкого сокращения эффективности воспроизводства.

Изучение ретроспективных данных о популяциях анадромных рыб Нижней Волги позволило выявить нелинейные особенности в формировании пополнения. Общей задачей наших исследований является анализ на основе вычислительных модельных сценариев антропогенно-

обусловленных процессов деградации биоресурсов экосистемы Волго-Каспия. К прогнозированию уловов осетровых методы динамической теории пополнения У. Рикера ранее не применялись, основным фактором считался гидрологический режим реки.

Среди разнообразной опубликованной информации интересны сведения о количестве молоди осетра, полученной от некоторой запаса, обработанные нами методом скользящей средней (рис.

2). Количество пропущенных на нерест в Волгу производителей за период наблюдений изменя-

550

Известия Самарского научного центра Российской академии наук. 2014. Т. 16, № 5(1)

лась существенно: от 2500 тыс. экз. в 1979-81 гг. до 25 тыс. в 2000-х гг. [8]. В настоящее время мы не имеем достоверной информации о нынешнем состоянии популяции Acipenser gueldenstaedtii и численности нерестовой части. С 2010 г. волжские популяции осетровых внесены в «Красную книгу» и промысел запрещен. Имеющиеся сведения с 1968 г. об общей численности осетровых Каспия по данным тралового учета [9] многими специалистами оцениваются как сомнительные [10] и менее достоверные, чем информация о проходе производителей на нерестилища в створе реки.

Рис. 2. Данные о воспроизводстве осетра

Acipenser guldenstaedtii.

Кривая воспроизводства осетра (на рис. 2 среднее по трем ближайшим) характеризуется как изломанная куполообразная. Мы можем восстановить отчетливо гладкую кривую, а не сгущение точек. Следовательно, популяционная динамика определялась промысловой смертностью более, чем зависимостью f(S). Стрелкой отмечен момент: при пропуске менее 350-103шт. непропорционально резко уменьшается количество скатывающейся в море молоди. Выявлен статистически достоверный пороговый эффект скачкообразного перехода эффективности воспроизводства на некий отличный от нуля минимальный уровень.

Подобное явление оказалось характерно не только для единственной волжской популяции. Сам У. Рикер не рассматривал (1) как единственно возможную нелинейную форму зависимости между запасом и пополнением, описав в работе [11] неожиданный для промысла эффект резкого падения воспроизводства, нарисовав для его объяснения кривую воспроизводства как подобие ступенчатой функции. Рикер назвал эффект «лезвием ножа», т.к. даже небольшой перелов 5% может подорвать запасы такой популяции и потому оптимальная стратегия промысла, изымающая весь излишний прирост нерестового стада, крайне опасна для практики в плане стабильной долговременной эксплуатации.

Зоологами неоднократно подтверждалось для не заботящихся о потомстве видов, что существует оптимальный для воспроизводства диапазон

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

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

Разработка новой репродуктивной модели.

Наличие ограниченных пищевых ресурсов в неявном виде учитывалось известными моделями, но общим недостатком функциональных итераций можно считать игнорирование резкого проявления декомпенсационного фактора смертности и стремительного изменения пищевых потребностей по мере развития молоди. Декомпенсационные факторы - увеличивающие смертность при уменьшении плотности особей снижают эффективность нереста, уменьшая количество икры реально вступившей в репродуктивный процесс при S = R . Прямой учет такой зависимости становится очевидной необходимостью при моделировании динамики популяций подвергающихся перелову. Задачей разработки новой модели запаспополнение стало создание гибкого математического аппарата для согласования качественного характера поведения имитационной модели промысловой популяции с ретроспективными данными о численности молоди волжских анадром-ных рыб.

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

Формализуем увеличение пищевых потребностей молоди исходя из факта, что скорость роста находится в обратной зависимости от численности поколения, но не в обратно пропорциональной. Это согласуется с данными наблюдений биологов в частности, с результатами экспериментов над ростом рыб при различной их плотности. Согласно наблюдениям, при увеличении плотности возникает асимметричное распределение размерной структуры популяции в сторону преобладания особей с меньшими размерами. Описывать убыль численности поколения N(t) на интервале модельного времени [0,7] будут следующие объединённые в систему дифференциальные уравнения (а, Р, с - константы):

dN

dt

dw

dt

-(a w(t)N(t) + ®(S)$)N(t)

f—r—------, w(°) = w0, JV(0) = XS.

551

Водные экосистемы

где: S— величина нерестового запаса; wit) - отражает уровень размерного развития поколения, влияющий на увеличение пищевых потребностей; g - параметр, учитывающий ограниченность количества доступных для молоди кормовых объектов. Z _ параметр ограничивающий темп развития вне зависимости от N(t) ; 1 - средняя плодовитость особей; а - коэффициент компенсационной смертности; (3 - коэффициент декомпенсационной смертности; к £ (1,2); интервал уязвимости t £ [О, Т\ должен быть специфичным для вида.

У бывающая функция

0(A) = 1 + ехр(-сУ2), lims^oo 0(A) ® 1 отражает

действие важного в экологии «эффекта Олли» или «эффекта агрегированной группы», когда для животных ведущих стадный образ жизни отмечается увеличение выживаемости молоди, хотя скученность усиливает пищевую конкуренцию взрослых. Для промысловых рыб подобное явление актуально [12] и достоверно наблюдалось [13], когда численность нерестового стада низка. Протяженность русловых нерестилищ осетровых обуславливает снижение плотности.

Функция эффекта Олли ©(А) введена в уравнение для N(t), но не в переопределение начальных условий N(0) = v(XS + 5) т.к. даже не вступившие в репродуктивный процесс кладки икры после их гибели ухудшает токсикологическую обстановку на нерестилище. с -характеризует степень действия эффекта, который не проявляется если численность взрослой части популяции осетровых достаточно велика. Введя поправку 5 для учета незначительной сейчас доли повторно нерестующих рыб, численно получаем S = N(T).

Графиком исследованной с применением численного решения в инструментальной среде последовательности задачи Коши (2) для вычисления численности выживших особей поколения к моменту времени Т новой модели является унимодальная кривая с уменьшающимся наклоном ниспадающей правой ветви. Кривая имеет ненулевую горизонтальную асимптоту, две нетривиальные точки пересечения с биссектрисой координатного угла R=S (рис. 3). Фазовое пространство в случае динамической системы \|/” (Rtl) на

основе системы (2) разделяется на две области притяжения: ГГ и ГК Границей между ними является неустойчивая стационарная точка R] первого пересечения кривой с биссектрисой R=S.

__ ^

Второе равновесие /у гибридной динамической системы устойчиво.

Анализ устойчивости неподвижных точек дискретно-непрерывной системы в инструментальной среде можно проводить с использовани-

ем свойства второй итерации y(y(x)). Необходимым и достаточным условием устойчивости неподвижной точки одномерного отображения

х*является неравенство: y2(x) > X при х<х*,

У О) < X при х>х .

=

темы (2). (2) исследовалась в среде моделирования AnyLogic-6, предоставляющей выбор численных методов с изменяющимся шагом интегрирования.

ЗАКЛЮЧЕНИЕ

Поведение траектории динамической системы, с использованием в качестве оператора эволюции решение задачи Коши (2), качественно отличается от итераций (1) возможностью притяжения к двум аттракторам. Соответственно присутствуют две области притяжения, границей которых служит репеллер, неустойчивая особая точка первого пересечения кривой с биссектрисой координатного угла. Траектории с начальными условиями, разделёнными репеллером, покидают его окрестность и приближаются к разным аттракторам. Один из аттракторов точка (0,0) на плоскости R/S. Если численность популяции соответствует области притяжения этого аттрактора, произойдёт вымирание. Данное свойство отражает наличие минимально допустимой численности для существования вида.

Кривая вида рис. 3 означает, что популяция долго способна выдерживать высокий процент изъятия, но потом произойдет резкое снижение уловов, которое сложно предвидеть экспертам и ввести мораторий на промысел. Деградация может растянуться во времени, и численность популяции будет колебаться ниже неустойчивого равновесия. В подобной ситуации необходимо своевременно оценивать качественный характер динамики запасов.

Работа выполнена при поддержке РФФИ (Проект 14-07-00066).

СПИСОК ЛИТЕРАТУРЫ

1. Риске Я., Шулдалъ Х.Р., Шлагстаг Д. Экологическое

моделирование для рыболовства// Модели многовидового управления [под ред. Т. Редсета]. М.: Изд-во

ВНИРО, 2002. С. 17-79.

2. ОдумЮ. Экология. М.: Мир, 1986. 376 с.

3. 1Reave F. Principles affecting the size of pink and chum

salmon population in British Columbia // Journal of the

552

Известия Самарского научного центра Российской академии наук. 2014. Т. 16, № 5(1)

fisheries research board of Canada, 1953. Vol. 9. № 9. P. 450—491.

4. Ricker W. Stock and recruitment // J. Fisheries research

board of Canada, 1954. №.11. C. 559-623.

5. Переварюха А.Ю. Об определении фрактальных объек-

тов в динамике моделей управления биоресурсами // Труды СПИИРАН, 2013. Т.24. С. 211-221

6. Feigenbaum М. Universal behavior in nonlinear systems //

PhysicaD., 1983. Vol. 7. №1-3. P. 16-39.

7. Шарковский A.H. Сосуществование циклов у непрерывного преобразования прямой в себя // Украинский математический журнал. 1964. Т.26. №1. С. 61-71.

8. Вещее П.В., Гутенева Г.И. Эффективность естественно-

го воспроизводства осетровых в Низовьях Волги в современных условиях // Экология. 2012. № 2. С. 123127.

9. Панин Г.Н., Мамедов РМ., Митрофанов И.В. Совре-

менное состояние Каспийского моря. М.: Наука, 2005. С. 281.

10. Легеза М.И., Маилян Р.А. К методике исследования численности осетровых в Каспийском море // Тезисы отчетной сессии ЦНИОРХ. Астрахань, 1974. С. 82-83.

11. Ricker W.E. Big effect from small causes: two examples from fish population dynamics // Journal Fisheries research board of Canada. 1963. Vol. 20. №2. P. 257-264.

12. Stoner A.W., Davis M.H., Booker C.J. Negative consequences of Allee effect are compounded by fishing pressure: comparison of queen conch reproduction in fishing grounds and a marine protected area // Bulletin of marine science, 2012. Vol. 88. №1. P. 89-104.

13. Myers R.A., Barrowman N.J., Hutchings J.A., Rosenberg A.A. Population-dynamics of exploited fish stocks at low population levels // Science, 1995. Vol. 269. P. 1106— 1108.

SIMULATION OF THE THRESHOLD EFFECT IN EFFICIENCY OF REPRODUCTION OF RUSSIAN STURGEON IN THE VOLGA RIVER

© 2014 A.U. Perevaryukha

St. Petersburg Institute for Informatics and Automation of RAS

Acceptable balance between reproduction and mortality from various factors could dramatically disturbed because of the impermanence in reproductive potential, that showed a retrospective analysis of data on spawning sturgeon in the Volga. Article considers development of the model of population process with the existence of different levels of reproductive efficiency. We proposed a system of ODE, which formalizes the influence of the growth rate on mortality in juvenile period and take into account the loss of unfertilized eggs. In a computing environment was received discrete-continuous dynamical system which allowed assessing the effects of a sharp decrease in reproductive efficiency.

This effect leads to the presence minimum acceptable size for the survival of Acipenser gueldenstaed-tii spawning population.

Key words, threshold effects in the dynamics of populations, degradation of the Volga sturgeon, discrete-continuous dynamical systems

Perevarukha Andrew Y., Ph.D., senior researcher at the Laboratory of Applied Informatics and problems of information society SPIIRAS, madelf@pisem.net

553