УДК 577.95:519.233.2
ВРЕМЕННАЯ ИЗМЕНЧИВОСТЬ В КОЭФФИЦИЕНТАХ МОДЕЛИ:
ПРИМЕР BUPALUS PINIARIUS
© 2013 г. Л.В. Недорезов
Центр междисциплинарных исследований по проблемам окружающей среды РАН 191187 С.-Петербург, наб. Кутузова, 14 e-mail l.v.nedorezov@gmail.com Поступила в редакцию 24.10.2011 г.
Если задана модель динамики численности популяции и достаточно длинный временной ряд данных, то параметры модели можно оценивать на разных отрезках этого ряда, получая соответственно разные оценки и модель с переменными во времени параметрами. Вариации параметров модели могут быть вызваны либо вариациями условий среды, влияющих на численности популяции, либо просто демографической стохастичностью. Возможный подход к решению этого вопроса рассмотрен на примере данных по колебаниям численности сосновой пяденицы (Bupalus piniarius L.) в Германии. Для аппроксимации данных используется известная модель Морана - Риккера, которая обладает достаточно богатым набором динамических режимов. Оценка параметров проводилась по методу наименьших квадратов (по 12 значениям из 58). Анализ трендов в двух рядах из коэффициентов модели - максимального коэффициента размножения и коэффициента саморегуляции - показывает, что изменения обоих показателей невелики за довольно долгий промежуток времени (60 лет) и гипотеза об отсутствии направленного тренда в изменениях среды оказывается статистически правдоподобной.
Климатические факторы оказывают существенное влияние на динамику численности насекомых (Исаев и др., 1984, 2001; Воронцов, 1978; Варли и др., Беггушап, 1981, и др.). Это влияние имеет комплексный характер и проявляется как непосредственно, так и через другие компоненты биогеоценоза. Поэтому искать какие-либо прямые коррелятивные связи изменения численности (или плотности) популяции с одним каким-либо климатическим показателем или же группой показателей, по-видимому, бесперспективно для объяснения наблюдаемых колебаний, но целесообразно для построения прогнозов различной заблаговременности (Кондаков, 1974; ^аеу ег а1., 2009; Топпа^ ег а1., 2010).
Изменение условий обитания вида (в том числе и изменение климатических показателей) приводит к соответствующим изменениям популяционных показателей (плодовитости, смертности особей, интенсивности саморегуляции). Но если непосредственное оценивание изменений в условиях обитания является крайне сложной задачей (сюда входят и оценки изменения кормовых условий, и оценки состояний различных регуляторных механизмов, и многое другое), то оценивание популяционных показателей, которые опосредованно выражают влия-
ние условий обитания, может выявить и степень этого влияния, и его направленность. Например, снижение коэффициента внутрипопуляционной конкуренции может свидетельствовать об улучшении кормовых условий, снижение величины максимального коэффициента размножения -об ухудшении (для особей) климатических показателей и/или кормовых условий, и так далее.
Если имеется модель, дающая хорошую аппроксимацию эмпирических данных, то для достаточно длинного временного ряда оценки параметров модели по начальным значениям ряда должны достоверно отличаться от их оценок, полученных для "хвоста" временного ряда (если, конечно, изменение условий оказывает заметное влияние на популяционные показатели). Проведение многократных оценок параметров модели по отдельным частям временного ряда даёт новые временные ряды (оценок параметров модели), в динамике которых можно пытаться выявить тенденцию.
Если исходный временной ряд - это последовательность значений плотности популяции (определяемой, например, ежегодно), то временные ряды коэффициентов модели имеют качест-
ВРЕМЕННАЯ изменчивость в коэффициентах модели
35
венно иную природу. В модели Морана - Риккера (Moran, 1950; Ricker, 1954):
e byi /14
Ук+i = ayke , (1)
где yk - численность (или плотность) популяции в k-й год, параметр a равен максимальному коэффициенту размножения (определяется как отношение численностей двух смежных генераций), а параметр b трактуют как коэффициент саморегуляции. Некоторые исследователи полагают (Исаев и др., 1984), что параметр a определяется климатическими и кормовыми условиями. Поскольку исходная выборка состоит из значений величины, которую естественно считать случайной, оценки параметров модели, определяемые по частям исходной выборки, также являются значениями некоторых случайных величин. Поэтому к этим двум выборкам значений a и b также применимы статистические методы анализа.
Отметим, что выполнение такой программы: оценка параметров модели по подмножествам исходной выборки и оценка тенденций изменения популяционных показателей, встречает ряд серьезных трудностей. Во-первых, как показывает практика использования нелинейных математических моделей для аппроксимации эмпирических данных (Недорезов, Садыкова, 2005; Nedorezov, et al., 2008; Tonnang et al., 2009, 2010; Nedorezov, 2011), даже для относительно коротких временных рядов (10-15 значений) модели типа (1) могут не дать приемлемой аппроксимации, что зависит от поведения отклонений эмпирических данных от теоретических кривых (Дрейпер, Смит, 1986, 1987). Во-вторых, нередко бывает так, что аппроксимация коротких временных рядов с помощью моделей типа (1) требует больших вычислительных затрат (например, вследствие неудачно выбранных начальных значений). В-третьих, если для одной какой-либо части временного ряда модель (1) оказывается приемлемой, а для другой части нет, то возникает естественный вопрос, стоит ли использовать все оценки параметров как индикаторы изменений внешних факторов или же только ту часть ряда, где модель дает удовлетворительную аппроксимацию данных?
В настоящей работе для аппроксимации известных данных Швердтфегера (Schwerdtfeger, 1957, 1968) по колебаниям численности сосновой пяденицы (Bupalus piniarius L.) используется модель Морана - Риккера (1). Оценивание параметров модели проводили по 12-летним временным интервалам с использованием известного метода наименьших квадратов. Временные ряды для коэффициентов модели Морана - Риккера
получались следующим образом: оценивались коэффициенты линейного тренда и для коэффициента углового наклона регрессионной прямой проверялась гипотеза о равенстве этого коэффициента нулю с помощью непараметрического критерия Тейла (Theil, 1950; Холлендер, Вулф, 1983). Показано, что с 5%-ным уровнем значимости гипотеза о равенстве коэффициента углового наклона нулю не может быть отклонена. Это означает, что изменения условий обитания пяденицы за период с 1881 по 1940 год не привели к существенным изменениям популяционных показателей.
МЕТОД АНАЛИЗА
Пусть имеется временной ряд x1, x2,..., xN, N -число лет (объем выборки), xk - плотность популяции в k-й год. Для каждого отрезка xr, xr+1,., xr+m, r, m > 1 (при обработке реальных данных полагалось m = 11), исходного временного ряда параметры модели Морана - Риккера (1) оценивались общеизвестным методом "наименьших квадратов": среди множества всех траекторий модели выбиралась такая, на которой следующий функционал имел минимум:
r+m
Q(a *, b *, у *; r, m) = min V(xy - /(/, a, b, y))2, (2)
a, b, у 1
/ = r
где y1 - начальное значение для плотности популяции в модели (1), /(j, a, b, у1) - соответствующее значение, получаемое по модели (1) при заданных значениях параметров a, b, и начальном зна-ченииу1:/(r, a, у^ =у 1,fr 1, а, b, у^ = у2, и так далее. В (2) через a*, b* и у* обозначены оценки параметров (начальное значение также считается неизвестным параметром, подлежащим оценке по имеющейся выборке), при которых достигается минимум функционала.
Аппроксимируя все указанные подмножества временного ряда, получаем N - m значений параметров a и b. Фактически мы получаем два новых временных ряда: a1,., aN-m и b1,., bN-m, где в качестве индексов r, r = 1,., N - m выбраны номера первых элементов соответствующих отрезков исходного ряда: оценки параметров модели a1 и b1 получены при минимизации функционала (2) по подмножеству x1, x2,..., x12; оценки a2 и b2 получены при минимизации функционала (2) по подмножеству x2, x3,..., x13 и так далее. Таким образом, индекс r соответствует моменту времени и задаёт определённый (при выбранном значении m) отрезок времени, для которого и оцениваются параметры модели.
Однако использовать все полученные значения этих новых временных рядов, видимо, нельзя. Это связано с тем, что для одного отрезка временного ряда модель (1) может дать удовлетворительную аппроксимацию, а для другого может оказаться и непригодной. Последнее зависит от того, как именно ведут себя последовательности отклонений эмпирических данных от теоретических (модельных), полученных при наилучших значениях параметров модели (1).
В стандартном случае отклонения эмпирических от теоретических данных должны быть нормально распределены со средним арифметическим равным нулю (точнее, соответствующие гипотезы не должны отклоняться при выбранном уровне значимости). Для проверки использовались критерии Колмогорова - Смирнова, Лил-лиефорса и Шапиро - Вилка (Большев, Смирнов, 1983; Lilliefors, 1967; Shapiro et al., 1968). Кроме того, в последовательности отклонений не должно быть сериальной корреляции (Дрейпер, Смит, 1986, 1987). Если какой-либо из используемых статистических критериев позволял отклонить соответствующую гипотезу (в том числе гипотезу о равенстве среднего нулю, об отсутствии сериальной корреляции и так далее), то модель (1) признавалась непригодной для аппроксимации соответствующего подмножества временного ряда. В противном случае модель признавалась пригодной для аппроксимации.
Представляет интерес анализ тенденций в изменениях значений параметров a и b. Значения элементов во временных рядах a1, ..., aN-m и b1, ..., bN-m являются случайными, поскольку определяются по исходному случайному временному ряду. Высказать какую-либо правдоподобную гипотезу об отклонениях элементов этих рядов от "истинных значений" не представляется возможным, и потому для выявления тенденций в этих рядах использовался непараметрический критерий Тейла (Theil, 1950; Холлендер, Вулф, 1983). Если этот
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.