«Динамика популяции камчатского краба (Paralithodes camtschaticus) в Баренцевом море (опыт моделирования) ...»
ПОЛЯРНЫЙ НАУЧНО-ИССЛЕДОВАТЕЛЬСКИЙ ИНСТИТУТ МОРСКОГО
РЫБНОГО ХОЗЯЙСТВА И ОКЕАНОГРАФИИ ИМ. Н.М. КНИПОВИЧА
(ПИНРО)
На правах рукописи
Баканев Сергей Викторович
Динамика популяции камчатского краба (Paralithodes
camtschaticus) в Баренцевом море
(опыт моделирования)
Специальность 03.00.18 – Гидробиология
Диссертация на соискание ученой степени кандидата биологических наук
Научный руководитель – доктор биологических наук, профессор А. В. Коросов Мурманск – 2009 Содержание Введение …………………………………………………………………………………...
Глава 1. Биология камчатского краба и рациональное использование его запасов (обзор литературы) ………………………………………………………………………. 1.1. Интродукция ………………………………………………………………… 1.2. Биология …………………………………………………………………….. 1.3. Промысел ……………………………………………………………………. 1.4. Опыт моделирования динамики популяции и оценка параметров …….. 1.5. Принципы рационального управления запасами камчатского краба ….. Глава 2. Материал и методы ………………………………………………………….. 2.1. Данные по траловым съемкам ……………………………………………… 2.2. Промысловые данные ……………………………………………………….. 2.3. Моделирование динамики численности …………………………………… 2.4. Оценка параметров моделей ……………………………………………….. Глава 3. Популяционные характеристики камчатского краба в Баренцевом море … 3.1. Распределение и размерно-возрастной состав самцов …………………… 3.2. Репродуктивные параметры ………………………………………………… 3.3. Групповой рост камчатского краба ………………..……………………….. Глава 4. Оценка численности популяции камчатского краба в Баренцевом море … 4.1. Продукционная модель……………………………………………………... 4.2. Когортная модель LBA ……………………………………………………… 4.3. Когортная модель CSA ……………………………………………………… 4.4. Сравнительный анализ результатов моделирования ……………………… Глава 5. Динамика численности камчатского краба в Баренцевом море …………… Глава 6. Трофическое воздействие камчатского краба на биоту Баренцева моря….. Глава 7. Регулирование промысла камчатского краба в Баренцевом море ………… Выводы………………………………………………………………………………….
Список литературы……………………………………………………………………
СПИСОК СОКРАЩЕНИЙ
ОДУ (общий допустимый улов) – величина годового промыслового изъятия части популяции, рассчитанная с учетом биологических особенностей данной популяции (продуктивности, динамики численности и др.) и целей ее эксплуатации: соответствует оптимальной, с точки зрения выбранного критерия регулирования, интенсивности промысла.BMSY – численность (экз.), при которой популяция достигает своей наибольшей продуктивности.
Bt – численность животных (экз.) в году t;
Cl,t – улов особей (экз.) размерной группы l в году t.
CSA (catch survey analysis) – разновидность модели динамики популяции, основанная на изменениях численности более крупных размерных групп, чем LBA.
F – коэффициент мгновенной промысловой смертности, отношение числа особей, изъятых промыслом за единицу времени, к численности популяции в течении этого времени, при условии, что численность популяции за это время не меняется.
Gl – средний прирост (мм) особи камчатского краба для размерной группы l.
K (емкость среды) – максимальная биомасса (численность) популяции, которую может поддерживать среда обитания данной популяции в отсутствии промысла.
LBA (length-based analysis) – модель динамики популяции ракообразных, основанная на изменениях численности размерных групп.
MSY (maximum sustainable yield) – максимальный устойчивый улов, соответствующий уровню максимальной продуктивности популяции; наибольший уравновешенный улов, который теоретически можно получить от данной популяции;
широко используется в качестве критерия регулирования промысла.
M – коэффициент мгновенной естественной смертности, отношение числа особей, погибших за единицу времени по естественным причинам, к численности популяции в течении этого времени, при условии, что численность популяции за это время не меняется.
ml,t – вероятность линьки для краба размерной группы l в году t.
Nl,t – численность краба (экз.) с новым карапаксом размерной группы l в году t.
Ol,t – численность краба (экз.) со старым карапаксом размерной группы l в году t.
POt (пострекруты) – численность крабов (экз.), пополнивших промысловый запас в году t–1 и ранее. Пострекруты это все промысловые крабы со старым карапаксом, а также особи с новым панцирем, имеющие размеры большие, чем рекруты.
PRt (пререкруты) – пополнение (экз.) промыслового запаса в году t+1. Пререкруты принадлежат размерной группе, конец интервала которой граничит с минимальным промысловым размером, а длина интервала соответствует средней длине прироста пререкрута за одну линьку.
q – коэффициент улавливаемости или доля популяции, которая вылавливается единицей промыслового усилия.
REt (рекруты) – численность (экз.) крабов с новым карапаксом, пополнивших промысловый запас в году t. Рекруты принадлежат размерной группе, начало которой равно минимальному промысловому размеру, а длина интервала размерной группы соответствует средней длине прироста рекрута за одну линьку.
ВВЕДЕНИЕ
Камчатский краб (Paralithodes camtschaticus Tilesius, 1815) – один из важнейших промысловых объектов Баренцева моря. Интродуцированный в 60-х годах прошлого столетия, он успешно акклиматизировался в водах северной Атлантики и к настоящему моменту образовал самовоспроизводящуюся популяцию на акватории от Лофотенских островов на западе до архипелага Новая Земля на востоке (Беренбойм, 2003). В 2004 г., благодаря значительному увеличению численности популяции, российский флот начал коммерческую эксплуатацию запаса. По нашей экспертной оценке в 2006 г. вылов камчатского краба в этом регионе достиг 30% от мирового улова этого вида.На протяжении последних 15 лет Полярный научно-исследовательский институт морского рыбного хозяйства (ПИНРО) совместно с норвежским Институтом морских исследований (IMR) ведут регулярные исследования этого вида в новой для него среде обитания. Результаты исследований показали, что популяция камчатского краба в Баренцевом море достигла численности, при которой возможна его успешная промысловая эксплуатация (Камчатский краб…, 2001). Однако, по мнению российских и норвежских ученых акклиматизация краба в Баренцевом море к настоящему моменту не завершена и в разных частях современного баренцевоморского ареала находится на разных ее этапах (Анисимова, 2003; Nilssen, 2003).
При взрывном характере роста численности, который часто наблюдается у акклиматизированных животных, камчатский краб может стать пищевым конкурентом для некоторых бентофагов, в том числе промысловых рыб. Качественный и количественный состав пищи камчатского краба в новом районе обитания в целом сходен с таковым в местах его естественного распространения. Отличительной особенностью спектра питания краба в Баренцевом море является сравнительно высокая массовая доля иглокожих, в то время как в тихоокеанском регионе как по частоте встречаемости так и по массе доминируют моллюски (Герасимова, 1997). Являясь пищевым конкурентом, популяция вселенца значительно влияет на биологическую структуру экосистемы и ее продуктивность. Кроме того, распространение этого вида несет определенные факторы риска, которые условно можно подразделить на экологические (последствия для окружающей среды и традиционно населяющих ее организмов) и экономические (последствия для традиционных промыслов и марикультуры) (Левин, 2001).
Несмотря на весьма длительную историю акклиматизации и биомониторинга камчатского краба, многие вопросы его биологии важные для промысла остаются малоизученными. В литературе отсутствуют оценки абсолютной численности и продуктивности популяции, а использующиеся характеристики носят относительный характер (Камчатский краб.., 2003). Такие традиционные представления о динамики численности, пространственной и функциональной структуре популяции основаны на большом объеме эмпирической гидробиологической информации. Между тем, математические методы анализа накопленных данных могут дать существенно более детальное описание процессов адаптации популяции к новой среде обитания в абсолютных оценках. На основе таких моделей возможно организовать рациональную эксплуатацию промысловой популяции интродуцента.
Приемы стохастического моделирования динамики, основанные на аналитическом выражении элементарных характеристик систем, в настоящее время являются эффективным инструментом исследования. Применение этих методов позволяет получать не только количественные оценки, но и выявлять внутренние и внешние взаимосвязи, управляющие ходом популяционных изменений (Deriso et al., 1985; Васильев, 2001).
Успех такого рода исследований в значительной степени зависит не от объема и качества информационной базы, а от интеграции и последующей унификации всей имеющейся в распоряжении разнородной информации, что позволяет получить достоверные знания о количественных характеристиках объекта и факторах, влияющих на динамику численности.
К сожалению, изучение формирования баренцевоморской популяции камчатского краба связано с определенными трудностями. Оценка динамики численности ракообразных, в отличие, например, от рыб, сопряжена с трудностью определения точного возраста. Камчатский краб растет в течение короткого периода после линьки, т.е.
индивидуальный рост особей имеет сложный прерывистый характер. Частота линьки зависит от пола, возраста и условий обитания краба. Традиционные модели с использованием численности поколений становятся трудно применимыми (Smith, Addison, 2003).
Одной из главных проблем в определении продуктивности баренцевоморской популяции камчатского краба является его незаконный промысел, который не поддается точной оценке. По экспертным оценкам российских исследователей масштабы браконьерского промысла позволяют сделать вывод об отсутствии реального контроля над выловом крабов в водах России (Цыгир, 2006).
В настоящей работе осуществлена попытка учесть все вышеизложенные аспекты и адаптировать новейшие вероятностно-статистические методы к существующим подходам исследования популяционных характеристик: численности, размерно-возрастной структуры, репродукции, смертности, кривой роста, плотности распределения и пр. Цель работы: на основе математического моделирования оценить динамику популяции камчатского краба в процессе освоения им нового ареала и разработать рекомендации по рациональной эксплуатации его запаса.
Решались следующие задачи:
1) Сформировать репрезентативные массивы эмпирических данных по биологии и промыслу камчатского краба в Баренцевом море на основе обобщения всей имеющейся информации (в том числе собранной автором) о состоянии его популяции.
2) Оценить популяционные характеристики по данным натурных съемок и литературным данным (параметры размерно-возрастной структуры, скорость роста, смертности, репродукции и пр.).
3) Разработать серию имитационных моделей и с их помощью оценить популяционные параметры и реконструировать популяционную динамику;
4) В серии модельных экспериментов оценить устойчивость популяции краба к интенсивности промысловой эксплуатации.
5) Оценить кормовую обеспеченность популяции камчатского краба как фактора, влияющего на динамику его численности.
6) Разработать рекомендации по рациональной эксплуатации запаса камчатского краба в Баренцевом море.
Автор благодарен всем сотрудникам Полярного института, помогавшим ему при сборе данных в ходе морских экспедиций и при обработке полученных материалов.
Особую благодарность автор выражает сотрудникам лаборатории промысловых беспозвоночных ПИНРО: Б. И. Беренбойму, Н. А. Анисимовой, П. А. Любину, М. А.
Пинчукову, В. А. Павлову, И. Е. Манушину и Ю. И. Жаку. Автор искренне признателен коллегам по лаборатории математического обеспечения оценки запасов гидробионтов:
Ю.А. Ковалеву, В. А. Коржеву, В. Л. Третьяку. Особую признательность автор выражает своему научному руководителю, Коросову Андрею Викторовичу, вложившему много труда и оказавшему неоценимую помощь на всех этапах подготовки работы.
Биология камчатского краба и рациональное использование его Идея вселения камчатского краба в Баренцево море возникла в конце 20-х годов прошлого столетия в Государственном океанографическом институте и в Тихоокеанском институте рыбного хозяйства. В 30-е годы была начата разработка метода практического осуществления эксперимента и в это же время предприняты попытки вселения этого вида.
После первых неудачных экспериментов опыты по перевозке были продолжены в 50-х годах, но и они не послужили основой для будущей интродукции краба в баренцевоморские воды.
Работы по вселению краба возобновились в 60-е годы, когда после обращения Мурманского Совнархоза к отечественным специалистам было подготовлено биологическое обоснование на вселение животного в Баренцево море. Руководителем и организатором этой работы стала доктор биологических наук А.Ф. Карпевич.
Практические работы по перевозке крабов на Кольский полуостров начались в 1961 г.
Крабов отлавливали в районе залива Петра Великого Японского моря и доставляли к Баренцеву морю самолетами в каннах (рис. 1). Начиная с 1966 г. было решено перевозить крабов в живорыбных вагонах. За весь период интродукции (1961–1969 и 1977–1978 гг.) в Баренцево море было выпущено 1,525 млн. экз. личинок, 10,7 тыс. экз. молоди и около тыс. экз. взрослых крабов (Камчатский краб…, 2003).
Сведения о расселении этого вида в Баренцевом море стали поступать в ПИНРО с середины 70-х годов. На основании фактов неоднократной поимки крабов специалисты сделали вывод о том, что интродуценты могут жить в условиях Баренцева моря (Козлов, Строгонова, 1977; Орлов, 1978; Сенников, 1977; Orlov, Ivanov, 1978). С конца 70-х стала поступать информация о попадании крабов в орудия лова не только у берегов СССР, но и в водах Норвегии (Орлов, 1997).
Рис. 1. Тихоокеанский и баренцевоморский ареалы камчатского краба (пунктирная стрелка – маршрут транспортировки краба в 1961–1969 гг.) В 1991 г. на российско-норвежской встрече ученых в г. Мурманске было решено организовать систему регулярных наблюдений, включающую фиксацию границ распределения камчатского краба, оценку его численности, определение полового, размерно-возрастного состава и других биологических показателей популяции. К середине 90-х годов прошлого столетия по результатам экспедиций 1992–1993 гг.
учеными ПИНРО был выдвинут тезис об окончании периода акклиматизации камчатского краба в Баренцевом море и формировании самостоятельной самовоспроизводящейся популяции (Герасимова, Кузьмин, 1994).
До формирования баренцевоморской популяции камчатский краб считался эндемиком северной части Тихого океана. В Охотском море обитают две крупные популяции краба: западно-камчатская и аяно-шантарская. В Беринговом море обитает крупная бристольская, а в заливе Аляска – аляскинская популяции (рис. 1). Кроме них известны сравнительно малочисленные популяции Японского моря: приморская, западносахалинская, южно-курильская, хоккайдская популяции и др. (Родин, 1985).
Расширение баренцевоморского ареала началось регулярно фиксироваться российскими и норвежскими учеными с начала 90-х годов прошлого столетия (рис. 2, Герасимова, Кузьмин, 1994). До этого периода краб встречался в основном в прибрежной полосе от Варангерфьорда до мыса Святой Нос. С началом регулярных исследований краб значительно расширил свой ареал как на запад, так и на восток от мест выпуска, причем темпы увеличения ареала вида в Норвежском море до 2003 г. превышали таковые в восточной части Баренцева моря (Беренбойм, 2003). Основываясь на устных сообщениях руководителя исследований камчатского краба в Институте Морских Исследований (г.
Тромсё, Норвегия) Яна Сундета, с каждым годом ученые Норвегии расширяют акваторию исследовательской съемки, включая в зону наблюдения все новые и новые фьорды.
Рис. 2. Расселение камчатского краба (A – границы ареала до 1977 г.; B – до г.; C – до 1994 г.; D – до 1997 г.; E – до 2007 г.; звездочка – район вселения краба) и основные районы его промысла (1 – Варагерфьрд; 2 – Рыбачья банка; 3 – Мотовский залив; 4 – Кильдинская банка; 5 – Западный Прибрежный район; 6 – Восточный Прибрежный район; 7 – Мурманское мелководье; 8 – Канинская банка) в Баренцевом море (Кузьмин, Беренбойм, 2000; Sundet, Pinchukov, 2009).
В результате исследований камчатского краба в Баренцевом море в 90-е годы установлено, что характер его сезонных миграций очень схож с таковыми в нативном регионе обитания (Камчатский краб…, 2003). В Баренцевом море максимальная глубина, на которой были отмечены камчатские крабы во время нагульных миграций, составляет 300–350 м. Это самые глубокие места у Мурманского побережья, где происходит, очевидно, и зимовка крабов, после чего они совершают нерестовые миграции в мелководные части губ и заливов (Беренбойм, 2003). Зоной основного воспроизводства Дальнезеленецкой. Формирование плотных скоплений крупных самцов в юго-восточных районах Баренцева моря происходит преимущественно благодаря активной миграции особей из западных районов (Соколов, Милютин, 2008).
Нерест протекает в период с января по май, как правило, на глубинах не превышающих 100 м. Пик нереста приходится на март-апрель и зависит от температурного режима конкретного года. Сравнительный анализ основных репродуктивных параметров показал, что в настоящее время баренцевоморская популяция камчатского краба, являясь самой высокоширотной, тем не менее успешно воспроизводится и не испытывает заметного угнетения со стороны факторов внешней среды. По ряду репродуктивных параметров (максимальные размеры, максимальная абсолютная индивидуальная плодовитость, размер 50% созревания) самки камчатского краба Баренцева моря несколько превосходят самок тихоокеанского региона (Баканев, 1997, Баканев и др., 1997; Баканев, Матьков, 1997; Баканев, 2001; Матюшкин, 2003).
Массовый выклев личинок камчатского краба в прибрежных районах Баренцева моря отмечается в апреле. Личинки краба могут находиться в планктоне около 3– месяцев, мигрируя в Прибрежной ветви Мурманского течения в восточном направлении.
Завершая планктонный цикл, мальки оседают в июле-августе, преимущественно в прибрежной зоне (Матюшкин, 2003; Баканев, 2003).
Зона массового расселения ранней молоди камчатского краба ограничивается прибрежными мелководьями с глубиной 0–20 м. Молодь обитает в тех же биотопах, что и многие мелкие формы баренцевоморских десятиногих раков, предпочитая участки распространения твердых грунтов и зарослей макрофитов, где она находит богатую кормовую базу и укрытие от хищников. В распределении молоди по мере ее роста появляется тенденция расширения зоны обитания. Достигнув 4–6 летнего возраста животные начинают миграцию в сторону открытого моря плотными скоплениями, пополняя промысловые группировки (Переладов, 2003).
С выходом в открытую часть моря камчатские крабы совершают сезонные миграции, нагуливаясь, в основном, на глубинах 150–200 м. Для взрослых особей характерна значительная пищевая пластичность и довольно низкая избирательность пищевых объектов. В целом спектр питания крабов определяется набором наиболее массовых и доступных в районе откорма бентосных организмов и других пищевых объектов (мелкие моллюски, усоногие раки, раки-отшельники, офиуры, морские ежи, полихеты, остатки рыбного промысла). Переход к факультативной некрофагии и включении в основной рацион отходов рыбного промысла в виде остатков рыб и беспозвоночных происходят лишь при недостатке его основного естественного корма – донных беспозвоночных (Герасимова, 1997; Герасимова, Кочанов 1997; Манушин, Анисимова, 2003).
Камчатский краб, как и все ракообразные, растет в течение короткого периода после линьки, до образования твердых покровов. Молодь камчатского краба линяет несколько раз в год; с возрастом частота линьки сокращается. Результаты исследований роста камчатского краба в Баренцевом море показали, что величина прироста у крабов с длиной карапакса 60–180 мм практически одинакова и в среднем составляет 17 мм с диапазоном значений от 10 до 23 мм по длине карапакса (рис. 3).
Рис. 3. Данные о приросте длины карапакса (А, сплошная линия – линия тренда) и вероятности линьки самцов камчатского краба (Б, сплошная линия – моделируемая логистическая кривая, пунктирная – порог 50%-й вероятности) (Nilssen, Sundet, 2006).
В возрасте 10–12 лет самцы камчатского краба достигают промыслового размера (150 мм по ширине карапакса и 132 мм по длине). В этот период частота линьки самцов краба сокращается и в среднем составляет 1 раз в два года. Достигая размеров свыше мм по ДК, самцы линяют раз в 3–4 года. Максимальный возраст, до которого доживает камчатских краб, около 20 лет, но обычно с учетом промысла, не более 15 лет (Matsuura, Takeshita, 1990). Принимая во внимание, что смертность не может быть одинакова в разные годы и у разных возрастов, американскими исследователями неоднократно совершались попытки оценить ее с помощью данных по мечению (Hirschhorn, 1966;
Balsiger, 1974) или по данным траловых съемок (табл. 1, Reeves, Marasco, 1980; Reeves, 1988; Greenberg et al., 1991). Различия в коэффициентах естественной смертности связаны, по-видимому, с различными методами их оценки и фактическими данными. Оценки настолько вариабельны, что достаточно сложно сделать какие-либо выводы о надежности этих величин. Тем не менее, могут быть сделаны два заключения: 1) естественная смертность заметно выше у особей младших и старших возрастных групп, чем у крабов средних возрастов; 2) в разные годы исследований наблюдаются значительные колебания естественной смертности.
Коэффициенты мгновенной естественной смертности (М)* для самцов камчатского краба тихоокеанского региона Reeves, Marasco, 1980(данные Reeves, Marasco, 1980 (данные Reeves, 1988 (данные 1969–1980) 0,58 0,70 0,71 0,34 0, Reeves, 1988 (данные 1981–1986) 1,21 0,72 0,49 0,65 1, Greenberg et al., 1991 (данные Greenberg et al., 1991 (данные * при отсутствии промысловой смертности равен M= ln( Nt 1 / Nt ), где Nt – численность поколения в году t, Nt 1 – численность в году t+1.
С конца 90-х годов прошлого столетия помимо естественных факторов на выживаемость популяции камчатского краба в Баренцевом море существенно влияет антропогенное воздействие (Беренбойм и др., 2006; Баканев, 2009; Bakanev, Berenboim, 2006). В соответствии с решением 22 сессии Совместной российско-норвежской Комиссии по рыболовству с 1994 г. начался экспериментальный промысел камчатского краба в Баренцевом море. С 2004 г. в российских водах Баренцева моря начат промышленный лов. Историю российского промысла краба в Баренцевом море можно разделить на три периода:
1994–1997 гг. – начальный период;
1998–2003 гг. – переходный период;
с 2004 г. – период промышленного лова.
Промысел в начальный период осуществлялся с помощью японских конических ловушек и донных тралов с сентября по декабрь. Добычу вели от 2 до 6 судов разных типов, включая береговые катера, а также средние рыболовные траулеры. Годовой вылов не превышал 20 тыс. экземпляров. Основные районы промысла – Западный Мурман (от Варангерфьорда до Кильдинской банки).
С 1998 г. при промысле камчатского краба начали использовать американские краболовы-процессоры. Добычу осуществляли в основном японскими коническими ловушками и американскими прямоугольными ловушками. Ежегодный вылов возрос с в 1998 г. до 600 тыс. экз. в 2003 г. Район промысла значительно расширился за счет акваторий Восточного Мурмана (Восточный и Западный Прибрежные районы, Мурманское мелководье). Впервые в 2002 г. значительная часть рыболовного флота ( судов) дислоцировалась в российской экономической зоне за пределами территориальных вод. В российской прессе начинает появляться обеспокоенность о фактах бесконтрольного промысла баренцевоморского камчатского краба. По устным сообщениям наблюдателей количество неучтенного краба на судах составляло 50–80 % от выгруженного груза.
Период промышленного лова характеризуется значительным увеличением промысловых нагрузок в 2005–2006 гг. Количество судов возросло до 30. Основная акватория промысла – районы Восточного Мурмана. Оцененный фактический вылов составил около 3 и 5 млн. экз. в 2005 и 2006 гг., соответственно. Импорт камчатского краба в США возрос в 2 раза. Американские аналитики и журналисты, а также российские ученые (Цыгир, 2006) поддерживают гипотезу о бесконтрольном лове камчатского краба в российских водах Баренцева моря.
При организации промысла краба в Баренцевом море был использован опыт управления в его нативном ареале. Схема такого управления совершенствовалась в течении нескольких десятилетий. К началу экспериментального промысла в Баренцевом море были подготовлены рекомендации по его регламентации, которые сводились к главным регулирующим подходам (Герасимова, Кузьмин, 1997):
2) регулирование вылова на основании общего допустимого улова Регулирование по размеру/полу сводится к набору простых технических мер, направленных на сохранение нерестового запаса и пополнения. Это полный запрет лова самок, а также самцов менее 150 мм по ширине карапакса. Лов должен производиться ловушками, поскольку в этом случае в каждом улове краб сортируется по полу, размерам, и часть выловленных особей можно легко отпустить в целях сохранения пополнения и репродуктивного потенциала. Кроме того, ограничение на лов в прибрежных районах связан с охраной молоди и самок, чья плотность на акваториях с глубинами менее 100 м выше, чем самцов промыслового размера. Сезонные ограничения связаны с товарным качеством краба в период его линьки и направлены на сохранение промыслового запаса.
Регулирование на основании общего допустимого улова было введено с началом экспериментального промысла в 1994 г. Первоначальный выбор уровня рекомендованного ОДУ (25% от оцененного промыслового запаса) для баренцевоморской популяции камчатского краба в начальный период ее исследований был сделан с учетом дальневосточного опыта. Предполагалось планомерно изменять степень промысловой эксплуатации новой популяции, в зависимости от биологических и экономических откликов, с тем чтобы найти наилучший количественный выбор эмпирическим методом проб и ошибок. Так в 2001 г., учитывая высокую неопределенность в динамике интродуцированного краба, отсутствие какого-либо опыта в эксплуатации такого рода биоресурсов в мире и начало коммерческого промысла, был рекомендован более осторожный 20% уровень изъятия. Одновременно в ПИНРО была начата разработка количественных моделей с целью получения наилучших возможных альтернатив для имеющихся в распоряжении рядов наблюдений. Предполагалось совершенствовать алгоритмы расчетов, модифицировать их по мере поступления новой информации и на основании этих моделей улучшать рекомендации по выбору ОДУ.
1.4. Опыт моделирования динамики популяции и оценка параметров Для оценки биологических показателей популяции камчатского краба в Баренцевом море в 2003 г. была осуществлена первая попытка моделирования его динамики численности. Модель описывала элементарные биологические процессы и базировалась на традиционных биостатистических методах, разработанных для оценки эксплуатируемых популяций рыб. Рассматривая динамику численности одного поколения, профессор Ф.И. Баранов вывел в 1918 г. соотношение, на которое опирается в настоящее время большинство биостатистических моделей (Баранов, 1918):
где Nt – численность поколения в начальный момент времени t, Nt 1 – численность в момент времени t+1, e – основание натуральных логарифмов. Показатель степени Z называется коэффициентом мгновенной общей смертности (Рикер, 1979; Шибаев, 2007).
Смысл коэффициента: за элементарный (очень маленький) промежуток времени dt численность животных уменьшается на величину dN, равную Z-той части от фактической численности N или dN/dt=–ZN. Для эксплуатируемых популяций гидробионтов Z определяется изъятием животных промыслом (коэффициентом мгновенной промысловой смертности, F) и всеми остальными причинами, кроме промысла (то есть коэффициентом Спустя 40 лет Gulland (1964) и Murphy (1965) предложили метод анализа структуры уловов эксплуатируемых популяций рыб на основе подхода Ф. И. Баранова.
Использовались два уравнения. Первое отражает взаимосвязь численностей смежных возрастных групп:
Второе уравнение позволяет рассчитать величину улова в поштучном выражении, получаемого от каждой возрастной группы:
следующем t+1 году, когда возраст рыбы увеличится на единицу; CFt,a – улов возрастной группы a в году t;
смертности возрастной группы a в году t. Они могут быть зависимыми или независимыми от возраста, а величина F может изменяться по годам в связи с изменениями интенсивности промысла. Уравнения служат базовыми для группы методов широко используемых в настоящее время в промысловой ихтиологии, которая называется виртуально-популяционным анализом (Рикер, 1979; Jones, 1974; Murphy, 1965; Pope, 1972).
При виртуально-популяционном анализе, как правило, начальная численность Nt, определяется величиной пополнения Rt,1, то есть количеством особей пополнивших в году t промысловый запас (Рикер, 1979). Обозначим количество животных, погибших по естественным причинам, символом согласно уравнениям [2] и [3], можно представить графически (рис. 4).
Рис. 4. Блок-схема моделирования динамики численности эксплуатируемой популяции (сплошные стрелки – переход выживших особей в следующую возрастную группу, пунктирные стрелки – убыль животных от промысловой и естественной смертности).
ракообразных на основе возрастных когорт не получил широкого применения в силу сложности определения точного возраста животного (ICES, 2001a; ICES, 2001б). Однако с начала 90-х годов прошлого столетия, основные принципы изменения численности поколений, разработанные Ф. И. Барановым, стали с успехом применяться в моделировании динамики численности крабов и креветок (Collie, 1991; Kruse and Collie, 1991; Zheng et al., 1995; Quinn et al., 1998; Cadrin, 2000). Модели этого типа базировались на изменениях численности не поколений, а размерных групп (LBA, length-based analysis) и учитывали линьку и стохастический рост животных.
Рассмотрим принцип моделирования динамики численности эксплуатируемой популяции камчатского краба на основе модели LBA, описанной в работе J. Zheng и др.
(1995). Численность краба моделировалась раздельно по размерным группам l=1,2,3,…i и по полу. Каждая размерная группа l в году t состояла из особей с новым ( Nl,t ) и старым ( Ol,t ) карапаксом. В течение года t крабы, как с новым, так и старым панцирем могут линять с вероятностью ml,t или с вероятностью 1– ml,t не линяя, оставаться в этой же размерной группе.
Вероятность линьки для краба ml,t была представлена логистической функцией (Balsiger, 1974):
Нелиняющие крабы Nl,t и Ol,t в году t+1 становились Ol,t 1 в таком соотношении:
Линяющие крабы в t году – ml,t ( Nl,t Ol,t ), в году t+1 в зависимости от длины прироста попадают в последующие размерные группы и становятся крабами с новым карапаксом, принадлежащими размерным группам: 1, 2, 3,…, L. Вероятность перехода в определенную группу зависит от параметров прироста краба за одну линьку. Рост на индивидуальном уровне характеризуется длиной прироста особи. Длина прироста особи на групповом уровне варьирует как внутри размерной группы, так и между группами. Средний прирост особи камчатского краба для размерной группы l ( Gl ) обычно выражали линейной функцией от средней длины ( ) размерной группы l перед линькой краба (Weber, Miahara, 1962):
где a и b – параметры линейной функции.
Отклонения прироста особей размерной группы l за линьку от уравнения регрессии описывали через вероятностное распределение. Плотность вероятности величины прироста выражали через гамма-распределение, форма которого, в зависимости от распределении прироста (Sullivan et al., 1990):
двух параметров G и после линьки в последующую группу l, равна интегралу функции гамма-распределения на интервале группы l :
– средняя длина группы l, l1 – начало интервала группы l, l2 – конец интервала группы l. Для последней размерной группы L, PL,L 1.
Пополнение (R) моделируемой популяции камчатского краба американские исследователи описали с помощью двух параметров: 1) количество молодых крабов, пополняющих популяцию в году t, Rt и 2) параметр Ul, определяющий стохастическое изменение процесса пополнения (Sullivan et al., 1990; Zheng et al., 1995). Пополнение каждой размерной группы Rl,t выражалось через уравнение:
где Ul описывалось гамма-распределением, таким же как в уравнениях [7] и [8] с поэтому требовалось оценить только параметр r, при условии что был известен средний размер пополнения.
Учитывая описанные выше вероятностные процессы в динамике численности камчатского краба, моделируемую популяцию удобно представить в каждый момент времени вектором, состоящим из численности размерных групп:
где Nt – вектор L размерных групп, показывающий общее количество животных с новым карапаксом в году t. Расчет Nt 1 обеспечивает матрица переходов (transition matrix), содержащая параметры стохастического роста (m, P) и пополнения (R):
где Al,t – временная переменная равная:
где Nl,t и Ol,t численности крабов с новым и старым карапаксом в размерной группе l в году t, соответственно; Cl,t – вылов крабов (экз.) размерной группы l в году t;
Ml,t –коэффициент мгновенной естественной смертности для крабов размерной группы l в году t и yt – времення задержка от момента съемки до средней даты периода промысла.
Расчет переменной Al,t основывается на подходе Ф. И. Баранова, описанном уравнениями [2] и [3]. Он также схематично представлен на рисунке 4.
На основе модели LBA в 90-е годы прошлого столетия был разработан целый ряд моделей, успешно использующих возможности анализа размерных рядов (Collie, 1991;
Kruse and Collie, 1991; Zheng et al., 1995; Quinn et al., 1998; Cadrin, 2000). Одна из них, CSA (catch survey analysis), является частным случаем LBA, когда численность популяции мала для получения качественных данных по размерному составу. Модель была разработана Зенгом и др. (Zheng et al., 1997) для небольших запасов краба с относительно низкими уловами в ходе ежегодных съемок. Обе модели принципиально сходны. Отличие состоит лишь в том, что в модели CSA используются более укрупненные размерные группы, а в расчетах приемлемо использовать данные по размерному составу лишь половозрелых самцов. В модели CSA деление на размерные группы, как правило, основывается на особенностях биологии и промысла дальневосточных крабов семейства Lithodidae, в которое входит и камчатский краб. Промысловая часть популяции в году t, состоящая из крупных самцов, делится на две группы: рекруты ( REt ) и пострекруты ( POt ). REt – это численность крабов с новым карапаксом, пополнивших промысловый запас в году t.
Рекруты принадлежат размерной группе, начало которой равно минимальному промысловому размеру, а длина интервала размерной группы соответствует средней длине прироста рекрута за одну линьку. POt – численность крабов, пополнивших промысловый запас в году t–1 и ранее. Пострекруты это все промысловые крабы со старым карапаксом, а также особи с новым панцирем, имеющие размеры большие, чем рекруты. Пополнением промыслового запаса считаются пререкруты ( PRt ), особи, которые в следующем году станут рекрутами. Пререкруты принадлежат размерной группе, конец интервала которой граничит с минимальным промысловым размером, а длина интервала соответствует средней длине прироста пререкрута за одну линьку. Динамику моделируемой популяции можно представить в виде уравнений:
где mPR – вероятность линьки для пререкрутов; GPR,RE, GPR,PO – параметры роста, соответствующие доле крабов, которые переходят из группы PRy в группу REt и из группы PRt в группу POt 1 ; M – коэффициент естественной смертности; Ct – вылов крабов в году t; yt – времення задержка от момента съемки до средней даты периода промысла.
При оценках численности беспозвоночных наряду с моделями LBA и CSA, основанными на размерных группах, в 70-х годах прошлого столетия начали широко использоваться продукционные модели, базирующиеся на динамике промыслового запаса без учета возрастной или размерной структуры популяции.
При дефиците биологической и промысловой информации моделируемая численность может быть выражена упрощенно через продукционные способности популяции. Основная идея продукционного подхода сводится к следующей схеме (Шибаев, 2007):
1) в отсутствии промысла результаты роста биомассы популяции (пополнение плюс весовой рост) уравновешиваются потерями из-за естественной смертности. Численность запаса стабилизируется на некотором уровне (K), определяемом экологической емкостью среды (Одум, 1986). Промысел, как дополнительная причина смертности, сокращает численность популяции и тем самым нарушает сложившееся равновесие. Продуктивность популяции (запаса) возрастает благодаря высвободившимся кормовым ресурсам и стремится вернуть запас в утраченное равновесное состояние (Баранов, 1918, 1925);
2) если прибавочную продукцию уровнять величиной промыслового изъятия, то запас сохранится в этом новом состоянии равновесия, которое будет соответствовать данной интенсивности промысла (Russel, 1931).
Здесь под состоянием запаса (популяции) мы понимаем его комплексную характеристику, которая включает оценки важнейших параметров популяции – численности, смертности, пополнения и др. и тенденций их изменения (Левин, Коробков, 1998). Классической продукционной моделью, ставшей прообразом всех современных моделей этого класса, считается линейная модель Шефера (Schaefer, 1954, 1957). В качестве продукционной функции запаса в модели используется логистическое уравнение популяционного роста (Verhulst, 1938). Дискретная форма записи модели при отсутствии промысла имеет вид:
где Bt – биомасса животных в году t, K – максимальная возможная биомасса при отсутствии промысла (емкость среды), r – коэффициент мгновенного популяционного роста (внутренняя скорость роста) в отсутствие плотностной регуляции. Согласно свойству этого логистического уравнения максимальные продукционные способности запаса наблюдаются при биомассе равной K/2. Годовая прибавочная продукция запаса ( Yt, прирост биомассы запаса за вычетом потерь по естественным причинам) выражается как зависимость его биомассы:
Причем считается, что запас находится в состоянии статического равновесия с окружающей средой (рис. 5) Рис. 5. Равновесная продукционная кривая при K = 1000 т и r = 0,5.
Такая форма продукционной кривой обусловлена влиянием лимитирующих природных факторов, действие которых на запас возрастает с увеличением биомассы (плотности) запаса. Это явление известно как плотностная регуляция запаса (Одум, 1986).
В качестве примера приведем воображаемый эксперимент, поясняющий механизм плотностной регуляции, описанный в книге В.К. Бабаяна «Предосторожный подход к оценке общего допустимого улова» (2000). Поместим несколько пар половозрелых особей рыб в замкнутый водоем, где имеются все условия для их роста и размножения. Первое время, пока особей в водоеме еще немного, совокупный прирост биомассы будет незначительным (рис. 6). Однако с каждым нерестом темпы увеличения численности, а следовательно и биомассы, будут возрастать. Процесс ускоренного роста будет продолжаться до момента t * (на рис. 6 он соответствует 11-ому году), когда его начнет сдерживать ограниченная емкость среды обитания, которая характеризуется, прежде всего, величиной кормовой базы, а также другими лимитирующими факторами, например объемом жизненного пространства. На графике динамики роста экспериментального логистической кривой B B(t ) (рис. 6). Очевидно, что в условиях равновесия биомасса, при которой запас достигает своей наибольшей продуктивности, то есть генерирует максимальную прибавочную продукцию, соответствует 500 т. В промысловой биологии максимальную прибавочную продукцию (maximum sustainable yield) принято обозначать аббревиатурой MSY, а соответствующую ей биомассу запаса – BMSY.
Рис. 6. Динамика роста экспериментального запаса.
После того как начнет сказываться влияние лимитирующих плотностных факторов, темпы роста численности (биомассы) запаса постепенно уменьшатся до 0, а биомасса достигнет своего максимального уровня (К), определяемого емкостью среды.
Возвращаясь к рисунку 5, где MSY = 125 т, BMSY = 500 т, представим, что рассматриваемый запас стал объектом промыслового использования. Сохраняя сделанное ранее допущение о равновесном состоянии запаса, будем считать промысел фактором окружающей среды, с которым у запаса мгновенно устанавливается равновесие. При постепенном увеличении уровня эксплуатации биомасса запаса будет уменьшаться, вызывая соответствующий рост его продуктивности, пока величина прибавочной продукции не достигнет своего наибольшего значения. При дальнейшем увеличении промысловой эксплуатации биомасса, а вместе с ней и продуктивность запаса снизится до 0, и запас потеряет свое промысловое значение.
Выполненный анализ продукционной зависимости приводит к трем важным выводам, на которых строятся общие принципы управления промысловыми запасами (Шибаев, 2007):
1. Запасом (точнее, его продуктивностью) можно управлять с помощью регулирования интенсивности промысла.
2. Поддерживая промысловое усилие на определенном уровне, можно получать устойчивые уловы, равные соответствующему значению прибавочной продукции. В этом контексте термины «максимальный устойчивый улов» (MSY) и «максимальная прибавочная продукция» полностью тождественны.
3. Один и тот же устойчивый улов можно получить при двух разных значениях биомассы (ниже и выше уровня BMSY, см. рисунок 5), поэтому экономически выгодней не допускать сокращения биомассы запаса ниже уровня BMSY.
Кривая прибавочной продукции (рис. 5) количественно описывает продукционную способность популяции. Поэтому ее часто используют для обоснования стратегии долговременного промыслового использования запаса, а также мер регулирования, направленных на реализацию этой стратегии. В настоящее время, подавляющее большинство таких стратегий основывается на допущениях о равновесном состоянии запаса и концепции прибавочной продукции, разработанной в середине прошлого века (Рикер, 1979; Бивертон-Холт, 1969; Засосов, 1970; Бабаян, 2000).
Процесс оценки продукционных свойств эксплуатируемой популяции можно представить в виде двух взаимосвязанных основных этапов:
1) построение модели исследуемого явления;
наблюдения.
Основы первого этапа, описанные в предыдущей главе, заложены к середине прошлого века и не претерпели существенных изменений (Рикер, 1979, Одум, 1986).
Второй этап опирается на методы математической статистики, которые в настоящее время продолжают активно развиваться (Punt and Hilborn, 1997; Schnute et al., 1998).
Стохастическая природа оцениваемых параметров может быть описана с использованием различных алгоритмов идентификации, наиболее распространенным из которых является метод максимального правдоподобия (Hilborn and Walters, 1992).
Точность получаемых оценок параметров зависит от полноты и качества исходной информации. Фрагментарность данных, ошибки измерений, шум природных процессов в совокупности с короткими рядами наблюдений могут значительно исказить результаты оценок и, соответственно, понимание реальных биологических явлений (Васильев, 2001).
В условиях плохой информационной обеспеченности оценки параметров зачастую находят, используя алгоритмы на основе формулы Байеса, когда в качестве исходной информации берутся не только данные наблюдений, но и априорное (предварительное) знание о параметрах модели (Bayes, 1763). Подход основан на попытке начать статистический вывод с некоторых исходных предположений (догадок) о вероятностном распределении неизвестных параметров. Значения модельных параметров можно установить априорно на основе оценок, полученных, например, для одних и тех же видов, но из разных районов, или для схожих видов из одного района. Используя теорему Байеса, окончательные (апостериорные) значения параметров оцениваются с учетом как данных наблюдений, так и предварительно заданных значений параметров:
Здесь распределение P(H ) называется априорным распределением вероятностей возможных значений параметра или параметров H (это распределение принимается прежде, чем получены статистические данные); data – статистические данные, полученные в ходе эксперимента (наблюдений, опытов), которые используются в модели;
P(H | data) – условная апостериорная вероятность, количественно оценивает насколько модель соответствует гипотезе H в условиях, характеризуемых данными наблюдений статистическом оценивании параметров, которая определяет вероятность того, что полученные данные в ходе эксперимента соответствуют гипотезе H.
В рамках байесовского подхода, гипотезы (параметры) трактуются как случайные величины и описываются с использованием вектора непрерывных значений параметров вместо набора дискретных значений. Вероятности в этом случае будут выражаться через плотности вероятностей, а знаменатель P(data) в уравнении [17] согласно формуле полной вероятности будет равен Функцию Байеса можно представить иначе, в виде трех составляющих:
где 1/P(data) некая постоянная величина по отношению к гипотезе H. Поэтому апостериорное распределение пропорционально произведению функции правдоподобия на априорное распределение или:
где символ означает пропорциональность.
Рассмотрим методику применения байесовского подхода на примере оценки параметров продукционной модели Шефера, выше описанной уравнением [15]. С добавлением показателя вылова модель приобретает вид (Schaefer, 1954):
где Bt – численность в году t, Ct – вылов в году t, K – максимально возможная численность популяции и r – коэффициент мгновенного популяционного роста.
Ненаблюдаемая переменная Bt может быть выражена через наблюдаемый показатель исследовательской съемки:
где q – коэффициент пропорциональности, в данном случае коэффициент улавливаемости, a e – остаточная погрешность, имеющая логнормальное распределение (Haddon, 2001).
также начальное значение численности B0. Традиционно в популяционной биологии оценка параметров осуществляется либо методом наименьших квадратов, либо методом максимального правдоподобия (Hilborn and Walters, 1992). Как показано выше, функция правдоподобия P(data | H ) также входит в состав формулы Байеса [17] и выглядит так:
Функция правдоподобия определяет вероятность получения данных наблюдений (data) при определенном наборе параметров, значения которых соответствую гипотезе H, то есть r, K, q, и B0. Метод максимального правдоподобия заключается в подборе таких максимальной. Чем выше вероятность, тем лучше модель описывает данные наблюдений, а, следовательно, динамику численности популяции. При дефиците исходной информации результаты настройки модели могут привести к неверным выводам о численности животных. Например, если дать нереально высокую величину K или q, и тем самым значительно исказить расчетную численность. Для предотвращения таких артефактов необходимо внести дополнительную информацию, то есть некое априорное (доопытное) знание.
Обычно до начала оценки эксперт имеет, по крайней мере, приблизительное представление о возможных величинах параметров. Эти представления могут базироваться на информации, полученной из литературных источников или личного опыта эксперта. Рассмотрим ход рассуждений о возможной величине, например, коэффициента улавливаемости q.
Отображая соотношение реальной величины численности и наблюдаемого индекса относительной численности, он может быть предварительно оценен, посредством полевых экспериментов. Являясь коэффициентом пропорциональности для индекса численности, оцененной в ходе исследовательской съемки, этот показатель определяется, главным образом, уловистостью орудий лова. Доля пойманных животных из числа находящихся в зоне облова, естественно, может принимать значения от 0 (никто не пойман) до (отловлены все особи). На основании этой информации эксперт может допустить, что действительное значение коэффициента q равновероятно может находится в диапазоне от 0 до 1.
Допустим, что у эксперта имеются некая дополнительная информация в виде литературных данных о коэффициенте уловистости аналогичного трала по отношению к животным исследуемого нами вида. В результате он может прийти к выводу о том, что интервальная оценка коэффициента уловистости трала должна лежать, например, в пределах 0,4–0,6. Кроме того, из многолетнего опыта эксперта и литературных данных известно, что коэффициент уловистости варьирует от траления к тралению и зависит от множества факторов (параметров раскрытия трала, типа грунта, размеров животных и пр.). На основании этой информации, а также используя соображения о характере варьирования этого показателя, исследователь имеет некое мысленное представление о возможной величине этого параметра до начала моделирования процесса. Ожидается, что, используя традиционный метод максимального правдоподобия и данные наблюдений, величина q окажется, скорее всего, близкой к 0,4–0,6. В тоже время маловероятно, что значение этого коэффициента будет близким, например, к 0 или к 1.
Чтобы избежать некорректной оценки q, очевидно, необходимо задать ее пределы, в которых будет проходить поиск оптимального решения. Такие пределы задаются, как правило, с учетом предыдущей накопленной информации, то есть априори (в данном случае на основании литературных данных и личного опыта эксперта). Для того чтобы использовать эту информацию, наше мысленное представление о возможной величине q необходимо выразить количественно. Задача заключается в построение частотного распределения вероятностей этой величины, которое называется распределением параметра q.
Построение априорного распределения целесообразно начинать с выбора типа распределения. Однако на практике не всегда имеется исчерпывающая информация о том, какому закону подчиняется вероятностное распределение того или иного параметра. В настоящее время проблема выбора априорных распределений для параметров моделируемых биологических процессов освещена подробно (Kass, Wasserman, 1996; Punt and Hilborn, 1997). В данном случае, учитывая литературные данные, предполагаем, что распределение подчиняется нормальному закону с модой в интервале 0,4–0,6 и границами от 0 до 1. Графически такое распределение изображено на рисунке 7, а математически принято записывать как:
где dnorm – нормальное распределение, а символ “~” обозначает “распределяется как”. Параметры нормального распределения mu (среднее) и tau (точность, показатель обратно пропорциональный дисперсии, то есть ) указаны в скобках. При включении априорного распределения q в формулу Байеса [17] рассчитывается апостериорное распределение этого параметра. Иными словами наша предварительная или априорная оценка параметра q вновь оценивается с учетом данных наблюдений, априорное распределение корректируется и строится окончательное апостериорное распределение (рис. 7). На этом процесс оптимизации или настройки параметра q заканчивается. В данном случае, была получена уточненная оценка q с медианой равной 0,38.
Рис. 7. Распределение вероятности параметра q (черная линия – априорное, синяя – апостериорное).
Апостериорные распределения всегда находятся в границах своих априорных показывает, как много информации о значениях неизвестного параметра содержится в статистических данных наблюдений при условии, что априорное распределение было задано в широком диапазоне возможных значений. Чем выше и уже апостериорный пик, тем параметр точнее оценивается с помощью данных наблюдений. Если апостериорная вероятность слабо отличается от априорной, то, возможно, данные не несут полезной информации для корректировки первоначальной оценки параметра (рис. 8).
Рис. 8. Распределение вероятности параметра q (черная линия – априорное, синяя – апостериорное) при информативных (А) и малоинформативных (Б) данных наблюдений.
Наиболее важным и одновременно сложным является вопрос выбора априорного распределения параметров. Для определенных параметров выбор может основываться в некоторых случаях на субъективном мнении эксперта. В идеале это мнение должно быть основано на фактах. Эти факты могут включать все те сведения, которые обычно не входят в формальный анализ количественных данных, а используются как информация, полученная вне эксперимента. Такая «внешняя» информация может включать:
1) результаты исследований по той же тематике;
2) результаты исследований, в которых изучались сходные биологические механизмы;
3) результаты лабораторных исследований, посвященные природе изучаемого явления;
4) результаты исследований, если эти результаты могут иметь ту же природу;
5) сведения о промежуточных исходах, которые наблюдались в данном эксперименте и свидетельствуют в пользу предложенной гипотезы;
6) информацию, полученную у других популяций со схожими биологическими и промысловыми параметрами.
Эти виды внешней информации предполагают ту или иную форму экстраполяции причинно-следственных связей. При этом параметр должен имеет априорное распределение, включающее значения искомого апостериорного распределения. На практике иногда желательно обрезать или укоротить апостериорные распределения с целью избежать растягивания процесса поиска решения. В таких случаях лимиты распределения должны быть выбраны в широких, биологически обоснованных пределах, чтобы не перекрывать апостериорное распределение. Целесообразно рассматривать несколько возможных априорных распределений различных параметров и анализировать результаты расчетов (Баканев, 2007; Баканев, 2008).
Согласно формулам [20] и [21] для настройки всего моделируемого процесса должны быть заданы априорные распределения не только для параметра q, но и для r, K, и B1. Далее на основе априорных распределений пяти параметров моделируется общее априорное распределение. Затем, используя формулу Байеса и данные наблюдений, рассчитывается целевое или общее апостериорное распределение, статистические показатели которого являются искомыми оценками параметров модели.
Рассмотренный пример байесовского подхода на примере продукционной модели представляет собой принципиальную схему оценки параметров модели. Реализация формулы Байеса, когда все параметры и расчетные переменные трактуются как случайные величины, основана на поиске решения системы интегральных уравнений. Процедура поиска решения связана с большим количеством вычислений и базируется на итеративном методе с применением генератора случайных чисел. В настоящей работе используется подход, в англоязычной литературе известный как Markov Chain Monte Carlo (MCMC, Congdon, 2001). Этот алгоритм вычисления основан на построении цепи Маркова.
Процесс, протекающий в системе, называется марковским, если для каждого момента времени которой вероятность любого состояния системы в будущем определяется только состоянием системы в настоящий момент и не зависит от того, каким образом система пришла в это состояние (Вентцель, 2003). Реализация цепи Маркова в численном виде происходит с помощью модели Гиббса (Gibbs sampler). Она обеспечивает способ выборки из совместных распределений многомерных переменных с помощью применения такого положения: для получения выборки из совместного распределения делаются многократно выборки из его представленных одномерных условных распределений. Процесс реализации MCMC с помощью модели Гиббса подробно описан в литературе (Gilks et al., 1996; Meyer, Millar, 2000).
История применения подходов на основе теоремы Байеса в популяционной биологии насчитывает порядка 15 лет (Burgner, 1991; Schnute, 1994; McAllister, Ianelli, 1997; Punt, Butterworth, 1996; Angel et al., 1994; Kinas, 1996; Meyer and Millar, 1999;
Patterson, 1999). В настоящее время такие подходы широко используются для описания состояния и прогнозирования динамики запасов морских млекопитающих, рыб и беспозвоночных. Модели, использующие теорему Байеса, были успешно адаптированы для оценки запасов северной креветки, синего краба, канадского лосося, южноафриканского анчоуса, усатых китов, тихоокеанских тунцов, палтуса и сельди (Михеев, 2003; Meyer, Millar, 1999; Butterworth, Punt, 1995; Schnute et al., 1998; Millar and Meyer, 2000; Bakanev, Berenboim, 2007). Оценки состояния запасов и прогнозы, выполненные такими моделями, используются в различных международных организациях, участвующих в формировании рекомендаций по эксплуатации гидробионтов: ФАО, АНТКОМ, НАФО, ИКЕС, НЕАФК (ICES, 2008; Punt, Hilborn, 2001;
McAllister and Kirkwood, 1998; Hvingel and Kingsley, 2006; Adkison and Peterman, 1996;
Annala, 1993; Beddington and Cooke, 1983; Bergh and Butterworth, 1987; Meyer and Millar, 1999; Bakanev, 2008).
1.5. Принципы рационального управления запасами камчатского краба эксплуатируемой популяции, которые во взаимодействии с рыболовством образуют систему «запас-промысел», являются информационной основой рационального рыболовства. Исторический опыт показывает, что неконтролируемый промысел рано или поздно переходит в фазу экономического перелова, а в худшем варианте и к перелову биологическому, который может иметь необратимый характер (Хилборн, Уолтерс, 2001).
Следовательно, обеспечение рационального рыболовства невозможно без применения соответствующих мер регулирования, направленных на поддержание системы запаспромысел на некотором оптимальном уровне.
Для обоснования стратегии долговременного промыслового использования популяции, а также мер регулирования, направленных на реализацию этой стратегии, в настоящее время широко используют концепцию прибавочной продукции, описанную в предыдущем разделе (рис. 5). Если известна цель эксплуатации данного запаса, например максимизация среднемноголетнего вылова или среднемноголетней прибыли, эта цель выражается в биологических и промысловых терминах, например BMSY. Параметры, характеризующие выбранную цель, называются критериями регулирования, или целевыми ориентирами управления (FAO, 1993). В зависимости от состояния запаса, при котором начата реализация стратегии его рационального использования, управление может быть направлено, на его восстановление (если B< BMSY ), снижение (если B> BMSY ) или поддержание на исходном уровне (если B BMSY ). Целенаправленное воздействие на запас осуществляется с помощью научно обоснованных мер регулирования рыболовства, к которым относятся ограничения на различные составляющие промысловой деятельности: вылов, сроки и районы лова, селективность орудий лова и др. Критерии регулирования при этом выполняют роль своеобразных точек отсчета, относительно которых определяется удаление текущего состояния запаса от целевого и оценивается величина промыслового воздействия для внесения необходимых корректив в динамику популяции (Шибаев, 2007; Бабаян; 2000).
Одной из основных эффективных мер регулирования промысла является установление общего допустимого улова (ОДУ). Общий допустимый улов – это биологически приемлемая для запаса величина годового вылова, соответствующая долговременной стратегии рационального промыслового использования данного запаса (Бабаян, 2000).
Процедура оценки ОДУ может включать в себя различные методы, объединенные рамками избранного алгоритма расчетов, конечным этапом которого является соотношение, связывающее искомую оценку с рекомендуемым значением интенсивности промысла и прогнозом величины промысловой части запаса.
Общий принцип оценки ОДУ крайне прост и заключается в следующем:
где t – индекс года промысла; I rec – рекомендуемое значение интенсивности промысла; FSB – биомасса промысловой части запаса. Интенсивность промысла может выражаться через коэффициент мгновенной промысловой смертности ( Frec ), эффективное промысловое усилие или коэффициент промысловой убыли (Бабаян, 2000).
Оценка ОДУ предусматривает решение двух самостоятельных задач: оценку биомассы запаса и обоснование величины управляющего воздействия на запас (Бабаян, 1982). Оценка, а точнее прогноз состояния запаса в год t получается на основе анализа особенностей динамики запаса по данным за имеющийся период наблюдений и экстраполяции выявленных тенденций на заданную перспективу. Оценка рекомендуемого промыслового воздействия на запас учитывает продукционные возможности запаса и напрямую связана с заранее принятыми целями и стратегией его эксплуатации. При этом рекомендуемая интенсивность промысла не должна превышать некоторый уровень, который отвечает продукционному потенциалу запаса и обычно ассоциируется с максимальным устойчивым выловом (MSY).
В общем случае продукционная способность зависит от множества факторов, учитывающих различные аспекты процессов размножения, весового роста и смертности, однако биомасса запаса в известном смысле является интегральным показателем ожидаемой суммарной результативности всех этих процессов. Таким образом, текущее состояние запаса (популяции) можно оценить путем сопоставления фактической величины биомассы с ее теоретическим значением, характеризующим некоторую важную особенность продукционного процесса в данном запасе. К таким значениям можно, в частности, отнести BMSY – величину биомассы, которая обеспечивает максимальную продуктивность запаса, или Bcrash – значение биомассы, ниже которго запас теряет устойчивость, то есть способность к расширенному воспроизводству. Часто состояние популяции связывают с нерестовой биомассой, учитывая ключевую роль, которую нерестовый потенциал играет в формировании динамики численности запаса.
Значения биомассы и других биологических показателей, прямо или косвенно характеризующих особые (или пороговые) состояния запаса, принято называть биологическими ориентирами (FAO, 1993). Биологические ориентиры выполняют роль точек отсчета при оценке текущего (или прогнозируемого) состояния запаса по известным фактическим значениям соответствующих параметров (например, биомассы нерестового запаса). Для выполнения этой функции биологические ориентиры должны соответствовать достаточно устойчивым продукционным характеристикам запаса.
промысловыми показателями (уловом, промысловым усилием и уловом на промысловое усилие), причем некоторые из них допускают двойную трактовку (табл. 2).
Различия в интерпретации некоторых параметров системы запас-промысел В настоящее время в основу современной теории рационального рыболовства положена концепция устойчивого развития и принцип предосторожности, выраженные четырьмя составляющими (FAO, 1993; 1995; UN, 1995; Бабаян, 2000):
1. Предупреждение или минимизация рисков, связанных с возможностью причинения ущерба эксплуатируемому запасу и соответствующей водной экосистеме.
2. Незамедлительное принятие заранее согласованных мер в случае реальной опасности для состояния эксплуатируемого запаса или среды его обитания.
3. Учет неопределенности (неполноты знаний о запасах промысловых видов и среде их обитания) в качестве объективно неизбежного фактора регулирования рыболовства.
4. Восстановление эксплуатируемых запасов, а также запасов ассоциированных или зависимых видов до высоких уровней продуктивности и поддержание их на этом уровне в течение всего периода эксплуатации.
В процессе адаптации теоретических принципов устойчивого развития и предосторожного подхода к регулированию конкретных промыслов, появилась необходимость расширить перечень традиционных ориентиров управления за счет опорных точек, учитывающих неопределенность в оценках используемых параметров.
Такие ориентиры помогают существенно уменьшить риск подрыва запаса за счет некоторого ограничения промысла.
Одной из главных особенностей предосторожного подхода является зональный принцип регулирования рыболовства, когда весь диапазон возможных состояний запаса (0, BK ) разбивается на отрезки, для каждого из которых устанавливается особый режим регулирования (рис. 9).
Рис. 9. Зональное представление области управления при предосторожном подходе: I – показатель интенсивности промысла; B – показатель состояния запаса; lim, buf и tr – индексы соответственно граничных, буферных и целевых ориентиров управления.
В общем случае при предосторожном подходе к регулированию рыболовства в зоне B Blim на промышленный лов вводится полный запрет, а допустимая область управления (промысловой эксплуатации) делится на 3 зоны:
– зону устойчивого состояния запаса без учета неопределенности (определяется граничными ориентирами I lim, Blim );
– зону устойчивого состояния запаса с учетом неопределенности (определяется буферными ориентирами I buf, Bbuf );
– зону стабилизации промысла на целевом уровне (определяется ориентирами управления Itr, Btr ).
Таким образом, в зависимости от уровней целевых и буферных ориентиров управления, зона стабилизации промысла может иметь самостоятельное значение при определении режимов регулирования или совпадать с зоной устойчивого состояния запаса с учетом неопределенностей.
Существует несколько версий интерпретаций рекомендаций по осуществлению предосторожного подхода, изложенных в Приложении II “Соглашения ООН по сохранению трансграничных запасов и запасов далеко мигрирующих видов рыб и управлению ими” (UN, 1995; Бабаян, 2000). В настоящей работе для оценки ориентиров по продукционной модели рассматривается интерпретация НАФО (Североатлантическая рыболовная организация), разработанная на исследовательской группе по граничным ориентирам в 2004 г. (Report…, 2004). В докладе группы есть рекомендации по использованию ориентиров, где в основе расчетов лежит продукционная модель Шефера.
По версии НАФО граничный ориентир управления по промысловой смертности Flim равен FMSY. Если этот показатель превышен, то считается, что запас находится в состоянии перелова. Blim – граничный ориентир по биомассе, ниже которого запас не должен опускаться. Считается, что при биомассе ниже уровня Blim величина пополнения, продуцируемая запасом, «низка» или неизвестна. При использовании продукционной модели Шефера рекомендуется Blim устанавливать на уровне 30% от BMSY. Буферные ориентиры Bbuf и Fbuf – это расчетные значения, при которых запас или его эксплуатация с высокой степенью вероятности не опускается до своих граничных ориентиров. Расстояние между граничными и буферными ориентирами выражает степень неопределенности оценки состояния запаса, то есть чем выше неопределенность, тем буферные ориентиры больше граничных. В промысловой биологии невозможность получения точных оценок, неполнота или заведомая недостоверность нашего знания об объекте исследования, получила название – неопределённость (Punt, Hilborn, 1997; Haddon, 2001). Она обусловлена рядом факторов.
Важнейшими источниками неопределённости являются нерепрезентативность имеющихся данных и изменчивость процессов, которые пытаются прогнозировать. Кроме того, существенное влияние на точность оценки имеет правильность аппроксимации процессов описываемых моделью, которая использована для оценки запаса. Хорошо известно, что оценки, получаемые по модели, могут существенно изменяться при добавлении новых данных (Improving…, 1998). Причины подобных изменений не всегда легко определить, но одной из основных является неправильная аппроксимация данных моделью и недостоверное, как правило, слишком упрощённое, описание процессов.
Уменьшение уровня неопределённости, через улучшение сбора данных и совершенствование самих моделей, имеет естественный предел. Поэтому одной из основных задач моделирования динамики численности является оценка уровня неопределённости, доверительного интервала получаемых оценок, или уровня риска, а также их учёт при выработке рекомендаций по эксплуатации промыслового запаса.
Под уровнем риска здесь понимается вероятность какого-либо нежелательного явления. Например, вероятность того, что биомасса нерестового запаса или уровень промысловой смертности достигнут некоторой критической величины или определенного ориентира управления. Частично эту проблему решает предосторожный подход, при котором уровень неопределённости учитывается при определении значений буферных ориентиров. Кроме того, в последнее время, при оценке запаса и прогнозировании его динамики, всё чаще используются вероятностные модели, в которых оценки и прогнозные величины представлены не точечными величинами, а в виде диапазона их возможных значений. Современные модели оценки запаса включают в свой состав возможность описания некоторых элементов имеющейся неопределённости. При использовании байесовских стохастических подходов в выработке управленческих решений буферные ориентиры могут быть опущены, так как расчет граничных ориентиров связан так же с расчетом их неопределенности и рисков.
На протяжении 10 лет ориентирами управления при промысле камчатского краба в заливах полуострова Аляска служат численность минимального зарегистрированного нерестового запаса, а также принимается в расчет средний вылов (mean yield), его стабильность, промысловые возможности (harvest opportunity), стабильность нерестового запаса. Пороговый уровень управленческой стратегии при промысле камчатского краба в Бристольском заливе служит для защиты репродуктивного потенциала запаса и достижения максимального долговременного ежегодного вылова (Zheng et. al., 1997). При выборе этого ориентира руководствуются соображением, что есть некий уровень продукционного потенциала, ниже которого риск перелова повышается при сохранении усилий.
В последние годы в ПИНРО активно ведутся поиски приемлемых биологических ориентиров управления баренцевоморским запасом камчатского краба (Баканев, 2006).
Выбор таких ориентиров также должен быть связан с концепцией прибавочной продукции и определением продукционной кривой. Однако в настоящее время нет данных о величине продукционного потенциала баренцевоморской популяции, а также о его теоретически возможной максимальной величине ( BK ) после завершения натурализации краба в Баренцевом море. Скорее всего, расчеты ориентиров и рисков их превышения должны быть сопряжены с некими условными уровнями численности, ниже которых падение численности будет нежелательным. Такие пороговые уровни численности могут быть оценены без учета допущения о равновесном состоянии запаса и, возможно, опираться как на оценку продукционных показателей, так и на опыт эксплуатации запасов в его нативном ареале.
Материалом для диссертационной работы послужили результаты количественного учета камчатского краба в Баренцевом море, выполненного в 1994–2008 гг. в ходе осенних траловых съемок (табл. 3), а также 21 научно-промысловых рейсов. Автор участвовал в 12 рейсах в качестве инженера-гидробиолога, научного сотрудника и начальника экспедиций. Основой для оценки популяционной численности послужили уловы 1481 траловых станций исследовательских съемок, в которых 37345 особей было взято на биологический анализ. Из научно-промысловых рейсов было проанализировано свыше 60 тыс. особей камчатского краба, из них биологическому анализу подверглось тыс. экземпляров. Автором было промерено около 7000 особей, из них на биологический анализ взято 5500 животных.
Характеристика использованного материала, собранного во время Стратифицированные траловые съемки выполнялись в конце лета – начале осени, обычно с 15 августа по 15 сентября. Период съемки выбран таким образом, чтобы зафиксировать максимальное количество крабов, которые в это время распределяются в районах нагула. С расширением баренцевоморского ареала камчатского краба сроки стратифицированной съемки постепенно увеличивались. Принималось, что съемки охватывали большую часть ареала вида в баренцевоморском регионе. В случае возрастания плотности скоплений в районах его учета потребовались дополнительные затраты рабочего времени, так как при повышении плотности отдельных скоплений объекта точность его оценки в страте может уменьшаться и поэтому в таких районах выполнялись дополнительные траловые станции. Съемка проводилась в пределах следующих промысловых районов Баренцева моря (рис. 10): Рыбачья (11), Кильдинская (12) и Канинская (3) банки, Мурманское мелководье (7), Западный и Восточный Прибрежные районы (13, 14).
Рис. 10. Карта-схема районов траловой съемки камчатского краба в РЭЗ Баренцева моря (мелкими цифрами обозначены номера страт, крупными – номера стандартных промысловых районов).
В 1994 г. акватория районов была разделена на страты, границы которых соответствуют изобатам 100, 200 и 300 м (Герасимова, Кузьмин, 1997). Такое разделение связано с фактическим распределением камчатского краба на момент начала учетных траловых съемок. Общее количество станций при проведении съемки распределялось по районам и стратам пропорционально их площади с учетом накопленных знаний о вариациях распределения краба в каждом районе и страте. Расчеты осуществлялись согласно методике разработанной для стратифицированных съемок гидробионтов (Cochran, 1963; Баканев, Пинчуков; 2006). Количество площадок Nk в страте k, равных Среднеарифметический индекс численности краба на обследованной акватории рассчитывается по формуле:
где nk – количество тралений в страте k, yi,k – улов краба (экземпляры) за i-траление в страте k. Дисперсия s y среднеарифметического индекса численности в страте k определяется по формуле:
По району m или по всей площади съемки в целом рассчитывается средневзвешенное значение индекса численности ym :
где Nm – количество площадок в районе m, равных площади одного стандартного численности в районе m определяется по следующей формуле:
Ошибка среднего по формуле:
Коэффициент вариации по формуле:
Индекс численности Ym краба по району m:
где t (P ) – значение t-критерия, зависящее от выбора доверительной вероятности Рс и числа степеней свободы.
Выбор величины Рс – скорее результат договоренности, чем строгого анализа. При оценке ориентиров управления Международный совет по изучению морей (ICES), например, рекомендует Рс=95%, тогда как другие, не менее авторитетные источники, считают возможным уменьшить доверительную вероятность до Рс 50% (Anon, 1992, Anon, 2000). В нашем случае мы выбрали Рс=95%, тогда:
Согласно методике рандомизированно-стратифицированной съемки (Doubleday, 1981) для получения более точных результатов количество тралений в стратах с максимальными уловами могло быть увеличено при наличии резерва рабочего времени.
доверительный интервал находился в границах ±25-50% при доверительной вероятности 95% (ICNAF Redbook, 1978).
В качестве учетного орудия лова использовался донный трал (тип ДТ-1625) с мелкоячейной вставкой в кутке (внутренний размер ячеи – 20 мм). Продолжительность учетных тралений составляла 1 час, скорость тралений – 3,0–3,2 узла. Маршрут съемки планировали таким образом, чтобы в кратчайшие сроки и максимально полно обследовать районы распределения краба. Съемки проводились круглосуточно. Переходы между станциями совершались по кратчайшему расстоянию.
Сбор и обработка биологического материала в съемках выполнялась в соответствии с методиками, принятыми в ПИНРО (Инструкции и методические рекомендации…, 2001). Расчеты индексов численности по стратам и районам производились раздельно по полу, стадиям линьки и размерным группам с интервалом в 10 мм, а также по укрупненным размерным группам: промысловые (при длине карапакса (ДК) больше 132 мм) и непромысловые самцы (ДК < 132 мм), икроносные самки и самки без икры. Для расчетов по модели CSA в качестве исходных данных брали количество половозрелых самцов, пойманных за период съемки. Половозрелых самцов разделяли на четыре группы: пререкруты-2, пререкруты (PR), рекруты (RE) и пострекруты (PO) (Баканев, 2003). Средний годовой прирост половозрелых самцов по длине карапакса в Баренцевом море составляет 17 мм (Nilssen, Sundet, 2006). В соответствии с годовым приростом пререкрутами-2 было условлено считать самцов с длиной карапакса 96– 113 мм, пререкрутами – 114–131 мм, рекрутами – самцов на 1–2-й стадиях линьки с размерами 132–149 мм, пострекрутами – самцов на 1–2-й стадиях линьки с ДК более 150 мм и самцов на 3–4-й стадиях линьки с размерами 132 мм и более. Расчеты осуществлялись согласно стандартной методике, принятой в ПИНРО и разработанной W.
G. Cochran (1963).
При отборе проб на плодовитость кладку икры вместе с плеоподами помещали в отдельную емкость и фиксировали 4 %-ным раствором формальдегида. При оценке плодовитости кладку икры отделяли от плеопод и обсушивали на фильтровальной бумаге до постоянной массы. Кладку взвешивали, просчитывали количество икры в двух навесках по 500 мг. Определяли массу одной икринки и общее количество икры в кладке по соотношению массы кладки и средней массы икринки.
Стадии полового созревания определялись путем построения логистической кривой по данным о доле икроносных самок в размерных группах с шагом 10 мм по длине карапакса. Кривые были построены методом взвешенной нелинейной регрессии (Otto, 1986). Сравнительный анализ скорости полового созревания проводили путем оценки и сравнения размеров, при которых 50 % самок данной популяции камчатского краба имеют наружную икру. Значимость наблюдаемых различий оценивали по критерию Стьюдента.
За весь период исследований было обработано 462 пробы на плодовитость.
При обработке исходной биологической информации и расчетов индексов численности используются пакеты прикладных программ BIOFOX FoxPro, Microsoft Excel, Microsoft Access, а также макросы к этим программам, написанные старшим научным сотрудником лаборатории промысловых беспозвоночных ПИНРО Любиным П.
Использовались данные вылова (масса и численность улова с учетом пола), размерного состава уловов и среднего улова на единицу усилия промыслового флота.
Объём вылова и производительность лова вычислялись по ежесуточным донесениям, поступающим со всех промысловых судов. Стандартизированные индексы численности (улов на единицу промыслового усилия) вычислялись для каждого типа ловушек.
Информация о половом и размерно-возрастном составе уловов собиралась автором и наблюдателями на борту промысловых судов (Инструкции и методические рекомендации по сбору..., 2001). В процессе подготовки к оценке запаса массовые промеры объединялись по районам.
Завершающим этапом подготовки массивов данных о промысле являлась оценка общего фактического вылова. В нашей работе использовались два ряда данных, характеризующие ежегодную интенсивность промыслового изъятия: официальный вылов и общий вылов с учетом нелегального. Нелегальный вылов учитывал браконьерские выгрузки на внутренний рынок и согласно литературным источникам составлял не менее 0,40 млн. экз. в 2003–2004 гг. (Соколов, 2005). Также была оценена величина экспорта краба в США российским флотом в 2005–2006 гг. по данным статистического отдела (http://www.st.nmfs.gov/). Средняя величина экспорта краба Россией в США в 1997– гг. составила около 10 тыс. тонн в год. В 2005 г. экспорт увеличился на 70% и превысил среднемноголетнюю величину на 7 тыс. тонн, что с учетом среднего веса промыслового краба (4,5 кг) и переходного технологического коэффициента (65%, Степаненко, Двинин, 2006) составляет 2,6 млн. крабов. Учитывая современное состояние ресурсов камчатского краба в дальневосточных морях (Долженков, Болдырев, 2006) можно допустить с достаточной степенью уверенности, что скачок в российском экспорте камчатского краба в США обусловлен поставками баренцевоморского краба с открытием его коммерческого промысла.
Для оценки запаса, динамики численности, расчета ОДУ и ориентиров управления использовалась продукционная модель на основе уравнения Шефера (Schaefer, 1954), а также две модели, основанные на динамике численности размерных групп: LBA, состоящую из 12 размерных групп (Zheng et al., 1995) и CSA, включающую 3 размерновозрастных группы (Zheng et al., 1997).
В качестве продукционной модели использовалась стандартная логистическая кривая, описываемая уравнениями [15, 20], часто называемая «продукционной моделью промысловую смертность, где скорость экспоненциального роста численности популяции (r) выражается через максимально устойчивый вылов (MSY), имеет вид (Richards, 1959;
Pella, Tomplinson, 1969):
где Bt – численность в году t, MSY – максимальный устойчивый вылов, K – максимальная возможная численность при отсутствии промысла (емкость среды), Ct – вылов.
Оценка абсолютной численности в большинстве моделей имеет большую неопределенность, если нет информации о соотношении индексов и абсолютной численности запаса. С целью уменьшения неопределенности в оценке параметра улавливаемости для управленческих целей обычно пользуются относительными величинами (Prager, 1994, 2002; Cardin, 2000; Improving…, 1998; Zheng et. al., 1997).
Относительная оценка становится возможна при делении уравнения [33] на BMSY, то есть на численность, соответствующую максимальной прибавочной продукции. С добавлением вычислительной ошибки, уравнение принимает следующий вид:
где Pt – уровень численности в году t относительно BMSY (Pt = Bt / BMSY).
Вычислительная ошибка оценивания моделируемой численности (v) имеет нормальное Ненаблюдаемая переменная Pt может быть выражена через наблюдаемый показатель относительной численности (survt), то есть индекса, рассчитанного по результатам исследовательской съемки:
Отношение индекса к реальной величине биомассы выражается через коэффициент улавливаемости qs, a ek – ошибка измерения индекса численности, имеющее нормальное Исходя из функциональной связи, выраженной уравнениями [34, 35], величина вылова в расчетах влияет на оценку численности и коэффициентов улавливаемости, то характеризующих ежегодную интенсивность промыслового изъятия: официальный вылов Coff, общий вылов с учетом нелегального Ctot (табл. 4). Требуется определить следующие параметры: m, MSY, BMSY, коэффициент улавливаемости, ошибки наблюдений и численность краба для первого года исследований.
Официальный и общий выловы, индекс промысловых самцов камчатского краба (тыс. экз.) в Баренцевом море в 1994–2006 гг., использующиеся в продукционной модели В качестве модели, описывающей динамику численности размерных групп, использовалась LBA (length-based analysis; Zheng et al., 1995), в основу, которой положены уравнения [4–12] (раздел 1.3). Моделирование численности самцов камчатского краба с новым карапаксом и со старым происходит раздельно, при этом допускается, что естественная смертность у крабов всех размерных групп одинаковая. Сокращенный вид модели можно выразить через два основных уравнения:
где Nl,t – численность крабов с новым карапаксом в размерной группе l и году t, Ol, y – численность крабов со старым карапаксом в размерной группе l и году t, M – коэффициент естественной смертности, ml, y – вероятность линьки крабов размерной группы l в году t, Rl, y – пополнение размерной группы l в году t, y – времення задержка от момента съемки до средней даты периода промысла, P,l – пропорция линяющих крабов группы l после одной линьки, достигающих группы l, Cl,t – вылов крабов группы l в году t, nl и ol – ошибки наблюдений, априорные распределения которых равны обратному гамма-распределению с модами коэффициентов вариации индексов Nl,t и Ol,t, полученных по съемкам. Учитывая, что съемки выполняются в конце года, допущение о межлиночном периоде крабов делалось на основании исследований российских и американских ученых (Левин, 2001; Zheng et al., 1995). Предполагалось, что крабы с экзоскелетом на стадии линьки 1 и 2 линяли в течение текущего года (крабы с новым экзоскелетом), а крабы с экзоскелетом на стадии линьки 3 и 4 – в течение предыдущих лет (крабы со старым экзоскелетом).
Для оценки абсолютной численности камчатского краба использованы серии индексов размерных групп, рассчитанные по результатам траловых съемок в 1994– гг. Были определены 12 размерных групп с интервалом 10 мм по длине карапакса (indCL90 – indCL200, табл. 5). Из-за низкой облавливаемости молодых особей для первой размерной группы был установлен минимальный размер 90 мм по длине карапакса. В последнюю группу входят самцы с размерами равными и более 200 мм.
Индексы численности камчатского краба (тыс. экз.), использующиеся в модели Коэффициент улавливаемости для разных размерных групп функционально зависит от размеров крабов, входящих в эти группы (Zheng et al., 1995, Баканев, 2003).
Отношение индексов к абсолютной численности групп обозначалось уравнением улавливаемости для каждой из 12 размерной группы:
где a и b – коэффициенты, CL90t,…, CL200t - численности групп по модели, рассчитываются как:
Принимается, что ошибки наблюдений CL90t : CL200t распределяются логнормально и распределение в уравнениях наблюдений данных выглядит так:
Данные по общему вылову были разделены на 8 размерных групп (табл. 6).
улавливаемости, ошибок наблюдений и численности крабов для первого года исследований по размерным группам. С тем чтобы избежать излишней параметризации модели оценка вероятности линьки ml,t и перехода из одной размерной группы в другую была получена с помощью модуля «Поиск решения», программы MS Excel с использованием метода поиска Ньютона.
Общий вылов камчатского краба по размерным группам (Сl,t, тыс. экз.) в Баренцевом море На основе модели LBA в 90-е годы прошлого столетия была разработана модель CSA (catch survey analysis). CSA является частным случаем LBA, когда численность популяции мала для получения качественных данных по размерному составу (см. раздел 1.4). Отличие модели CSA от LBA состоит лишь в том, что в ней используются более укрупненные размерные группы, а в расчетах приемлемо использовать данные по размерному составу лишь половозрелых самцов. Согласно уравнениям [13, 14], для оценки численности ненаблюдаемых переменных пререкрутов ( PRt ), рекрутов ( REt ) и пострекрутов ( POt ) использованы серии индексов соответствующих размерных групп, полученные в ходе траловых съемок в 1994–2006 гг. (survPRt, survREt, survPOt, табл. 7).
Индексы численности размерных групп камчатского краба (тыс. экз.) в Баренцевом дифференцированными коэффициентами улавливаемости для каждой из трех групп: qpr, qre, qpo. Принимается, что ошибки наблюдений pr, re, po распределяются логнормально и survPR ~ q pr PRt e pr, survRE ~ qreREt ere, survPO ~ q poPOt e po. В модели использовались два ряда данных, характеризующие ежегодную интенсивность промыслового изъятия:
официальный вылов Coff, общий вылов с учетом нелегального Ctot (табл. 4). Требуется определить следующие параметры: вероятность линьки для пререкрутов ( mPR ), параметры роста ( GPR,RE, GPR,PO ), коэффициент естественной смертности (M), коэффициенты улавливаемости, ошибки наблюдений и численность трех размерных групп краба для первого года исследований.
Для оценки параметров моделей и расчетных переменных используется байесовский подход. Он базируется на следующих четырех компонентах. Априорная (предварительная) информация (PRIORS) и данные (DATA), полученные в ходе эксперимента, объединяются с помощью некоего метода (method) для получения апостериорного (POSTERIORS) знания. Эти компоненты могут быть представлены в виде:
PRIORS DATA POSTERIORS
Данная формула – запись уравнения Байеса, принятая в статистике (см. раздел 1.3;уравнения [17–19]). Входной компонент DATA представляет собой индексы численности краба, оцененные по съемкам, и вылов. В качестве компоненты PRIORS берется предварительное (приблизительное) знание о параметрах модели, количественно представленное в виде распределений вероятностей их возможных значений. В качестве компоненты POSTERIORS выступает существенно уточненное знание о параметрах модели, также представленное в виде распределений вероятностей значений. Метод расчета POSTERIORS зависит от сложности моделируемых процессов. Для описания биологических процессов с высокой степенью неопределенности использовался алгоритм вычисления, основанный на построении цепи Маркова (рис. 11). Для этого применяется итерационная процедура, на каждом шаге которой рассчитываются модельные значения параметров и переменных. Каждый итерационный шаг включает в себя три этапа.
На первом этапе первой итерации происходит расчет численности животных (data) по биологической модели (BIO) с использованием стартовых (начальных) значений параметров и переменных (INITS). Биологическая модель в нашем случае представляет собой одну из трех, описанных ранее: продукционную, LBA или CSA. Стартовые значения генерируются встроенным модулем программы, а затем, если необходимо, корректируются с учетом биологически правдоподобных значений. На первом этапе второй и последующих итераций берутся величины параметров и переменных, рассчитанные в ходе предыдущей итерации.
На втором этапе модельные и эмпирические (DATA) значения численности включаются в функцию правдоподобия (LIKELIHOOD, уравнение [22]) и вычисляется (like) мера их согласованности, то есть определяется вероятность того, что эмпирические значения могут быть получены с использованием данной биологической модели с заданными стартовыми значениями параметров и переменных.
На третьем этапе по формуле Байеса (BAYES, уравнение [19]) количественно оценивается способность модели с заданными параметрами (PRIORS) генерировать эмпирические значения DATA. В процессе расчетов происходит настройка параметров (PRIORS) в условиях, характеризуемых данными наблюдения (DATA), и рассчитываются новые апостериорные значения модельных параметров и переменных. Оцененные показатели параметров и рассчитанные модельные переменные, характеризующие модель в данный момент, играют роль начальных условий для следующей итерации.
Данная процедура расчета является заключительным этапом одного итерационного шага или звена в цепи Маркова (MCMC) и реализована с использованием метода Гиббса.
В процессе итераций генерируются (stat) выборки распределений значений параметров и расчетных переменных. Статистический анализ таких выборок, в который входит расчет средних, ошибок средних, стандартных отклонений, медиан и доверительных интервалов параметров, является результатом прогона модели. Окончание итерационного процесса (прогона) происходит, когда ошибка среднего искомого параметра, рассчитанная в цепи Маркова, составит менее 5% его стандартного отклонения (Spiegelhalter et al., 2000).
BAYES LIKELIHOOD BIO
BAYES LIKELIHOOD BIO
POSTERIORS
Рис. 11. Принципиальная блок-схема оценки параметров биологических моделей (показаны два итерационных шага).POSTERIORS – апостериорные распределения параметров и BIO – биологическая модель (уравнение популяционной динамики);
data – расчет численности животных по биологической модели;
stat – расчет апостериорных распределений параметров и После итерационного прогона суммарная выходная статистика подвергалась анализу. Для оценки качества настройки модели выполнялось сравнение фактических значений и их апостериорных (рассчитанных в модели) распределений с использованием двух критериев согласия:
Во-первых, рассчитывались остатки, т.е. разности между логарифмами фактических индексов численности по съемке и рассчитанных индексов по модели. Чем лучше модель согласуется с данными, тем меньше величина остатков. Суммарная статистика отклонений (RSS, от residual sum of squares) показала степень смещения модельных значений от наблюдаемых всего временного ряда: RSS (ln Nt ln Nt )2, где Nt – наблюдаемый индекс численности по съемки, Nt – расчетный индекс численности по модели.
Во-вторых, адекватность модели анализировалась путем сравнения распределения каждой фактической величины со своим апостериорным прогнозным распределением (Gelman et al., 1992, 1995). При этом после настройки параметров в модели имитируется набор наблюдаемых данных. Вероятности имитируемых значений фактических величин datarep имеют распределение:
В выражении [42] показатель P(datarep | ) представляет собой выборку наблюдений из распределений, рассчитанных в модели – P( | dataobs).
Данный критерий проверяет, нет ли смещения наблюдённого индекса численности относительно медианы распределения имитируемых на каждой итерации модели индексов численности. Степень отклонения суммируется в векторе p-значений, рассчитываемых как пропорция N-итераций, в которых выборка имитируемых данных наблюдения из апостериорного распределения превышает истинные значения данных наблюдения:
где I(x) = 1, если x истинно, и I(x) = 0 если, x ложно. Величины близкие к 1 или 0 в векторе p-значений показывают, что набор фактических значений в значительной мере отклоняется от медианы апостериорного распределения. Величина p.value близкая к 0, показывает, что расчетная величина datarep приблизительно в 50% случаев больше входного значения и в 50% – меньше, то есть отклонения не имеют систематического смещения.
За этапом оценки качества настройки параметров производится исследование моделируемой популяционной динамики, выясняется степень адекватности полученной модели. При этом оценивается поведение популяции при различных исходных предпосылках, например, меняются показатели степени эксплуатации, величины пополнения, уровня естественной смертности, продукционной способности популяции.