CC BY
34Моделирование биологических популяций с учетом свойств отдельной особи
Сорокин П.А. ( pavel.sorokin@mobico.ru )
Центр Открытых Систем и Высоких Технологий МФТИ
Введение
Эколого-биологическая область - одна из наиболее перпективных, прогресс в этой области зависит от применения методов математического моделирования, использования вычислительных технологий. Одной из центральных задач биомоделирования является моделирование биологических популяций. Моделирование используется для прогнозирования их численности, изучения миграций, планирования полевых исследований.
На данный момент накоплен обширный опыт имитационного моделирования популяций. Обычно модель представляет собой систему уравнений феноменологического типа, описывающими связи между макропоказателями (например, численностью популяции) и обобщенными экологическими характеристиками. Характерным примером такого подхода является хорошо известная модель "хищник-жертва" Лотки-Вольтерра, где численность популяции хищников и жертв описывается системой двух обыкновенных дифференциальных уравнений, в которых наряду с численностью используются такие обобщенные показатели, как коэффициент размножения жертв, смертности хищников, обобщенный показатель мх взаимодействия. Во всех таких моделях описываются свойства популяции в целом.
В данной работе предложен другой подход, при котором в модели описываются свойства отдельной особи, а свойства популяции являются интегральными характеристиками модели. Это позволяет учитывать ряд дополнительных аспектов поведения популяции (эколого-физиологических, социально-территориальных), не описываемых в рамках традиционых моделей. В то же время процесс настройки и использования модели становится более наглядным, что дает возможность строить системы моделирования, которыми смогут пользоваться биологи и экологи, не являющиеся специалистами в области математики и моделирования. Постановка задачи
Основная потребность современного этапа развития биолого-экологического моделирования -возврат доверия биологов к исследованиям математиков. Поэтому рассматривалась задача построить модель, позволяющую наладить диалог с биологами, при обсуждении которой можно было оперировать терминами, понятными им ("разговаривать" на языке биологии, а не математики.)
Один из существенных недостатков аналитических моделей - сложность их приближения к биологическим реалиям, отсутствие прозрачности в установлении отношения между моделью и реальной экосистемой. В связи с этим в последнее время приобретает все большую популярность имитационное моделирование, которое стимулируется бурным развитием вычислительной техники. В то же время, как и во всякой новой области, традиции имитационного моделирования еще не устроялись, существует потребность в новых подходах.
Полученный опыт при моделировании тундровых популяций и сообществ, направленный на анализ конкретных эколого-биологических процессов во взаимодействии с биологами, привел к формализованной постановке задачи, которая позволяет учитывать специфику полевых исследований, дает возможность учесть экспертные оценки на стадии постановки задачи и обсудить результаты в привычной для специалистов-биологов терминологии.
Основной целью исследования ставилось построить модель, которая позволила бы анализировать биофизические механизмы, лежащие в основе формирования динамики численности популяции, была бы открыта для модификаций и формулировалась бы в терминах, позволяющих экспертам-
биологам участвовать в обсуждении и формировании модели. Ставилась задача преодолеть ограниченность аналитических моделей, не позволяющих эффективно учитывать данные о популяциях, полученные в результате биолого-экологических исследований. Возможность анализа ведущих механизмов явления дает ключ к обоснованию аналитических подходов, сложившихся в популяционном моделировании. Построенная модель должна была показать возможность такого обоснования и таким образом установить связь конкретных биофизических механизмов, действующих внутри популяций, с параметрами популяционных уравнений.
Данная работа являлась продолжением комплекса исследований по созданию набора математических моделей для описания эколого-биологических процессов. Схематически эволюция моделей показана на рис. 1. Развитие моделей шло по пути все большей детализации. От моделей, описывающих биосферные процессы, осуществлялся переход к моделям экосистем и популяций.
Обоснование по пуляцио нных
уравнений _
Особь
Рисунок 1. Эволюция моделей.
Логичным завершением процесса детализации описание является построение модели, учитывающей свойство отдельной особи.
Схематически комплекс факторов, существенных для анализа биофизических механизмов, формирующих динамику численности, показан на рис. 2.
В модели необходимо было учитывать:
• эколого-физиолигические свойства особей
• особенности взаимодействия особей
• влияние на их поведение окружающей среды (в том числе пространственных особенностей ареала)
• сезонные факторы
Рисунок 2. Свойства особи.
Объект моделирования
Для моделирования была выбрана экосистема тундры. Рассматривалась задача изучения динамики численности тундровых животных. Модель настраивалась на данных по численности популяции леммингов. Численность лемминговой популяции является фактором, определяющим численность основного промыслового животного тундры - песца. Кроме того она интересна с научной точки зрения сама по себе. Численность лемминговых популяций подвержена значительным колебаниям; она может меняться из года в год в несколько сот раз. Для большинства компонентов тундровых экосистем последствия этих вспышек весьма существенны. Между тем, механизмы регуляции численности леммингов остаются мало изученными. Поэтому исследование особенностей динамики численности леммингов может помочь раскрыть экологические закономерности, приводящие к колебаниям численности животных.
Лемминги - основное растительноядное животное тундры весом 80-100 грамм. Новорожденные лемминги слепые и почти голые, их вес порядка 5,0-5,5 грамм. Сроки прозревания около 2 недель, это период детеныши проводят в норе. Начало сезона активного размножения февраль-март, конец -август-сентябрь, с перерывом в весенний период (часть мая или июня). Лемминги живут один год. При размножении в помете обычно бывает в среднем 6 детенышей. Период полового созревания наступает в возрасте 1-4.5 месяца в зависимости от плотности популяции и условий внешней среды. Данные преведены по материалам В. А. Орлова [10]. Модели популяции леммингов рассматривались в работах [1,10,12]. Описание модели
Модель строилась для анализа динамики численности популяции леммингов. Для настройки модели использованы реальные данные о численности популяции, зарегистрированные на острове
Врангеля. Тем не менее предлагаемый подход к моделированию в достаточной степени независим от конкретного объекта моделирования и может использоваться для решения широкого класса задач. Поэтому при описании модели насколько это возможно применялась терминология ресурс-потребитель, которая может быть использована в самых различных предметных областях. При описании моментов, специфичных для популяции леммингов, использовались естественные биологические понятия.
Одним из существенных отличий от ранее строившихся моделей является наличие в модели пространственного распределения популяции. Это дает в руки исследователям принципиальную возможность анализировать влияние геометрии ареала на популяцию. Возможно так же изучать пространственную плотность популяции.
Краткое описание модели
В модель включены потребители - лемминги и ресурс - растительность. Пространственные аспекты в модели являются существенными, для моделирования используется разновидность клеточного автомата. Модель функционирует следующим образом.
Пространство моделируется прямоугольной таблицей. Состояние каждой клетки таблицы описывается величиной ресурса, и списком потребителей, находящихся в данной клетке. Каждый потребитель описывается набором свойств. К ним относятся возраст, жизненный потенциал (характеризует жизнеспособность особи), и способность к размножению.
В клетке, в которой он находится, он потребляет некоторое количество ресурса. Если потребляется большее количество ресурса, жизненный потенциал увеличивается. В противном случае потенциал уменьшается (но не может стать отрицательным). Увеличивается потенциал не безгранично, а до некоторого порогового значения, имеющего смысл максимальной жизнестойкости потребителя. Он зависит от возраста. В модели использовалась экспоненциальная зависимость от возраста. Если жизненный потенциал потребителя падает до 0, он исключается из модели (гибнет). С течением времени ресурс восстанавливается.
Между потребителями существует конкуренция за обладание ресурсом. Поэтому увеличение плотности популяции приводит к снижению жизнестойкости. Если два потребителя находятся в одной клетке, происходит уменьшение жизненного потенциала на некоторую величину.
На каждой итерации значение способности к размножению увеличивается на величину, зависящую от жизненного потенциала. Более жизнестойкие потребители быстрее наращивают способность к произведению потомства. Когда способность к размножению увеличивается до некоторого порогового значения, происходит появление потомства. При этом значение способности к размножению полагается равным 0, и процесс повторяется. После трех размножений потребитель гибнет.
На каждой итерации происходит перемещение потребителей. Они происходят с учетом благоприятности условий обитания в различных клетках. Степень благоприятности клетки описывается величиной штрафа. Штраф для каждой клетки начисляется в зависимости от плотности популяции и количества растительности.
Если возраст потребителя меньше некоторого значения, он не передвигается, не потребляет ресурс, не участвует в столкновениях, не размножается. Для леммингов это возраст достижения половой зрелости. Однако существование такого потребителя учитывается в плотностной составляющей штрафа. Тем самым в модели введено две возрастные группы потребителей.
Существенным обстоятельством в формировании колебаний численности является наличие сезонов (зимы и лета). В ряде построенных ранее популяционных моделей колебания численности были получены именно путем введения в модель двух различных режимов (соответствующих двум сезонам) с различными значениями коэффициентов модели. В предложенной модели различаются три сезона - лето, перезимовка без размножения, подснежное размножение. В летнем режиме модель функционирует, как описано выше.
Поведение леммингов в зимний период имеет некоторые отличия. Под снегом существует оптимальный для жизнедеятельности биотоп, в котором лемминги гарантированно выживают.
Выживание леммингов за пределами этого биотопа зависит от жизнестойкости особи (жизненного потенциала), которая определяется в первую очередь плотностью популяции.
Подробное описание модели
Таблица, моделирующая пространство, имеет размеры М х N клеток. Величина ресурса в клетке с координатами (т, п) - /тп (/тп е [0,1]). С течением времени ресурс восстанавливается согласно уравнению Ферхюльста
1тп ($ +1) = 1тп(г) * (Ь - 1тп (г)) (1)
Значение коэффициента Ь принималось равным 2 (при этом /тп = 1 является устойчивым значением). Здесь и далее г - целочисленная переменная, обозначающая номер итерации. Потребитель описывается возрастом а., жизненным потенциалом р., и способностью к
размножению di. В клетке, в которой он находится (обозначим ее (т, п)), он потребляет ресурс в количестве А/ { , причем выполняются соотношения
А/ =Аfo * (г) (2)
Л»(г +1) = /тп(г) - А/. (г) = (1 - А/0) * /тп (г) (3)
Здесь А/0 - некоторая постоянная. Поскольку А/0 принимает значения от 0 до 1, значение этой постоянной равно максимальной биомассе, выедаемой леммингом за одну итерацию. Жизненный потенциал потребителя изменяется согласно уравнению (4):
рг (г +1) = рг (г) + А/. - А/шт (4)
Значение pt принадлежит отрезку [0,1], в случае, если значение pt (t +1), рассчитанное по
формуле (4), выходит за пределы этого отрезка, устанавливается соответствующее граничное значение. Величина A/min имеет смысл минимального количества ресурса, необходимого потребителю для поддержания нормальной жизнедеятельности в течение одной итерации. Если потребляется большее количество ресурса, жизненный потенциал увеличивается (но не превосходит pmax. В противном случае потенциал уменьшается (но не может стать отрицательным).
Полагается А/0 > A/min, и при fmn, близких к 1, происходит увеличение жизненного потенциала. Таким образом, наличие ресурса обуславливает восстановление потребителей.
Величина pmax имеет смысл максимальной жизнестойкости потребителя и зависит от возраста. В модели использовалась экспоненциальная зависимость от возраста: 0 ln10 -
Pmax = eXP|- 2—Г * a I (5)
Здесь l - число итераций в периоде, соответствующем одному году, a - возраст потребителя, а
коэффициент ln10 выбран из того условия, чтобы в возрасте двух лет максимальный жизненный 2 * l
потенциал был равен 0.1.
Механизм влияния плотности популяции на жизнестойкость учитывается следующим образом. При попадании двух потребителей в одну клетку, происходит уменьшение жизненного потенциала на величину Ар0.
рг (t +1) = min(р^ (t) -Ap0,0) (6)
На каждой итерации значение способности к размножению увеличивается на величину Дdi, причем выполняются соотношения
Дdг (г) = Мо * р.(г) (7)
Дdг (г +1) = Ш8Х(^,. (г) + Дdг (г),1) (8)
Здесь Д^0 - некоторая постоянная. Поскольку р. принимает значения от 0 до 1, значение этой постоянной равно максимальному значению, на которое может увеличиться способность к размножению за одну итерацию. Начальное значение di полагается равным 0. При достижении
значения di = 1 происходит появление потомства количеством пс потребителей (для леммингов пс = 3 принималось как характерное число детенышей в помете). После появления потомства значение di полагается равным 0, и процесс повторяется. Штраф для каждой клетки вычисляется по формуле
Р = Р0 + Рп + Рр (9)
тп тп тп тп \ /
где РПп - плотностная составляющая штрафа, Р0Рп - составляющая учитывающая количество растительности, Р°п - случайное возмущение. Рр линейно зависит от биомассы растительности (Р0Рп =-/ тп), РОПп прямо пропорционально числу потребителей в соседних клетках (поскольку
потребитель на каждой итерации совершает движение, на следующей итерации в данной клетке могут оказаться потребители из соседних клеток).
Значение коэффициента пропорциональности - один из настраиваемых параметров модели. Его величина определяет влияние плотности потребителей на условия существования. Р 0п - малая
величина, постоянная для каждой клетки. Она вводится для того, чтобы сделать клетки изначально неравноправными, и обеспечить случайное перемещение потребителей при малых плотностях популяции.
Приведем наглядный пример функционирования модели. На рис. 3 показан один и тот же фрагмент таблицы для двух последовательных итераций модели. Кружками помечены потребители. Цветом помечено количество ресурса (чем темнее клетка, тем больше ресурс). Цифрами обозначены координаты клеток. Рассмотрим функционирование модели на примере потребителя в клетке (2,2)(на рисунке слева). Он может переместиться в одну из клеток: (1,2), (2,1), (2,3) или (3,2). Выбор клетки осуществляется путем определения значения штрафа для каждой из них согласно уравнению (10). Составляющая Р0Рп пропорциональна количеству растительности, Р0Пп - числу потребителей в соседних клетках, Р°п постоянна. Для каждой клетки штраф вычисляется по формуле
Ртп = Щ/тп +а 2итп + Р1п (10)
Клетка (1,2): Р12 = а 1 * 0.80 +а2 * 0 + 0.00024 Клетка (2,1): Р21 =а 1 * 0.83 +а2 * 2 + 0.00037 Клетка (3,2): Р32 = а1 * 0.90 +а2 * 1 + 0.00012 Клетка (2,3): Р23 = а1 * 0.90 +а2 * 2 + 0.00062
Если, например, а1 = 1.0 и а 2 = 0.3, потребитель переместится в клетку (1,2). Расположение потребителей на следующей итерации показано на рисунке справа.
• • •
• 1 1
• • 2 • • 2
• 3 • 3
1 2 3 1 2 3
а) б)
Рисунок 3. Пример функционирования модели.
Возраст достижения половой зрелости а0. Хотя неполовозрелый лемминг не участвует в
потреблении растительности, размножении и столкновениях, еуществование такого потребителя учитывается в плотностной составляющей штрафа.
В зимний период таблица разбивается на три части (соответствующие трем биотопам, в различной степени благоприятным для проживания). Для каждой части устанавливается порог жизненного потенциала, при котором лемминг выживает (порог выживания). При переходе в режим моделирования перезимовки без размножения переназначаются штрафы клеток. Появляется дополнительная составляющая штрафа Р^п, которая линейно зависит от порогов выживания и
убывает по направлению к более благоприятным частям таблицы. Таким образом, зимой штраф ячейки вычисляется по формуле
Р = р0 + рп + рр + р£ (11)
п п п п п
Дальнейшее поведение леммингов аналогично их поведению в летний период, за исключением того, что размножения не происходит (соответственно, значение способности к размножению не меняется). Отсутствует так же восстановление растительности.
В конце этого периода сравнивается жизненный потенциал каждого лемминга с пороговым значением потенциала выживания для ареала, в котором он находится. Если жизненный потенциал меньше порогового значения, лемминг гибнет.
В период подснежного размножения начинается размножение леммингов. В остальном модель функционирует, как в период перезимовки без размножения.
Результаты моделирования
Наиболее важной характеристикой популяции является ее численность. На рис. 4 расчетные значения численности (сплошная линия) соотнесены с наблюдаемыми в природе (помечены отдельными точками).
На приведенном графике нанесены только те точки, которые попадают в периоды, за которые имеются данные о природной численности популяции (июнь-август). Полный график численности, полученный с помощью модели, и из которого выбраны данные для построения графика рис. 4, приведен на графике рис. 5.
N зоо
250
200
100 ■■
50
о щ
■■ / ° /и°о
' ° X- п 3 X О I °
.^О ° оо I о О У о О —^
1 1 1 1 1 1972 I I I I 1974 I I I I I I 1976 I I I I I I I I 1978 1980
Рисунок 4. Результаты одного из имитационных экспериментов с моделью
(сплошная линия) и зарегистрированная на острове Врангеля динамика изменения численности копытного лемминга (отдельные точки на графике).
Рисунок 5. Полный график динамики численности.
В природе обычно наблюдаются 3-х либо 4-х летние циклы. Четырехлетний цикл четко выделяется на рис. 5. Трехлетние циклы также были получены в экспериментах с моделью.
Численность популяции является результатом двух процессов - размножения и вымирания. Интенсивность этих процессов является одной из характеристик популяции. Модель позволяет изучать интенсивность размножения и отмирания потребителей. Результаты такого исследования приведены на рис. 6. На графике показаны суммарные (от начала периода) числа отмерших и родившихся потребителей.
Для анализа процессов формирования динамики численности полезно рассмотреть изменение биомассы растительности и среднего жизненного потенциала в сопоставлении с численностью потребителей (рис. 7, рис. 8).
Можно также изучать влияние численности популяции на средний потенциал. Для этого достаточно данные, представленные на графике рис. 8 представить в виде графика зависимости среднего потенциала от численности. Результат такого преобразования показан на рис. 9 и рис. 10. На рис. 9 использованы все данные графика рис. 8. На рис. 10 взяты данные за период с мая 02 по май 03.
Год
Числетость Прирост Смертность
Рисунок 6. Графики интенсивности размножения и отмирания.
Рисунок 7. Динамика изменения биомассы растительности в сравнении с численностью популяции леммингов.
200 150 100 50
А
--1-1-1- -1-1-1-
Год
•Численность
■Потенциал
Рисунок 8. Динамика изменения среднего жизненного потенциала в сравнении с численностью популяции леммингов.
Рисунок 9. Зависимость среднего жизненного потенциала от численности популяции леммингов.
Рисунок 10. Зависимость среднего жизненного потенциала от численности популяции леммингов.
На графике рис. 8 прослеживается характерная зависимость среднего потенциала отдельной особи от численности популяции. На тех участках, где численность популяции невелика, средний
потенциал более высокий и изменяется без какой-либо четко выраженной тенденции (участки с мая 01 по май 02 и с мая 03 по май 04). Там же, где численность популяции существенно возрастает, средний потенциал имеет ярко выраженную тенденцию к уменьшению (участок с мая 02 по май 03). Отсутствие заметной корреляции среднего значения жизненного потенциала и численности при небольших значениях численности свидетельствует о том, что при малой плотности популяции взаимное влияние особей незначительно. При увеличении плотности популяции это влияние усиливается и приводит к заметному уменьшению среднего потенциала. Однако, даже на тех участках, где потенциал не имеет четко выраженной тенденции, его изменение нельзя назвать хаотическим. Здесь график состоит из достаточно продолжительных промежутков монотонности, что свидетельствует об инерционности жизненного потенциала. Вероятно, это связано с тем, что на очередную итерацию модели переходит значение потенциала, которое было на предыдущей итерации. На рис. 9 использованы все данные графика рис. 8. Наблюдается тенденция к уменьшению среднего потенциала при увеличении численности (на графике проведена линия тренда). Однако данные сильно зашумлены. Действительно, ведь при сравнительно небольших численностях популяции не наблюдается явной зависимости между численностью и потенциалом. Зависимость проявляется только при больших численностях, и это видно на графике. На рис. 10 взяты данные за период с мая 02 по май 03, когда зависимость между численностью и потенциалом хорошо прослеживается. Поэтому зашумленность графика заметно меньше.
Можно также изучать численность неполовозрелых леммингов. Это лемминги, возраст которых меньше а0 - они не передвигаются и влияют только на плотность популяции. На рис. 11 приведена
численность неполовозрелых особей в сравнении с общей численностью популяции.
Рисунок 11. Численность неполовозрелых особей.
Из графиков видно, что размножение происходит скачкообразно, то есть можно выделить промежутки времени, когда размножение наиболее интенсивно.
Одной из задач исследования было получение из свойств особи популяционного уравнения. Под популяционным уравнением понимается уравнение, которое позволило бы по текущей численности популяции получить прогноз численности через некоторый промежуток времени. Одно из возможных уравнений - это функция последования. Под функцией последования понимается зависимость численности популяции в конце некоторого периода времени от численности популяции начале этого периода. Имеет смысл исследовать несколько функций последования ¥х, отличающихся периодом времени т . Практическое значение имеют зимняя, летняя и годовая функции последования. Каждая из них позволяет получить численность в конце соответствующего периода в зависимости от численности в его начале.
На рис. 12 приведена характерная зимняя функция последования. На рис. 13 приведена характерная летняя функция последования. На рис. 14 приведена характерная годовая функция последования.
Рисунок 12. Зимняя функция последования.
Рисунок 13. Летняя функция последования.
Рисунок 14. Годовая функция последования.
Функция последования является удобным инструментом исследования численности популяции. В общем случае она не является однозначной, но если однозначность наблюдается, мы получаем простейший алгоритм прогнозирования численности популяции. Подставляя значение численности в начале периода, получаем значение численности в конце периода. Функция последования, таким образом, является "сверткой" сложных алгоритмов имитационного моделирования, требующих для своей работы значительных ресурсов. Процедура прогнозирования при этом значительно упрощается. Сначала, используя сложные ресурсоемкие алгоритмы, мы строим уравнение популяции (в нашем случае это функция последования). Затем используем простой алгоритм определения численности популяции с помощью функции последования.
Проводимые ранее исследования показали, что часто удается построить однозначную функцию последования, что стимулирует использование механизма функций последования для решения задачи прогнозирования численности популяций.
На рис. 12 хорошо видно, что до некоторого порогового значения начальной численности конечная численность популяции пропорциональна начальной, а после порогового значения остается постоянной, равной емкости оптимального биотопа.
На рис. 13 до некоторого предела зависимость так же линейная. При дальнейшем увеличении начальной численности наступает насыщение - конечная численность остается практически постоянной. Такое поведение функции последования обусловлено, очевидно, ограниченностью ареала обитания.
Годовая функция, подобная представленной на рис. 14 была получена и в других работах (см. [10,12]). Сходные свойства численности популяций наблюдаются в природе. Перспективы развития
Модель будет полезна при планировании полевых исследований (например, можно пытаться прогнозировать плотность популяции в различных точках ареала).
Развитие этой задачи требует привлечения методов, используемых в статистической механике, статистических экспериментов с разреженными газами. Привлекательна возможность использования предложенного подхода к проверке диффузионных подходов. Такие характеристики, как интенсивность столкновений, плотность потока особей могут быть получены из модели на основе клеточных автоматов естественным образом.
Интересна возможность учитывать обучение особей. Действительно, моделирование отдельной особи позволяет хранить информацию о ее индивидуальной истории, что значительно приближает модель к реальной дейтвительности.
Перспективным представляется использование подобных моделей для построения так называемого "биологического конструктора" - программного комплекса для моделирования реальных ареалов с возможностью наглядного представления исходных данных. При использовании моделей, в основе которых лежат уравнения того или иного вида, неизменно встает проблема обоснования и настройки коэффициентов модели, поскольку модель настраивается как на свойства популяции, так и на свойства конкретного ареала (в частности, его геометрию). Биологический конструктор же позволит описать конкретный ареал явно (буквально, нарисовать его) и производить моделирование для данного конкретного ареала без дополнительных усилий по настройке модели. Это возможно, так как модель настраивается только на свойства популяции, свойства же ареала задаются в модели явно. Важно так же отметить, что с биологическим конструктором могут работать не только квалифицированные математики, но и специалисты, не имеющие математического образования (в частности, биологи). Это возможно, поскольку взаимодействие человека с биологическим конструктором можно организовать на интуитивно понятном уровне, без оперирования математическими терминами и понятиями. При этом задача математиков - один единственный раз решить задачу настройки, настроив используемые модели на свойства конкретной популяции (вне зависимости от ареала обитания). Практическая же работа с построенными моделями производится непосредственно специалистами-биологами.
ЛИТЕРАТУРА
1. Байбиков Е.В., Белотелов Н.В., Завьялова С.В., Обридко И.В., Орлов В.А., Саранча Д.А., Шелепова О.В., Шиляева Л.М. О Моделировании тундровых популяций и сообществ. Математическое моделирование. Процессы в сложных экономических и экологических системах. // М.: Наука, 1986. с.207-219.
2. Н.Н. Боголюбов. Проблемы динамической теории в статистической физике. // М.-Л.: Гостехиздат, 1946, 119 c.
3. В. Вольтерра. Математическая теория борьбы за существование. // М.: Наука, 1976, 286 с.
4. Гаузе Г.Ф. Математический подход к проблемам борьбы за существование. // Зоол. журнал, 1933, т.12, с.170-177.
5. Гаузе Г.Ф., Витт А. А. О периодических колебаниях численности популяций. Математическая теория релаксационного взаимодействия между хищниками и жертвами и ее применение к популяции двух простейших. // Изв. АН СССР, Отд. мед., мат., и естеств. наук, 1934, N10, с.1555-1559.
6. Дородницын А.А. Математика и описательные науки. // М.: Знание, 1982, с.6-15
7. Колмогоров A.H Качественное изучение математических моделей динамики популяций. // М.: Наука, 1972, 254 p. (Сб. научных тр. // Проблемы кибернетики, вып. 25, с.100-106).
8. Костицын В. А. Эволюция атмосферы, биосферы и климата. // М.: Наука , 1984, 96 с.
9. Одум Ю. Основы экологии. // М.: Мир, 1975, 740 с.
10. Орлов В.А., Саранча Д.А., Шелепова О.А. Математическая модель численности популяции леммингов (Lemmus, Diorostonyx) и ее использование для описания популяций Восточного Таймыра. // "Экология", 1986. N2. С. 43-51.
11. Р. А. Полуэктов, Ю.А. Пых. Динамика однолокусных полиаллельных популяций. Динамика нейтральных популяций. // Динамическая теория биологических популяций." // М.: Наука, 1974, с.~241-349.
12. Саранча~Д.А. Материалы по количественной экологии. Математическое моделирование и биофизические аспекты. // М.: "Биомоделирование", ВЦ РАН, 1995.
13. Ю.М. Свирежев, Д.О. Логофет. Устойчивость биологических сообществ. // М.: Наука, 1978, 352 с.
14. Форрестер Дж. Мировая динамика. // М.: Наука , 1978, 268 с.
15. Чернявский Ф.Б., Ткачев А.В. Популяционные циклы леммингов в Арктике (экологические и эндокринные аспекты). // М.: Наука, 1982, 164 с.
16. А.М. Шарковский. Разностные уравнения и динамика численности популяций. // Киев: Ин-т математики АН УССР, 1982, 24 с.
17. А.М. Шарковский. Существование циклов непрерывного преобразования прямой в себя. // Укр. мат. журнал, 1964, т.16, н.1 с.61-65
18. C. S. Elton. Periodic fluctuations in the numbers of animals: their causes and effects. // BritJ.Exper.Biol., v.2, P.119-163, 1924.
19. G. B. Ermentrout, L. Edelstein-Keshet. Cellular automata approaches to biological modeling. // Journal of Theoretical Biology, 160:97-133, 1993.
20. E. Fredkin, T. Toffoli. Conservative logic. // International Journal of Theoretical Physics, 21:219253, 1982.
21. E.H. Kerner. A Statistical Mechanics of Interacting Biological Species. // Bull.Math.Biophys., 1957, v.19, no.2 p.121-143.
22. A. J. Lotka. Elements of physical biology. // Baltimore: Williams and Wilkins, 1925, 368 p.
23. A. J. Lotka. Elements of mathematical biology. // N.Y., 1956, 460 p.
24. J. von Neumann. Theory of Self-Reproducing Automata. // Illinois: University of Illinois Press, 1966. (Edited and completed by A.W. Burks).
25. A. Okubo. Diffusion and Ecological Problems: Mathematical Models. // N.Y. etc.: Springer-Verlag, 1980, 254 p. (Biomathematics; Vol. 10).
26. T. Toffoli, N. Margolus. Cellular Automata Machines. // Cambridge, Massachusetts: The MIT Press, 1987.
27. T. Toffoli. Cellular automata mechanics. Technical Report 208. // The University of Michigan: Comp. Comm. Sci. Dept., 1977.
28. S. Wolfram. Universality and complexity in cellular automata. // Physica D, 10:1-35, 1984.
Статья поступила в редакцию 11.07.2002 г.
