CC BY
91
УДК 57.087.1
ОБ ОДНОМ МЕТОДЕ ОЦЕНКИ ЗАПАСОВ И ПРОГНОЗА ОДУ МОРСКИХ ПРОМЫСЛОВЫХ РЫБ НА ОСНОВЕ НЕПРЕРЫВНОЙ МОДЕЛИ ДИНАМИКИ ВОЗРАСТНОЙ СТРУКТУРЫ ПОПУЛЯЦИИ
О. И. Ильин
Н. с., Камчатский научно-исследовательский институт рыбного хозяйства и океанографии 683000 Петропавловск-Камчатский, Набережная, 18 Тел., факс: (415-2) 41-27-01; 42-57-96 E-mail: ilin@kamniro-avacha. kamchatka. ru
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ, ОЦЕНКА ЗАПАСОВ, ОПТИМАЛЬНЫЙ РЕЖИМ ИЗЪЯТИЯ
Рассмотрен один метод оценки запасов и прогноза ОДУ морских промысловых рыб на основе непрерывной модели динамики возрастной структуры биологической популяции. Приведены результаты численного моделирования ретроспективной динамики запаса и оптимального режима изъятия восточноохотоморского минтая (Therarga Chalcogramma).
ON THE ONE METHOD FOR STOCK ABUNDANCE ESTIMATION AND FORECAST OF TAC OF MARINE COMMERCIAL FISH SPECIES BASED ON THE MODEL OF AGE-STRUCTURED POPULATION HAS ANALYZED О. I. Ilyin
Scientist, Kamchatka Research Institute of Fisheries and Oceanography 683000 Petropavlovsk-Kamchatsky, Naberejnaya, 18 Tel., fax: (415-2) 41-27-01; 42-57-96 E-mail: ilin@kamniro-avacha. kamchatka. ru
MATHEMATICAL SIMULATION, STOCK ABUNDANCE ASSESSMENT, OPTIMAL REMOVAL REGIME
One method for stock abundance estimation and forecast of TAC of marine commercial fish species based on the model of age-structured population has analyzed. The results of the numerical modelling of the age structure retrospective dynamics and the optimal regime of exploitation of the walleye pollock (Therarga Chalcogramma) stock in the east part of the Okhotsk Sea are provided.
Сегодня существует огромное количество математических моделей и методов теории рыболовства для оценки состояния запасов и прогнозирования величины изъятия морских промысловых объектов. При этом особую популярность получают те методы, для применения которых требуется наиболее доступный набор данных — возрастной состав уловов, а также масса, мгновенные коэффициенты естественной смертности и доля половозрелых рыб в зависимости от возраста. К числу таких моделей можно отнести и непрерывную модель биологической популяции с возрастной структурой. Под непрерывной понимается модель с непрерывными временем и возрастом.
Первая непрерывная модель динамики возрастной структуры биологической популяции построена, по-видимому, в работе фон Ферстера (Von Foerster, 1959). В дальнейшем эта модель с ее различными модификациями и обобщениями исследовалась различными авторами (Sinko et al., 1967; Динамическая теория ..., 1974; Абакумов и др., 1994). Учет возрастной (размерной) структуры в непрерывном случае приводит или к дифференци-
альным уравнениям в частных производных типа уравнения неразрывности (Von Foerster, 1959), или к эквивалентным им интегральным уравнениям (Перцев, 1999).
В работе А.И. Абакумова с соавторами (1994) подробно проанализирован один метод оценки состояния популяции (запаса) с учетом возрастной структуры, построена специальная численная схема расчетов по этому методу, учитывающая специфику рыбопромысловой информации. Эту схему можно рассматривать как некоторое обобщение процедуры виртуально-популяционного анализа (Gulland, 1969). Расширению методов оценки запаса на общий случай учета популяционных структур посвящена работа (Абакумов, 1991). Однако этот метод, как и ВПА, имеет принципиальный недостаток — сильную чувствительность к оценкам естественной смертности и ошибкам в значениях уловов по возрастным группам. В настоящей работе используется несколько другой подход, близкий по смыслу к статистическим методам и настроечным процедурам с привлечением дополнительной информации из теории рыболовства
(Deriso et al., 1985; Richard Methot, 1990; Kizner et al., 1997), основным преимуществом которых является возможность привлечения в качестве дополнительной информации данных по промысловому усилию, данных съемок и индексов численности с возрастной структурой. Ретроспективные оценки популяционных параметров (промысловой смертности, пополнения) в указанных процедурах находятся из условия наилучшего приближения модельного запаса к данным прямых наблюдений методом наименьших квадратов или максимального правдоподобия.
В настоящей работе наилучшее приближение обеспечивают оценки, доставляющие минимум некоторому квадратичному функционалу. Оптимальный режим изъятия определяется из решения задачи оптимальной эксплуатации (оптимального сбора урожая) популяции с функциональными ограничениями. Задача оптимальной эксплуатации для непрерывной модели с учетом внутрипопуля-ционных структур в стационарном состоянии рассматривалась различными авторами (Свирежев и др., 1980; Абакумов, 1991; Ильин, 2007) и в динамике — в работах (Ильин, 2006, 2007). Условие сохранения популяции на достаточно высоком уровне продуктивности и поддержания ее на этом уровне в течение всего периода эксплуатации приводит к задаче с функциональными ограничениями (Ильин, 2007). Ограничения на величину нерестового запаса популяции имеют тот же смысл, что и граничные ориентиры по биомассе в предо-сторожном подходе к управлению промыслом у
В.К. Бабаяна (2000).
МАТЕРИАЛ И МЕТОДИКА
Уравнение модели динамики возрастного состава биологической популяции с непрерывным временем и непрерывным возрастом имеет вид (см., например, Динамическая теория..., 1974):
Здесь t — время, т — возраст объектов, © — максимальный возраст, x(t,x) — плотность распределения численности популяции; u(t,x) — функция, характеризующая интенсивность промысла данной популяции, m(x) — коэффициент естественной смертности объектов.
Это уравнение необходимо дополнить краевыми и граничными условиями. В зависимости от целей моделирования и имеющейся в распоряжении исследователя информации об объек-
те, эти условия могут быть заданы в следующих комбинациях.
1) Краевое условие на левой стороне области G = {(7, т):0<КХ, 0<т<©} в виде x(t, 0) = F(t, s), где s — нерестовый запас моделируемой популяции, и начальное условие вида x(0, т) = х0(т). Некоторые задачи с такими условиями рассмотрены в работах (Von Foerster, 1959; Динамическая теория..., 1974; Ильин, 2007). Здесь краевое условие определяет интенсивность появления новых объектов в популяции в момент времени t, F(s) описывает процессы плотностной саморегуляции на ранних стадиях жизненного цикла (Бивертон и др., 1969; Рикер, 1979; Максименко и др., 2003). Нерестовый запас s определяется по формуле:
где a(r) — доля половозрелых объектов, w(t) — масса объектов. Начальное условие задает «стартовое» возрастное распределение численности. В случае F(t, s) = R(t), где R(t) — заданная интенсивность появления новых объектов в популяции, решение может быть найдено аналитически (например, Абакумов, 1991).
2) Краевое условие на левой стороне области G в виде x(t, 0) = F(s) + краевое условие на правой стороне в виде x(t, ©) = xx(t). Последнее краевое условие определяет плотность объектов «терминального» возраста. Примеры задач с такими условиями приведены в работах (Абакумов, 1991; Ильин, 2007). В качестве плотности объектов x{t, ©) «терминального» возраста © можно выбрать сумму выловленных объектов с возрастом, превосходящим ©.
3) Начальное условие x(T, т) = х2(т), с помощью которого задаётся возрастное распределение объектов в «терминальный» год + краевое условие для плотности объектов «терминального» возраста x(t, ©) = xx{t). В этом случае решение уравнения может быть найдено аналитически и не представляет особых проблем (см., например, Абакумов, 1991). Начальное распределение объектов может быть задано по результатам съемок.
Если необходимо получить ретроспективные оценки состояния популяции (запаса), используют условия 2 и 3. В случае, когда поставлена задача прогноза состояния запаса и определения оптимальной стратегии изъятия, задается условие 1. В настоящей работе используются первые два типа начальных и краевых условий.
Решение уравнений модели динамики возрастного состава биологической популяции с непрерыв-
ным временем и непрерывным возрастом осуществлялось численными методами на компьютере. Написана программа, реализующая алгоритм численного решения (Ильин, 2007).
Оценка запасов и прогноз ОДУ морских промысловых рыб на основе непрерывной модели динамики возрастной структуры биологической популяции осуществляется в два этапа. На начальном этапе находятся ретроспективные оценки состояния запаса из условия наилучшего приближения модельных оценок к данным по возрастному составу уловов и съемкам. Пусть имеются данные о возрастном составе уловов — плотность выловленных объектов С(£ т) и данные т съемок 8(£), 1 <\<т в моменты времени Ьх. Формулируется следующая оптимизационная задача (см., например, Ильин, 2007): на решениях системы уравнений динамики возрастного состава популяции
(1)
требуется минимизировать квадратичный функционал
(2)
при ограничениях на коэффициенты Рикера и интенсивность промысла
0 <а, 0<р, 0 <и<и (3),
где и — некоторое число. В настоящей работе для запаса восточноохотоморского минтая коэффициент интенсивности промысла в ретроспективе и($, т) ограничивается сверху величиной и = 1/год.
Первое слагаемое в (2) характеризует сумму квадратов отклонений модельных значений нерестового запаса от наблюдаемых 8(А) в моменты времени и, второе слагаемое представляет собой среднеквадратическое уклонение модельной плотности выловленных объектов и(/; т) х х(/; т) от на-
блюдаемой C (t, т). Здесь прослеживается связь задачи (1)-(2)-(3) с такими общеизвестными моделями, как «CAGEAN» (Deriso R.B. et al., 1985), «ISVPA» (Kizner et al., 1997) и «Synthesis» (Richard, Methot, 1990). В этих моделях минимизируется взвешенная сумма суммы квадратов отклонений логарифмов наблюдаемых уловов от логарифмов «теоретических» уловов и суммы квадратов отклонений логарифмов модельной и наблюдаемой нерестовой биомассы.
Задача оптимизации (1)-(2)-(3) может быть решена только численными методами. Автором настоящей работы построен алгоритм (Ильин, 2007) приближенного решения указанной задачи методом условного градиента (Васильев, 1980) и создана компьютерная программа, осуществляющая поиск решения по указанному алгоритму.
После того, как решена задача (1)-(2)-(3) и получены ретроспективные оценки плотности численности x(t, т), интенсивности промысла u(t, т) и коэффициентов Рикера в связи «запас-пополнение», формулируется задача определения оптимального режима изъятия (задача оптимальной эксплуатации). В качестве «стартового» распределения численности на прогноз выбирается «терминальное» возрастное распределение численности из ретроспективной оценки. Интенсивность промысла представляется в виде произведения промыслового усилия p(t) на коэффициент улавливаемости д(т) (долю изымаемой биомассы на единицу промыслового усилия). Из ретроспективной оценки u(t, т) определяется коэффициент улавливаемости д(т) по формуле (см., например, Ильин, 2007):
1
-idt, pit) Ф 0
(4),
тЪ р(п
где/?(/), 0<Х<Т— известное промысловое усилие в ретроспективный период. Интенсивность рождения новых объектов в популяции в прогнозные годы предполагается равным некоторому постоянному значению Я (например, равное среднему арифметическому от плотности пополнения за весь ретроспективный период). Из каких-либо соображений определяется критический уровень нерестовой биомассы 8, ниже которого не должна понижаться модельная биомасса нерестового запаса на начало каждого календарного прогнозного года. Критический уровень 8 соответствует граничному ориентиру по биомассе 88В11т в пре-досторожном подходе к регулированию промысла у В.К. Бабаяна (Бабаян, 2000). Задача оптимальной эксплуатации состоит в определении такого оптимального промыслового усилия р(/% ко-
торое на решениях уравнения модели возрастного состава популяции
(5)
максимизирует функционал
F=\ jP(t)'q(T)'x(t>T)dtdt —>sup (5^
выражающий суммарный вылов рыб с возрастом от 0 до 0 за период времени от 0 до Т\ При этом ограничения на искомое промысловое усилие и биомассу нерестового имеют вид:
^ ^ 0 < p(t)<P (7),
В ограничении (7) величина P представляет максимально возможную интенсивность промыслового усилия. Другими словами, P — это максимальное мгновенное промысловое усилие, которое может создать промысловая флотилия.
Необходимо отметить, что при определении оптимальной стратегии изъятия здесь не используется концепция максимального уравновешенного улова (Бабаян, 2000). Можно сказать, что оптимальная величина изъятия определяется из возможностей промыслового флота и из условия сохранения нерестовой биомассы популяции не ниже некоторого уровня.
Для приближенного решения задачи оптимального управления (5)-(6)-(7)-(8) в работе (Ильин, 2007) используется численный метод условного градиента. Алгоритм решения реализован в компьютерной программе.
Набор данных, необходимых для расчетов, составляют:
1) коэффициенты естественной смертности по возрастным группам;
2) масса объектов по возрастным группам;
3) доля половозрелых объектов по возрастным группам;
4) данные о количестве промысловых усилий;
5) данные съемок о величине нерестовой (или общей) биомассы запаса.
Такой же набор входных данных имеют такие модели, как «CAGEAN» и «Synthesis». Последняя модель базируется на уравнении Баранова (Баранов, 1918) и относится к классу стохастических се-
парабельных когортных методов. Автором настоящей работы создана компьютерная программа, в которой реализован алгоритм североамериканской «Synthesis». Она успешно использовалась для расчетов запасов западноберинговоморского минтая, корфо-карагинской сельди, трески и наваги юго-западной части Берингова моря.
Рассмотрен пример использования указанной схемы для оценки состояния запасов и прогнозирования оптимальной величины изъятия восточноохотоморского минтая. Ретроспективные оценки, полученные по непрерывной модели динамики возрастного состава биологической популяции, сравнивались с ретроспективными оценками по модели «Synthesis».
Все необходимые данные для расчетов предоставлены лабораторией морских промысловых рыб КамчатНИРО. В расчетах использовались данные об уловах восточноохотоморского минтая по возрастным группам 0+ - 11+ за период с 1964 по 2008 годы. Фактические значения среднемноголетней массы (кг), доли половозрелых рыб и мгновенных коэффициентов естественной смертности (1/год) представлены на рис. 1, по оси т — возраст (в годах). Значения мгновенных коэффициентов естественной смертности (МКЕС) определены по методу Гундерсона-Ди-герта (Gunderson et al., 1988). Кроме того, использованы оценки биомассы нерестового запаса по результатам икорных съемок КамчатНИРО в период с 1973 по 2008 гг.
Мгновенные коэффициенты естественной смертности m(r) аппроксимировались двумя ветвями парабол:
Среднемноголетняя масса w^) восточноохотоморского минтая аппроксимировалась кривой Берталанффи (Bertalanffy, L. von, І938):
Зависимость доли половозрелых рыб от возраста а(т) аппроксимировалась «кусочками» полиномов:
0.7 -і 0.6 -3 05
- 0.4 -ш 0.3 Н 2 0.2 -0.1 -0
0.9 -0.8 -0.7 -0.6 -оз 0.5 -У 0.4 -2 0.3 -0.2 -0.1 -0 --
♦ ♦ ♦
і----------1---------1---------1---------1—
2 4 6 8 10
Возраст рыб, лет
о 1.2 -і -□
а 1 Н
-О
5 0.8 -
Cl
8 п с j
CQ U.D -
о
| 04 ос § 0.2
♦ * *
-І-------------»-г
2 4
6 8 Возраст рыб, лет
10
12
“1 Рис. 1. Среднемноголетние значения МКЕС (1/год)
12 (а), доли половозрелых рыб (b) и массы (кг) (c) вос-
точноохотоморского минтая
РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЕ На рис. 2 и 3 в виде линий уровня представлены ретроспективные оценки плотности биомассы х(ї, т) (в тыс. тонн/год) и интенсивности промысла и(ї, т) (в 1/год) восточноохотоморского минтая по результатам применения непрерывной модели динамики возрастного состава биологической популяции. Коэффициенты Рикера оказались равными а= 22,94, 3 =9,3 х 10 4. При этом значение минимизируемого функционала (2) составило 5779200,018.
На рис. 4 представлена ретроспективная динамика биомассы общего (ТВ) и нерестового (88В) запаса восточноохотоморского минтая по непре-
рывной модели и модели «Synthesis» в сравнении с данными икорных съемок КамчатНИРО. Результаты двух моделей неплохо согласуются между собой и с данными съемок. В период времени с 1996 по 2008 гг. отклонение между значениями биомассы нерестового запаса по данным съемок и модельными значениями не превышает 8,13% для непрерывной модели и 13,7% для когортной модели «Synthesis».
Для определения оптимального режима эксплуатации восточноохотоморского минтая по непрерывной модели динамики возрастной структуры популяции необходимо задать коэффициент
Рис. 2. Ретроспективная динамика возрастного состава восточноохотоморского минтая (плотность биомассы, тыс. тонн/год)
1965 1970 1975 1980 1985 1990 1995 2000
Рис. 3. Ретроспективная интенсивность промысла (1/год) восточноохотоморского минтая
1960 1970 1980 1990 2000 2010
Рис. 4. Ретроспективная динамика биомассы общего (ТВ) и нерестового (ББВ) запаса восточноохотоморского минтая
Рис. 5. Коэффициент улавливаемости для восточноохотоморского минтая
улавливаемости д(т), интенсивность рождения новых организмов, ограничения на величину биомассы нерестового запаса и на количество промысловых усилий.
Из ретроспективных оценок интенсивности промысла и(/, т) и данных о количестве промысловых усилий по формуле (4) определяется коэффициент улавливаемости д(т) для восточноохотоморского минтая. График полученной зависимости коэффициента улавливаемости от возраста рыб представлен на рис. 5.
Интенсивность рождения для восточноохотоморского минтая (плотность пополнения в возрасте 0) в прогнозные годы предполагается равным среднему арифметическому от плотности по-
Возраст рыб, лет
Рис. 6. Плотность биомассы восточноохотоморского минтая при режиме оптимальной эксплуатации (тыс. тонн/год)
полнения за последние 10 предпрогнозных лет (х(ї ,0) = Я = 8,35 млрд рыб/год).
В качестве критического уровня по биомассе для восточноохотоморского минтая выбрано среднее арифметическое значение биомассы нерестового запаса по данным съемок за последние 10 лет (8=1525,6 тыс. тонн).
Ограничением на интенсивность промыслового усилия служит максимальное значение количества промысловых усилий за весь период ведения промысла (Р = 46,3 тыс. с.-с./год).
Был найден оптимальный режим эксплуатации восточноохотоморского минтая на 3 года вперёд (2009-2011 гг.). На рис. 6 представлена плотность биомассы х(ї, т) восточноохотоморского минтая (в тыс. тонн/год) при режиме оптимальной эксплуатации на 2009-2011 гг.
Модельная динамика годовых уловов, биомассы общего (ТВ) и нерестового (88В) запасов восточноохотоморского минтая при оптимальной эксплуатации представлена на рис. 7. Оптимальная интенсивность промыслового усилия и годовое промысловое усилие изображены на рис. 8.
Решение задачи оптимальной эксплуатации для восточноохотоморского минтая таково, что в каждый календарный год запас сначала эксплуатируется с максимальной интенсивностью, а затем изъятие прекращается с тем, чтобы к концу года биомасса нерестового запаса восстановилась до уровня 8. Сохраняя нерестовую биомассу восточноохотоморского минтая к концу каждого прогнозного года на достаточно высоком уровне продуктивности — 1525,6 тыс. тонн, за три года из запаса можно изъять 1279,012 тыс. тонн.
ЗАКЛЮЧЕНИЕ
Представленные результаты моделирования ретроспективной динамики возрастного состава и оптимального режима изъятия восточноохотоморского минтая показывают принципиальную возможность использования непрерывной модели динамики возрастной структуры биологической популяции для оценки текущего состояния запасов и прогнозирования величин ОДУ морских промысловых объектов.
2009 2010 2011 2012
Рис. 7. Вылов, биомасса общего (ТВ) и нерестового (88В) запаса восточноохотоморского минтая при режиме оптимальной эксплуатации
Рис. 8. Интенсивность промыслового усилия и годовое промысловое усилие при режиме оптимальной эксплуатации восточноохотоморского минтая
Алгоритмы нахождения ретроспективных оценок запаса и определения оптимального режима эксплуатации реализованы в виде компьютерной программы.
СПИСОК ЛИТЕРАТУРЫ
Абакумов А.И. 1991. Модельный анализ оптимальных режимов эксплуатации популяций // Управление и оптимизация. Владивосток: ДВО АН СССР. С.3-16.
Абакумов А.И. 1991. Обобщение математических методов оценки запаса для промысловых популяций. Препринт / Интститут прикладной математики ДВО АН СССР. Владивосток, 11 с.
Абакумов А.И., Кольев Н.В., Максименко В.П., Горр С.В. 1994. Матричный метод оценки запаса и прогнозирования вылова популяций морских организмов // Вопросы ихтиологии. Т. 34. № 3. С. 400-407.
Бабаян В.К. 2000. Предосторожный подход к оценке общего допустимого улова (ОДУ). Анализ и рекомендации по применению. М.: ВНИРО, 190 с.
Баранов Ф.И. 1918. К вопросу о биологических основаниях рыбного хозяйства // Изв. отдела рыбоводства и науч.-промысл. исслед. Т. 1. Вып. 1. С. 84-128.
Бивертон Р., Холт С. 1969. Динамика численности промысловых рыб. М.: Пищ. пром-сть, 248 с.
Васильев Ф.П. 1980. Численные методы решения экстремальных задач. М.: Наука, 520 с.
Динамическая теория биологических популяций / Под ред. Р. А. Полуэктова. 1974. М.: Наука, 456 с.
Ильин О.И. 2006. Об одной задаче оптимального управления биологическим сообществом // Сиб. журн. индустриальной математики. Т. 9. № 4 (28). С. 75-81.
Ильин О.И. 2007. Об оптимальной эксплуатации популяции несколькими типами орудий лова // Сб. науч. тр. Камчат. НИИ рыб. хоз-ва и океанографии. Вып. 9. С. 258-260.
Ильин О.И. 2007. Об оптимальной эксплуатации популяции рыб с возрастной структурой // Сиб. журн. индустриальной математики. Т. 10. № 3 (31). С.43-57.
Ильин О.И. 2007. Об одном методе определения оценки запасов и оптимального вылова популяции с возрастной структурой // Математическое моде-
лирование и информационные технологии в исследованиях биоресурсов Мирового океана: Тез. докл. Материалы отраслевого семинара. Владивосток: ОИТ ТИНРО-центра, 84 с.
Ильин О.И. 2007. О численном решении задачи оптимальной эксплуатации промысловой популяции биологических объектов // Тез. докл. математической шк.-семинара имени академика Е.В. Золотова. Владивосток. С. 70-71.
Максименко В.П., Антонов Н.П. 2003. Количественные методы оценки рыбных запасов. Петро-павловск-Камчатский: КамчатНИРО, 256 с.
Перцев Н.В. 1999. Исследование решений интегральной модели Лотки-Вольтерра // Сиб. журн. индустриальной математики. Т. 2. № 2 (4). С. 153-167.
Рикер У.Е. 1979. Методы оценки и интерпретация биологических показателей популяций рыб. М.: Пищ. пром-сть, 408 с.
Свирежев Ю.М., Тимофеев Н.Н. 1980. О регулировании численности популяции с возрастной структурой // Журн. общ. биологии. Т. 41. № 2.
С. 200-209.
Bertalanffy, L. von. 1938. A Quantitative Theory Of Organic Growth. // Human Biology. V. 10. № 1. P. 181-213.
Deriso R.B., P.R.Neal, and T.J. Quinn II. 1985. Catch-age analysis with auxiliary information. // Can. J. Fish. Aquat. Sci. 42. P. 815-824.
Gulland J.A. 1969. Manual of methods for fish assessment. Pt.1 // Fish population analysis. FAO. Nan., Fish. Sci. № 4. P. 1-154.
Gunderson D.R., Dygert PH. 1988. Reproductive effort as a predictor of natural mortality rate // J. Cons. Int. Explor. Mer. 44. P. 200-209.
Kizner Z.I., Vasilyev D.A. 1997. Instantaneous Separable VPA (ISVPA) // ICES Journal Of Marine Science. Vol. 54. № 3. P. 399-411.
Richard D. Methot. 1990. Synthesis model: an adaptable framework for analysis of diverse stock assessment data // International North Pacific Fisheries commission. Bul. № 50. Vancouver, Canada.
Sinko J. W., Streifer W. 1967. A new model for age-sized structure for population // Ecology. V. 48. P. 910-918.
Von Foerster H. 1959. Some remarks on changing populations // The Kinetics of cellular proliferation. N.Y. P. 382-407.
