CC BY
23
УДК 528.88
СПУТНИКОВЫЕ ДАННЫЕ В ПРОГНОЗЕ УРОЖАЙНОСТИ ПШЕНИЦЫ ДЛЯ АДМИНИСТРАТИВНЫХ РАЙОНОВ СИБИРИ
Сергей Митрофанович Кононенко
Сибирский институт управления - филиал Российской академии народного хозяйства и государственной службы при Президенте РФ, 630102, Россия, г. Новосибирск, ул. Нижегородская, 6, кандидат физико-математических наук, доцент, тел. (383)210-17-15, e-mail: sergej_k@cn.ru
Таисия Васильевна Старостина
Сибирский научно-исследовательский гидрометеорологический институт, 630090, Россия, Новосибирск, ул. Советская, 30, кандидат сельскохозяйственных наук, ведущий научный сотрудник, тел. (383)222-25-30, e-mail: starostina@sibnigmi.ru
Разработаны модели прогноза урожайности яровой пшеницы для некоторых административных районов Сибирского региона. Используются стандартные приземные агрометеорологические данные и спутниковая информация о вегетационном индексе NDVI. Проверено качество полученных моделей и проведены авторские испытания на независимой выборке. Лучшие модели переданы для производственных испытаний и для практического использования.
Ключевые слова: вегетационный индекс, агрометеорологические наблюдения, модель прогноза урожайности, травы.
SATELLITE DATA IN FORECAST FOR WHEAT YIELD FOR ADMINISTRATIVE AREAS OF SIBERIA
Sergey M. Kononenko
Siberian Institute of Management - a branch of the Russian Academy of National Economy and Public Administration under the President of the Russian Federation, 6, Nizhegorodskaya St., Novosibirsk, 630102, Russia, Ph. D., Associate Professor, phone: (383)210-17-15, e-mail: sergej_k@cn.ru
Taisia V. Starostina
Siberian Research Hydrometeorological Institute, 30, Soviet St., Novosibirsk, 630090, Russia, Ph. D., Leading Researcher, phone: (383)222-25-30, e-mail: starostina@sibnigmi.ru
The models of spring wheat yield forecast for some administrative regions of the Siberian region have been developed. Standard surface agrometeorological data and satellite information on the vegetative index NDVI are used. The quality of the received models was checked and author's tests were carried out on an independent sample. The best models are transferred for production tests and for practical use.
Key words: vegetation index, agrometeorological observations, yield forecast model, grasses.
1. Введение, обзор
Прогноз урожайности яровой пшеницы рассмотрен для некоторых административных районов Сибирского региона. Модели используют данные стандартных агрометеорологических наблюдений (осадки, температуру и дефицит
насыщения воздуха) и новую информацию дистанционного зондирования. Вовремя предоставленный прогноз администрациям районов помогает уменьшить экономические потери, связанные с рисками природного характера. В данной работе делается оценка влияния различных факторов на точность прогноза.
Появление новых модулей спутникового зондирования поверхности земли с возможностью измерять и обрабатывать различные участки спектра привело к новым показателям вегетационного периода роста растений. К ним относится нормализованный относительный индекс растительности NDVI (Normalized Difference Vegetation Index), показатель количества активной биомассы и самый распространенный среди подобных ему индексов. Он используется в различных работах, связанных с дистанционным зондированием поверхности земли и растительного покрова. Накопленная в Институте космических исследований РАН по всей территории России информация по этому индексу с реальной регулярностью 3-4 раза в неделю и разрешением 250 м архивирована с 2000г., дала начало новому направлению исследований.
В [1] дается краткий обзор теоретических и прикладных разработок этого направления, проводимых в Институте космических исследований РАН. Методом поиска аналогов временного хода вегетационных индексов разработаны автоматизированные мониторинг и анализ состояния посевов в любом регионе России, оперативная и объективная оценка последствий воздействия на посевы различных факторов. В работе описана технология построения для любого административного района карт отклонений вегетационного индекса от средне-многолетних значений. Графики временной изменчивости значений индекса можно строить в указанных точках земной поверхности. Мониторинг состояния посевов даётся для субъектов РФ.
В [2] вегетационный индекс используется вместе со стандартными данными агрометеорологических наблюдений для построения регрессионных моделей и делается ежедекадный расчет урожайность озимых культур по территории субъектов Северо-Кавказского, Приволжского и Центрально-Черноземного регионов. По оценке авторов ошибка прогноза не превышала 10 % , что даёт возможность использовать данную технологию в оперативном обеспечении сельского хозяйства УГМС Росгидромета.
В работе [3] авторы заменили в своей динамической модели прогноза урожайности рассчитанную относительную площадь листьев посевов на значения измеренного вегетационного индекса NDVI, получив усовершенствованную модель, с помощью которой сделали следующие выводы - величина NDVI достаточно адекватно отражает ход фотосинтеза всех сельскохозяйственных культур; использование NDVI в блоке расчета фотосинтеза в действующей модели прогнозирования урожайности яровой пшеницы дало положительные результаты; замена теоретической кривой сезонного хода фотосинтеза на измеренные значения вегетационного индекса, привела к повышению оправдывае-мости в годы с экстремальными условиями. Полученные результаты дают основание утверждать, что вегетационный индекс NDVI может быть использован
как во вновь создаваемых динамических моделях прогнозирования урожайности, так и в действующих моделях для их усовершенствования.
В качестве примера использования вегетационного индекса NDVI в ЕС можно привести материалы регулярного бюллетеня [4], издающегося по проекту Европейской комиссии MARS (Monitoring Agricultural ResourceS), в котором регулярно публикуются точные и своевременные прогнозы ожидаемой урожайности сельскохозяйственных культур по территории союза и других стратегических регионов мира. Технология системы AGRI4CAST, оперативно работающая с 1992г., в состоянии контролировать весь вегетационный период таких культур как, зерновые, масляничные, белковые, сахарной свеклы, картофеля, рис, пастбища, способна реагировать на эффекты краткосрочного метеорологического воздействия. Эффективность системы обеспечивают данные дистанционного зондирования и метеорологических наблюдений, агрометеорологическое моделирование и статистический анализ.
2. Описание данных, моделей и терминов
Модели включают данные стандартных агрометеорологических наблюдений (осадки, температуру и дефицит насыщения воздуха) и спутниковые данные дистанционного зондирования (вегетационный индекс NDVI). Данные по урожайности взяты из статистических сборников субъектов региона. Данные NDVI=(nir-red)/(nir+red), где nir - отражение в ближней инфракрасной области спектра (диапазон волн 0,7-1,0 мкм), red - отражение в красной области спектра (диапазон волн 0,6-0,7 мкм) из архива ИКИ РАН. Спектральные каналы выбирались из условия наименьшего отражения растительности в красном диапазоне спектра и, наоборот, самое высокое отражение по сравнению с другими природными объектами в ближнем ИК-диапазоне спектра. Падение спектральной кривой в красном диапазоне и резкий подъем в ближнем ИК-диапазоне характерно для растительности в хорошем состоянии. Фотосинтез высших растений приводит к максимуму отражения листьями в инфракрасной области и к максимальному поглощению хлорофиллом солнечной радиации в красной области спектра. Фазы вегетации яровых зерновых культур влияют на величину вегетационного индекса. Фазы кущения и выхода в трубку на 23 неделе делают максимальной скорость роста кривой функции вегетационного индекса. На 24 неделю приходится нижний узел соломины, на 27 неделю фазы колошения и цветения. Далее, 29 неделя, идут фазы молочной спелости, 32 неделя, восковой спелости и, 32 неделя, полной спелости и, наконец, на 37 неделе уборка урожая. На 29 неделю приходится максимум среднего многолетнего хода кривой ndvi. Максимальное поступление азота в растения яровой пшеницы (50 - 60 %), происходит между фазами кущения и колошения. В результате перед фазой коло -шения накапливается 70 - 80 % азота и формируют 50 - 60 % сухой массы растения. В этот период происходит формирование генеративных клеток, дифференциация органов колоса, в частности цветков в колосе. Потребление воды по фазам развития яровой пшеницы распределяется следующим образом. Макси-
мальное потребление воды (50-60 %) за весь вегетационный период происходит в фазе выхода растений в трубку и колошения, в фазе молочной спелости зерна 20-30 % и восковой спелости 3-5 %. При построении моделей прогноза урожайности фазы вегетации определяют выбор агрометеорологических параметров и интервалов их временного суммирования. Параметры регрессионных моделей находились с помощью метода наименьших квадратов, среда программирования МЛТЬЛБ. Как и в [5] оценка значимости коэффициента корреляции для парных моделей и в многомерном случае проводилась с помощью критерия Стьюдента, выборочное значение ^ которого, рассчитывались по формулам и сравнивались с *;кр (а; V) из таблицы, где уровень значимости а принимался равным 5%, а число степеней свободы v=N-m-1, где N=12 длина выборки, т число независимых переменных.
При оценке качества полученных моделей проверялось выполнение следующих условий:
1) Близость коэффициента детерминации г2 к значению 0.70, показывающую какой процент дисперсии исходного ряда урожайности у = {уь.. .,у№} описывает модель.
2) По критерию Стьюдента проверялась значимость коэффициентов регрессии. Критическое значение Ькр (а, V) определяется по уровню значимости а (принималось равным 5%) и числу степеней свободы V = N-2, где N - длина ряда.
3) Адекватность регрессионной модели исходным данным проверялась с помощью критерия Фишера Б.
4) Неравенство оу >а£, где оЕ - стандартное отклонение ошибок, оу - стандартное отклонение исходного ряда . Авторские испытания моделей проводились в соответствии с [6] на двухлетней контрольной выборке. По результатам этих проверок принималось решение о передаче построенных моделей на производственные испытания.
3. Построение моделей, расчеты
При построении моделей для одного из районов в качестве независимых параметров среди агрометеорологических данных выберем температуру (временной интервал I - 20У11, коэффициент корреляции г=-0,65) и дефицит насыщения воздуха (временной интервал 21.У1 - 20.У11, коэффициент корреляции г = -0,75). К ним добавим данные по вегетационному индексу иёу1.
Вначале построим однофакторную модель (независимый признак температура) и оценим ее качество. Проверку значимости г проведем с помощью критерия Стьюдента, выборочное значение 1;* которого, рассчитаем по формулам и сравним с *;кр( а, V), где уровень значимости а принимается равным 5 %, а число степеней свободы V = N-2, где N длина выборки. Так как
5.63 = 1* > = 2.20,
то коэффициент корреляции г значим, т.е. между двумя переменными существует статистически значимая обратная связь. Расчет коэффициентов модели дает следующую зависимость урожайности от дефицита насыщения
у = - 1.09х3+ 24.37.
Коэффициент детерминации модели равен 0.56, что указывает на недостаточность числа независимых переменных для описания ряда. При этом статистическую значимость или адекватность построенной модели подтверждает величина критерия Фишера ( 20.2 = Б > ^р = 4.8) и его уровень значимости 0.01 значительно меньший задаваемого при расчете коэффициентов модели 0.05. Как видим статистика Фишера получилась намного хуже чем для модели предыдущего района и это должно повлиять на среднюю относительную ошибку конечной многофакторной модели.
Далее усложним модель включением еще одного параметра. Расчет коэффициентов модели дает следующую зависимость урожайности от температуры и дефицита насыщения воздуха
у = - 0.50х2- 0.83х3 + 31.79.
Коэффициент детерминации модели равен 0.59, что ниже требуемой границей 0.70 и указывает на недостаточность числа независимых переменных. При этом статистическую значимость или адекватность построенной модели подтверждает величина критерия Фишера (7.2 = Б > Ркр = 4.1) и его уровень значимости 0.01 значительно меньший задаваемого при расчете коэффициентов модели 0.05. Как видим, включение еще одного параметра не улучшило статистику Фишера, она осталась такой же низкой как и в модели с одним параметром. Далее включим в модель еще один параметр интегральный индекс иёуь В этом случае зависимость урожайности от трех параметров температуры, дефицита и индекса иёу1 будет следующей
у = - 0.62х2 + 0.19х3 + 18.3х4 - 10.
Включение спутниковых данных по вегетационному индексу значительно повысило качество модели. Рассмотрим подробнее изменения. Коэффициент детерминации вырос до 0.82, а это доля дисперсии зависимой переменной, объясняемая рассматриваемой моделью зависимости, то есть объясняющими переменными. Статистическую значимость или адекватность построенной модели подтверждает величина критерия Фишера (18.6 = Б > Ркр= 3.9), причем уровень значимости в рассчитываемой статистике Фишера равен 0.0003. Статистика модели от трех независимых переменных стала значительно лучше по сравнению с моделью от двух переменных.
Проверим модель на двухлетней контрольной выборке. Прогноз (у/, х^} считается оправдавшимся, если его ошибка \у[ -у^ с контрольной выборкой
(yf, Xi} меньше 0,67g при заблаговременности не более двух месяцев (см. формулы [6]). Имеем для первого (2012г.) и второго (2013г.) года соответственно
\у( - У"! = |7.38 - 7.9| = 0.52 < 0.67g = 0.67*3.5 = 2.3 \у[ - ycj = \14.98 - 16.2\ = 1.22 < 0.67g= 0.67*3.2 = 2.2.
Тогда оправдываемо сть метода у (см. формулы [6]) для двух контроль -ных годов равна 100 %. И средняя относительная ошибка оправдавшихся прогнозов Р, считаемая по формулам [6] равна
1 0.52 1.22
Р = - (-100 +-100) = 6.4 %.
2 v13.6 13.6 '
Средняя относительная ошибка возросла на порядок (0.6% и 6.4%) по сравнению с результатом. Как видим, по уровню оценок статистик Фишера для одно и дух параметрических моделей можно, по крайней мере качественно, давать оценку относительной ошибки трехпараметрической модели. Для сравнения приведем еще два прогноза. Инерционный прогноз дает на эти два года 14.4 ц/га и 7.9 ц/га, не оправдывается и в результате получается средняя относительная ошибка 54 %. Такая ошибка объяснима с учетом аномального по агроклиматическим параметрам 2012 года. Прогноз по тренду (или климатический прогноз)
y = - 0.4143x + 846.333
дает на эти два года 12.8 ц/га и 12.3 ц/га, также не оправдывается и в результате получается средняя относительная ошибка 36 %.
Для другого района в качестве независимых параметров модели среди агрометеорологических данных выберем осадки (временной интервал I - 20VII, коэффициент корреляции г=0,78) и дефицит насыщения (временной интервал 11-20VII, коэффициент корреляции г=-0,61). К ним добавим данные по вегетационному индексу ndvi.
Вначале построим однофакторную модель (независимый признак осадки) и оценим ее качество. Проверку значимости г проведем с помощью критерия Стьюдента, выборочное значение t* которого, рассчитаем по формулам и сравним с tKp( a, v), где уровень значимости а принимается равным 5 %, а число степеней свободы v = N-2, где N длина выборки. Так как
6.7 = t* > tKp = 2.20,
то коэффициент корреляции г значим, т.е. между двумя переменными существует статистически значимая прямая связь. Расчет коэффициентов модели дает следующую зависимость урожайности от дефицита насыщения
y = 0.10x_2 + 8.32.
Коэффициент детерминации модели равен 0.61, что указывает на недостаточность числа независимых переменных для описания ряда. При этом статистическую значимость или адекватность построенной модели подтверждает величина критерия Фишера (17.5 = Б > ^р = 4.8) и его уровень значимости 0.001 значительно меньший задаваемого при расчете коэффициентов модели 0.05.
Далее усложним модель включением еще одного параметра. Расчет коэффициентов модели дает следующую зависимость урожайности от осадков и дефицита насыщения
Коэффициент детерминации модели равен 0.66, что ниже требуемой границей 0.70 и указывает на недостаточность числа независимых переменных. При этом статистическую значимость или адекватность построенной модели подтверждает величина критерия Фишера (9.6 = Б > Ркр = 4.1) и его уровень значимости 0.005 значительно меньший задаваемого при расчете коэффициентов модели 0.05. Далее включим в модель еще один параметр интегральный индекс иёу1. В этом случае зависимость урожайности от трех параметров температуры, дефицита насыщения и индекса иёу1 будет следующей
Включение спутниковых данных по вегетационному индексу значительно повысило качество модели. Рассмотрим подробнее изменения. Коэффициент детерминации вырос до 0.82, а это доля дисперсии зависимой переменной, объясняемая рассматриваемой моделью зависимости, то есть объясняющими переменными. Статистическую значимость или адекватность построенной модели подтверждает величина критерия Фишера (13.9 = Б > Ркр= 3.9), причем уровень значимости в рассчитываемой статистике Фишера равен 0.001.
Проверим модель на двухлетней контрольной выборке. Прогноз (у/, х^} считается оправдавшимся, если его ошибка \у[ -у^ с контрольной выборкой (у-% х{} меньше 0,67о при заблаговременности не более двух месяцев (см. формулы [13]). Имеем для первого (2012 года) и второго (2013года) года соответственно
Тогда оправдываемо сть метода у (см. формулы [6]) для двух контроль -ных годов равна 100 %. И средняя относительная ошибка оправдавшихся прогнозов Р, считаемая по формулам [6] равна
у = 0.08х2 - 0.21х3 + 10.98.
у = 0.02х2 - 0.045х3 + 16.15*4 - 19.05.
Ь/-У1С \ = \ 6.76 - 7.1\ = 0.34 < 0.67о = 0.67* 3.7 = 2.5 \у( -ус[\ = \ 15.85 - 17.5\ = 1.64 < 0.67о= 0.67* 3.7 = 2.5.
1 0.34
1.64
Р = " (
2 14.4
100 +
14.4
100) = 6.9 %.
Для сравнения приведем еще два прогноза. Инерционный прогноз дает на эти два года 1б.4 ц/га и 7.1 ц/га, не оправдывается и в результате получается средняя относительная ошибка б8 %. Такая ошибка объяснима с учетом аномального по агроклиматическим параметрам 2G12 года. Прогноз по тренду (или климатический прогноз)
y = G.G995x - 18б.143 дает на эти два года 14.G ц/га и 14.15 ц/га, также не оправдывается.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Савин И.Ю., Лупян Е.А., Барталев С. А. Оперативный спутниковый мониторинг состояния посевов сельскохозяйственных культур в России. Геоматика, 2G11, №2, стр.б9-7б.
2. Клещенко А.Д., Савицкая О.В. Оценка пространственно-временного распределения урожайности зерновых культур и стандартизованного индекса осадков (SPI) по спутниковой и наземной информации Тру-ды ГГО - СПб.: ООО «Д'АРТ», 2G14. - Вып. 571. - C. 147 -1б1.
3. Клещенко А.Д., Лебедева ВМ., Найдина Т.А., Савицкая О.В. Использование спутниковой информации Modis в оперативной агрометеорологии, Современные проблемы дистанционного зондирования Земли из космоса, 2G15. - Т.12. - С. 143 - 154.
4. SPOT Vegetation satellite analysis, MARS Bulletin., 2G11, Vol.19, No. 22 (Web Version).
5. Кононенко СМ., Старостина Т.В. Совместное использование спутниковой и наземной метеорологической информации для мониторинга посевов сельскохозяйственных культур // Интерэкспо ГЕО-Сибирь^П. XIII Mеждунар. науч. конгр. : Mеждунар. науч. конф. «Дистанционные методы зондирования Земли и фотограмметрия, мониторинг окружающей среды, геоэкология» : сб. материалов в 2 т. (Новосибирск, 17-21 апреля 2G17 г.). - Новосибирск : СГУГиТ, 2G17. Т. 2. - С. 28-33.
6. РД 52.27.284-91 Mетодические указания. Проведение производственных (оперативных) испытаний новых и усовершенствованных методов гидрометеорологических и гелиоге-офизических прогнозов.- M., Гидрометеоиздат, 1991. - С. 98-1G7.
О С. М. Кононенко, Т. В. Старостина, 2018
