Вопросы ихтиологии, 2022, T. 62, № 2, стр. 198-208

МАТЕРИАЛ И МЕТОДИКА

Для моделирования динамики запаса русского осетра Азовского моря была разработана модель DAP (Depleted artificial population model). Данная модель реализована в среде R (R Core Team, 2016) и состоит из двух подмоделей. В первой из них с помощью метода наименьших квадратов выполняется настройка модели. Используется итеративная процедура сопоставления модельного ряда оценок биомассы половозрелой части популяции (B_y) с оценками того же запаса, полученными с помощью учётных съёмок 1981–1999 гг. (B_prior(y)).

Предварительно проводится сглаживание оценок запаса по данным съёмок методом скользящего среднего по трём точкам для устранения межгодовых флуктуаций и неточностей оценивания прямым методом. В данной подмодели учитываются многолетние ряды уловов, включая как официальные (С_y), так и ряды уловов ННН-промысла (C_iuu); ежегодный выпуск молоди с ОРЗ (N_juv(y)); средняя масса половозрелых особей (W_M), рассчитанная по данным о возрастной динамике массы тела W(x) с учётом пола. Значения ННН-уловов для периода 1989–1999 гг. получены путём кусочно-линейной аппроксимации экспертных оценок за период 1994–1999 гг. (рис. 1); для периода 1981−1988 гг. условно приняли, что ННН-уловы составляли 10% официального вылова. Каких-либо оценок этого параметра в данный период в Азово-Донском бассейне, к сожалению, не имеется.

В процессе настройки модели выполняется оптимизация следующих параметров: φ_M– коэффициент годовой естественной убыли половозрелого стада; rr – коэффициент выживания выпускаемой молоди до возраста достижения половой зрелости (эквивалент коэффициента промыслового возврата для промысловых популяций); B₀– биомасса половозрелой части запаса в первый год выбранного периода настройки модели.

Расчёт биомассы пополнения половозрелой части запаса в год y выполняется по формуле:

где i – возраст достижения половой зрелости, Njuv_y−i – численность выпущенной молоди русского осетра в год (y–i).

В модели в качестве стартового значения принято rr = 0.01 (1%), что согласуется с данными литературы: 1.1–1.3% (Бойко, Калинкина, 1961; Макаров, 1964; Реков, Корнеев, 1987), 0.6% (Зайдинер и др., 2000).

Биомасса половозрелой части популяции на начало следующего года вычисляется по формуле:

где B_y – биомасса половозрелой части запаса на начало года y; C_total(y) = C_y + C_iuu−y, где C_y и C_iuu–y – временны́е ряды соответственно официальных и ННН-уловов в год y; φ_M –коэффициент годовой естественной убыли, W_M – средняя масса особей половозрелого стада, W_y – средний годовой прирост массы тела половозрелых особей.

Зависимости массы тела половозрелой особи от возраста (W1(x) − самцы в возрасте от 6 до 21 года, W2(x) − самки 10−25 лет) рассчитаны по среднемноголетним данным (Козлитина и др., 2005) с помощью линейных функций (в обоих случаях R² = 0.97): W1(x) = 0.8482x + 1.47 и W2(x) = 1.764x – 5.013. В расчётах далее принято среднее значение годового прироста самцов и самок (соответственно 0.848 и 1.764 кг): W_y = 1.306 кг.

Средний возраст наступления половой зрелости определён по кривой полового созревания русского осетра логит-методом (Михайлюк, Пятинский, 2020) и принят равным 14 лет (17 лет для самок, 12 для самцов), что согласуется с данными Чепурной и Рекова (2017). Средняя масса особи в возрасте созревания (W_M) оценена как 17 кг.

Отметим, что в модели DAP ведутся расчёты для двух периодов: первый период, выбранный для настройки модели (1981–1999 гг.), второй – прогнозный (2021–2050 гг.). Как правило, в условиях моратория на промысел между этими периодами существует промежуток лет, для которого запас моделью не оценивается.

Поскольку оценка запаса методом учётных съёмок проводится не в начале года, а в некотором месяце m, выполнен пересчёт модельного ряда B_y к моменту проведения съёмки; биомасса запаса в месяце m года y рассчитывается по уравнению:

Вторая подмодель рассматривает разные сценарии управления запасом для прогнозного периода. Расчёты биомассы половозрелой части запаса для прогнозного периода выполняются по тому же уравнению (2) с тем же набором оптимальных значений коэффициентов (φ_M, rr), но для него теперь требуется задать биомассу запаса для предыдущего 2020 г., для которого съёмка показала близкую к нулю величину, потому принято B_0(2020) = 0.

В прогностической подмодели исследованы 24 сценария восстановления запаса, рассматривающие разные сочетания двух управляющих параметров: Njuv_y−14 – численность выпуска молоди с ОРЗ с временны́м лагом, равным возрасту созревания 14 лет; biuu-coef − доля ННН-уловов от биомассы запаса в начале каждого года. До 2020 г. включительно объём выпуска молоди (Njuv) известен (данные официальной статистики ОРЗ). Начиная с 2021 г. для каждого сценария задаётся постоянное значение Njuv, которое начиная с 2035 г. будет определять величину биомассы рекрутов (B_rec) согласно формуле (1). В разных сценариях рассматриваются четыре уровня Njuv в прогнозный период: 3 млн экз. – наблюдаемый среднемноголетний уровень выпуска; 5 млн экз. – оптимистичный вариант сценария выпуска; 1 млн экз. – депрессивный вариант сценария; 7 млн экз. – максимально наблюдённый показатель в ретроспективный период.

Значения коэффициента biuu-coef задаются в диапазоне от 1 до 20% биомассы запаса на начало года – 1, 5, 7, 10, 15, 20%, что не превышает значений этого параметра в ретроспективный период (20%).

Для проверки надёжности построенной модели выполняли процедуру диагностики – тест чувствительности модели. Тест чувствительности позволяет оценить возможные отклонения локальных оптимумов от глобального оптимума оценок параметров модели (биомассы запаса, биомассы пополнения, промысловой убыли), что может приводить к ошибкам при выборе оптимального набора значений (φ_M, rr, B₀). Тест неопределённости позволяет оценить возможные отклонения оценок биомассы запаса из-за неопределённости расчётных оценок промысловой убыли в ретроспективный период.

При выполнении теста чувствительности из всех наборов значений трёх параметров (φ_M, rr, B₀) случайным образом выбирали 75% наборов. Далее для каждого выбранного набора значений (φ_M, rr, B₀) выполняется повторный пересчёт всех характеристик запаса (биомассы запаса и пополнения, промысловой убыли), после чего полученные оценки для каждого расчётного года (у) упорядочиваются по возрастанию и для каждого года определяются медианные оценки и доверительные 95%-ные интервалы. Полученные интервалы являются оценками возможной степени разброса истинной величины запаса из-за возможной неточности параметризации модели. Процедура теста чувствительности выполняется как итеративная процедура метода Монте-Карло. Полученные медианные оценки и доверительные интервалы сравниваются с модельными оценками, полученными при генеральном оптимуме (φ_M, rr, B₀) для определения степени надёжности оценок запаса.

Тест неопределённости динамической системы запас−промысел к фактору изменчивости промысловой убыли (этот тест особенно важен из-за неопределённости оценок ННН-изъятия) выполнен в соответствии с процедурой тестирования моделей в организациях ICES/FAO (FAO, 1996, 2002). Суть процедуры заключается в оценке возможной степени вариативности величины запаса в зависимости от того, какая изменчивость промысловой убыли наблюдалась в ретроспективный период. При помощи данного теста выполняется прогноз возможных изменений биомассы запаса в будущем в зависимости от того, что наблюдалось в прошлом для динамической системы запас−промысел. В качестве оценки вариативности используется значение стандартного отклонения из ряда оценок промысловой убыли в ретроспективный период, после чего с помощью коэффициента вариации эта величина пересчитывается на шкалу оценок биомассы запаса. На основе полученных оценок вариации и известных величин средних значений биомассы запаса в ретроспективный период генерируются псевдослучайные ряды оценок биомассы запаса. Полученные ряды псевдослучайных значений сортируются по возрастанию, после чего рассчитываются 95%-ные доверительные интервалы.

Хотя промысловое изъятие русского осетра в прогнозный период не предвидится, тем не менее в работе моделируется допустимое изъятие для целей научных исследований (НИР) и искусственного воспроизводства (далее – ОДУ НИР). Оценка ОДУ НИР выполнена из соображений поиска такой величины изъятия, которая не может оказать существенного влияния на половозрелую часть запаса в прогнозный период. Для этого, исходя из результатов теста неопределённости при базовом сценарии прогноза (biuu-coef = 0.1, Njuv = 3 млн экз), вычисляется определённая доля от ширины доверительного интервала оценки биомассы запаса. Этот диапазон описывает вероятную степень изменчивости величины запаса за каждый год. Полагаем, что небольшая доля изъятия от этого диапазона не должна значительно повлиять на биомассу запаса в прогнозном периоде.

Для определения незначимой доли ОДУ НИР тестировали несколько уровней доли от диапазона неопределённости (1, 5, 10, 50, 100%) и вывод о незначимости влияния ОДУ НИР делали на основе исследования его влияния на период достижения целевого ориентира управления (B_tr)_. Незначимым считался такой уровень воздействия, при котором период достижения целевого ориентира увеличивался не более чем на 1 год.

РЕЗУЛЬТАТЫ

В результате выполнения процедуры поиска оптимального решения по входным данным за период с 1981 по 1999 гг. получены оценки трёх основных параметров модели (табл. 1). На основе найденных оптимальных значений параметров модели выполнен расчёт биомассы рекрутов и запаса (рис. 2) по уравнениям (1)–(3) для ретроспективного периода. Рост коэффициента промысловой убыли (рис. 2в) в конце периода настройки модели связан как с высокими оценками нелегального изъятия (рис. 1), так и с резким снижением биомассы запаса. Для диагностики результатов моделирования выполнены тесты чувствительности (рис. 3) и неопределённости (рис. 4).

Анализ динамики запаса русского осетра, полученной в ретроспективный период настройки модели (рис. 2а), показал, что в течение благоприятного для запаса периода (1988–1998 гг.) условия среды позволяли поддерживать биомассу половозрелой части популяции на уровне не менее 10 000 т, поэтому в прогнозных сценариях решено рассматривать задачу восстановления запаса до этого уровня (B_tr = 10 000 т), который и считаем целевым ориентиром управления запасом и промыслом русского осетра. Таким образом, при моделировании возможных прогнозных сценариев запас считается восстановленным, если его биомасса достигла величины 10 000 т.

Для прогнозного периода рассмотрены 24 сценария при разном сочетании управляющих параметров: численность выпущенной с ОРЗ молоди 14 лет назад (Njuv) и коэффициент biuu-coef, постоянный для всего прогнозного периода. В качестве базового (и наиболее реального) выбран сценарий, при котором рыбоводные заводы ежегодно будут выпускать 3 млн особей (на уровне среднемноголетнего выпуска в ретро периоде), а biuu-coef = 0.1 (10% биомассы запаса на начало года). На рис. 5 представлены оценки половозрелой части запаса и биомассы пополнения для базового сценария. Возобновление коммерческого промысла русского осетра не рассматривается, однако на нужды научного изъятия и для формирования ремонтно-маточного стада осетровых заводов закладывается незначительное изъятие ОДУ НИР на уровне 1% от ширины доверительных интервалов, полученных в тесте неопределённости (рис. 6).

При базовом сценарии биомасса запаса достигнет целевого ориентира 10 тыс. т в 2048 г. и продолжит расти, если не изменится режим промысла (рис. 5). Биомасса пополнения начиная с 2035 г. стабилизируется на уровне 931 т при фиксированном годовом выпуске с ОРЗ 3 млн экз. сеголеток.

На рис. 6 представлена динамика запаса с доверительными интервалами, рассчитанными при тех же условиях, что и в ретроспективной подмодели, на основе вариативности коэффициентов промысловой убыли согласно тесту неопределённости (рис. 6б) и тесту чувствительности к стартовым параметрам на основе процедуры Монте-Карло (рис. 6а) при уровне значимости p = 0.95.

Уравнение для общего улова в прогнозный период после 2034 г. имеет вид: C_total(y) = C_y + C_iuu−y = = C_y + biuu-coef × B_y, где C_y – допустимое изъятие на цели НИР и воспроизводства на уровне, не оказывающем значимое влияние на биомассу половозрелой части запаса (рис. 7); C_iuu−y – объём ННН-добычи, представленный как доля биомассы половозрелой части запаса на начало года (C_iuu−y = biuu-coef × B_y).

Результаты проверки различных уровней ОДУ НИР в течение прогнозного периода показали, что незначимое смещение года достижения целевого ориентира выполняется при уровне 1% изъятия от ширины доверительного интервала оценок запаса в тесте неопределённости. Для предотвращения нарушения процесса восстановления запаса из-за неопределённости уровня ННН-изъятия было решено оценивать величину лимитирующего (максимально допустимого) значения ОДУ НИР для разных сценариев с учётом возможности реализации самого худшего, депрессивного, сценария. Максимальная величина ОДУ НИР для всех сценариев зафиксирована на её максимальном вычисленном для депрессивного сценария уровне: 6.13 т (biuu-coef = 0.2 и Njuv = 3 млн экз.). Эта величина соответствует максимальному значению нелинейной функции, описывающей зависимость год–вылов в прогнозный период (рис. 7б). Полученные для базового сценария значения возможного ОДУ НИР заложены во все моделируемые сценарии в качестве величины официального изъятия (C_y).

Три других сценария восстановления запаса при иных величинах Njuv, но с тем же значением biuu-coef = 0.1, что для базового сценария, представлены на рис. 8 вместе с базовым. Для этих четырёх сценариев биомасса запаса меняется одинаково вплоть до 2034 г., пока Njuv определяется реальными объёмами выпуска молоди.

Результаты прогностической подмодели для всех 24 сценариев динамики запаса: четыре варианта численности ежегодного выпуска молоди (1, 3, 5 и 7 млн экз.) при шести вариантах значений biuu-coef (0.01, 0.05, 0.07, 0.10, 0.15 и 0.20) представлены в табл. 2

ОБСУЖДЕНИЕ

Результаты моделирования в ретроспективный период (рис. 2а) свидетельствуют о плавном росте биомассы половозрелой части запаса в 1985–1995 гг., что обусловлено высокочисленными выпусками ОРЗ в предшествующий период (1971–1981 гг.). Начиная с 1995 г. отмечается стремительное сокращение биомассы запаса на фоне резкого роста промысловой убыли, которая в большей степени была обусловлена объёмом ННН-изъятия в этот период.

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

Результат теста неопределённости (рис. 4) также указывает на увеличение разброса вероятностных оценок после 1990 г. По-видимому, снижение надёжности настройки ретроспективной подмодели после 1990 г. обусловлено недостаточной информационной обеспеченностью величин ННН-промысла и многократным усилением с годами его воздействия на популяцию русского осетра. Также, вероятно, снижение надёжности получаемых оценок биомассы запаса, рекрутов и промысловой убыли после 1990-х гг. обусловлено неопределённостью и (или) временнóй изменчивостью параметров rr, φ_M в разные модельные периоды. Параметры rr и φ_M в модели полагаются неизменными. Это условие может не выполняться, однако при ограниченности данных считаем такое допущение приемлемым.

Практически все прогнозные сценарии, за исключением самых депрессивных, указывают на возможность восстановления запаса русского осетра в долгосрочной перспективе. Так, наиболее вероятный сценарий, при котором ННН-изъятие не превысит 10% биомассы запаса на начало года, а объём выпуска молоди ОРЗ останется на среднемноголетнем уровне − 3 млн экз., предполагает восстановление биомассы половозрелой части запаса до уровня 10 тыс. т к 2048 г., а при более высоком выпуске молоди (5 млн экз.) – уже к 2037 г.

При наиболее оптимистичных сценариях (при более строгом регулировании ННН-промысла и наращивании мощностей ОРЗ) возможно восстановление запаса до целевого уровня уже к 2028 г. Оптимистичными выглядят и сценарии, при которых объём выпуска молоди ОРЗ составит 7 млн экз/год – даже при величине ННН-изъятия 15% биомассы запаса возможно достижение целевого уровня B_tr.

Рассмотрены и более депрессивные сценарии. Так, в случае ослабления регулирования ННН-добычи (вариант biuu-coef = 0.20) независимо от объёма выпуска молоди ОРЗ достичь целевого уровня восстановления запаса не удастся вплоть до 2050 г. При низком уровне выпуска молоди с ОРЗ (1 млн экз/год) прогноз также неблагоприятный – целевой ориентир не будет достигнут, а биомасса запаса будет постепенно снижаться.

Результаты данной работы согласуются с выводами других авторов (Чепурная и др., 2008; Чепурная, Реков, 2017; Реков, Чепурная, 2018), которые отмечают критическое состояние запаса русского осетра в Азовском море в период запрета промысла и главенствующее негативное влияние ННН-промысла, а для восстановления запаса требуется выпуск не менее 2.8 млн экз/год. Данный объём выпуска использован в базовом сценарии, рассмотренном в нашей работе. Шляхов с соавторами (2005) в качестве негативного фактора, оказывающего главенствующее влияние на биомассу запаса русского осетра, также указывают ННН-промысел.

Моделирование возможного ОДУ на цели НИР и воспроизводства указывает на возможность изъятия по мере увеличения запаса при базовом сценарии от 4.86 т в 2022 г. до 6.13 т в 2050 г. (рис. 7а). При депрессивном сценарии изъятие для НИР и воспроизводства может составить соответственно 3.46 и 5.90 т (рис. 7б). Как показали результаты моделирования, такое изъятие не окажет значимого влияния на биомассу половозрелой части запаса в рассматриваемый прогнозный период. Эта рекомендация согласуется с рекомендациями протокола “XXXII сессии Российско-украинской комиссии по вопросам рыболовства в Азовском море”, в котором промысел русского осетра запрещён, а на цели НИР и воспроизводства для Российской Федерации и Украины предлагается ОДУ соответственно 3.81 и 0.20 т. Следует заметить, что представленные рекомендации по ОДУ НИР не будут оказывать значимого влияния на динамику запаса только в рамках моделируемых сценариев.

ЗАКЛЮЧЕНИЕ

Предлагаемая модель DAP может быть применима в случае, когда воспроизводство запаса обеспечивается только выпусками молоди, получаемой в условиях ОРЗ. Другим важным условием для использования модели является наличие оценок запаса и промысла в благоприятный период развития запаса до закрытия промысла. Модель DAP может быть применена не только для запасов осетровых, но и для других видов рыб с искусственным воспроизводством.

Следует иметь в виду, что на практике сложно оценить реальную промысловую убыль из-за неопределённости уровня ННН-промысла, а тем более его изменчивости. Фактор ННН-промысла, существенно превышающий официальные величины изъятия, вносит в результаты моделирования высокую долю неопределённости. В случае, если данные по уловам представлены надёжным временны́м рядом, при моделировании можно избежать большой неопределённости, с которой пришлось столкнуться авторам настоящей работы.

Полученные результаты свидетельствуют о депрессивном состоянии запаса русского осетра в бассейне Азовского моря в настоящее время. Рассмотренные долгосрочные сценарии включают как оптимистическое развитие событий, сопровождающееся ростом биомассы запаса, так и депрессивное, которое может повлечь за собой длительный коллапс популяции. Имеющиеся работы других авторов (Шляхов и др., 2005; Чепурная и др., 2008; Чепурная, Реков, 2017; Реков, Чепурная, 2018) подтверждают высокую вероятность депрессивного сценария.

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