CC BY
54
Югорский НИИ информационных систем ул. Мира, 151, Ханты-Мансийск, 628011, Россия
E-mail: bryksin@uriit.ru
применение адаптированной модели биопродуктивности epic и космоснимков MODIS для прогнозирования урожайности зерновых культур на территории западной сибири
Введение
За последний год проблема развития сельского хозяйства в России вышла на передний план, что связано с реализацией Национального проекта развития агропромышленного комплекса. При этом дотационные бюджетные и коммерческие вливания в отрасль возросли в несколько раз. В данной ситуации особенно актуальной становится задача поэтапной оценки и уточнения урожайности зерновых культур, как в разрезе отдельных хозяйств, так и на уровне административных районов и субъектов РФ. Заблаговременный прогноз урожайности используется как коммерческими (для оценки кредитоспособности сельхозпроизводителя в рамках кредитования под будущий урожай), так и государственными структурами (для контроля над нецелевым использованием дотационных средств). Несмотря на важность проблемы, для осуществления прогноза в основном все еще используются традиционные приемы оценки агроклиматических ресурсов и районирования территорий с помощью комплексных агроклиматических показателей, что связано с существенными материальными и временными затратами [Хворова, Брыксин, 2002]. Ситуация осложняется еще и тем, что в период недофинансирования произошло ухудшение экологического состояния сельскохозяйственных земель, вызванное существенными нарушениями агротехники возделывания, наблюдавшимися в 90-е гг. прошлого века, что приводит к уменьшению продуктивности сельскохозяйственных культур.
Практически оперативную оценку и мониторинг биопродуктивности зерновых культур с достаточно высокой точностью можно осуществлять используя шесть основных метеорологических характеристик (максимальная и минимальная температура воздуха, приходящая солнечная радиация, количество осадков, относительная влажность и скорость ветра), которые измеряются на российских метеостанциях, и космические сканерные изображения среднего разрешения.
Использование имитационных моделей биопродуктивности могло бы служить решением описанной проблемы, однако, настройка моделей и последующая корректировка по фактическим данным требует существенных материальных и временных затрат [Хворова, Брык-син, 2002]. Использование только данных дистанционного зондирования существенно снижает производственные затраты, но, как было показано в [Евтюшкин, Юшаков, 2001], также не дает достоверных результатов.
Выходом из сложившейся ситуации может служить применение имитационной модели биопродуктивности в совокупности с данными дистанционного зондирования Земли (ДЗЗ), что позволяет осуществлять оценку урожайности без проведения дорогостоящих полевых измерений. Для данных целей авторами предлагается использование модели биопродуктивности EPIC, разработанной J. R. Williams в сельскохозяйственной службе США [Williams. 1997], и космоснимки сканера MODIS со спутников TERRA и AQUA, принимаемых в Центре ДЗЗ Югорского НИИ информационных технологий (ЮНИИ ИТ). Сканер MODIS имеет
ISSN 1818-7900. Вестник НГУ. Серия: Информационные технологии. 2007. Том 5, выпуск 2 © В. М. Брыксин, 2007
полосу захвата 2330 км, 36 спектральных каналов в интервале от 0,4 до 14,38 мкм. Пространственное разрешение в каналах 1 и 2 равно 250 м, в каналах с 3 по 7 - 500 метров, в остальных - 1000 м.
По орбитальным изображениям MODIS вычислялся нормализованный вегетационный индекс (NDVI), позволяющий уменьшить влияние условий съемки в разное время суток:
NDVI = (L2 - Ll)/(L2 + Ll), где L1 и L2 - яркости 1 и 2 каналов сканера MODIS. Использование сканера MODIS, позволяет отслеживать динамику развития культур на уровне региона до четырех раз в сутки, что повышает частоту наблюдений при разреженной динамической облачности.
Ранее в результате подспутниковых экспериментов [Евтюшкин, Рычкова, 2004] было показано, что значения NDVI хорошо коррелируют с такими характеристиками растительного покрова, как индекс площади листьев (LAI), коэффициент вегетационного покрытия и зеленая биомасса (чистая первичная продуктивность). С увеличением значений перечисленных параметров растут и значения NDVI, но зависимость эта линейна только в определенном интервале. Уже при содержании хлорофилла выше 3-5 мг/см2 и листовом индексе выше 4-5 наблюдается эффект насыщения NDVI.
Для перехода от индекса NDVI, рассчитанного по данным радиометра MODIS, к листовому индексу (LAI), который является одним из наиболее существенных параметров в модели EPIC, используется кусочно-линейная зависимость, предложенная авторами алгоритма MOD15 (NASA).
Для прогноза урожайности используется модель биопродуктивности MIDC, которая представляет собой адаптацию модели развития растений EPIC, разработанной в США, к условиям Западной Сибири и Урала. EPIC [Williams, 1997] является непрерывной имитационной моделью блочного типа, которая может использоваться для определения стратегий управления в сельском хозяйстве. Данная модель позволяет рассчитывать основные параметры культуры (биомассу, листовой индекс, высоту растения и др.) на каждый день, что позволяет сравнивать их с данными наземных измерений и дистанционного зондирования и в случае необходимости производить корректировку. В модели EPIC фенологическое развитие культуры основывается на сумме эффективных температур. Потенциальное ежедневное изменение в сухой биомассе (ЛБР) зависит от листового индекса (LAI) и солнечной радиации (КА):
ЛБР = 0,0005 • BE • RA • (l - e^0'65 ’LA), где BE - коэффициент преобразования энергии к биомассе.
Реальное ежедневное приращение биомассы корректируется коэффициентом стресса (REG), который определяется как минимальная величина из водного, температурного и азотного стрессов. Моделирование приращения листового индекса (LAI) разделено на два этапа: от появления листьев до начала увядания и от начала увядания до конца сельскохозяйственного сезона. Урожайность (YLD) рассчитывается пропорционально итоговой сухой биомассе (BIOM): YLD = HI • max(BIOM), где HI - индекс урожайности.
В ходе тестирования модели были выявлены некоторые существенные недостатки [Хво-рова и др., 2001]. В частности, в модели присутствует существенное число метеорологических и почвенных характеристик, которые в России не измеряются. Блоки радиации, температуры почвы и азотного баланса реагируют не адекватно. С другой стороны, к основным можно отнести всего шесть метеорологических характеристик (максимальная и минимальная температура воздуха, приходящая солнечная радиация, количество осадков, относительная влажность и скорость ветра), которые полностью или частично измеряются на российских метеостанциях. Из выше сказанного можно сделать вывод о целесообразности применения модели EPIC после проведения ее адаптации к условиям Западной Сибири.
При анализе данных, предоставленных Западно-Сибирским Управлением гидрометс-лужбы в рамках договора о научном сотрудничестве, было выявлено, что приходящая солнечная радиация не измеряется на большинстве метеостанций, а блок радиации в модели EPIC позволяет найти только максимальную суммарную приходящую радиацию. Поэтому было принято решение заменить блок радиации в модели EPIC. Для расчета суточного хода приходящей солнечной радиации предлагается воспользоваться формулой Т. Г. Берлянда, используемой в модели радиационного режима АФИ РАСХН [Полуэктов, 1991]:
п = .^т к,
1 + / • coseck/
где , = 1,98 кал • см2 • мин1 - солнечная постоянная, / - коэффициент, зависящий от географической зоны и времени года, к, - высота Солнца.
Суммарная интенсивность радиации на уровне Земли зависит от высоты Солнца, прозрачности атмосферы, степени облачности, продолжительности солнечного сияния. Так как в предоставляемых гидрометслужбой данных отсутствуют продолжительность солнечного сияния и облачность, то для учета этих факторов предлагается использовать расчет приходящей солнечной радиации (0 через относительную влажность воздуха:
Q = — 1,318 • ~ • Qo + 0,922 • ~ + 1,357 • Qo + 1,012, для мая,
Q = — 1,448 • ~ • Qo — 0,035 • ~ + 1,682 • Qo + 0,077, для июня-августа, где ~ - относительная влажность воздуха. При расчете предполагалось, что даже при наличии сплошной облачности определенный процент радиации достигает поверхности почвы. Ниже приведена оценка точности расчетов приходящей радиации:
средняя абсолютная процентная ошибка
Модель / город Барнаул Тобольск Коченево
Расчет по EPIC, % 98,66 97,77 111,28
АФИ, % 22,18 25,51 27,86
На рис. 1 приведен пример динамики расчетной и экспериментальной суммарной солнечной радиации.
30
Q.
МДж/м2
25
20
15
10
0
13 14 15 16 17 18 19 20 21 22 23 24
ШЗ Эксперементальная ♦ Расчетная Номер декады
Рис. 1. Изменения расчетной и экспериментальной суммарной солнечной радиации Q (МДж/(м2-сут.))
за вегетационный период в 2004 г. (Барнаул).
В ходе анализа блока температурного стресса модели EPIC было выявлено, что он не адекватно реагирует на изменение температуры почвы (заморозки и перегревание) [Хворова и др., 2001]. Таким образом, возникла проблема разработки блока расчета температуры почвы, зависящей от минимального числа параметров.
Из литературы известно, что днем поверхность почвы нагревается прямой солнечной радиацией, а ночью поверхность охлаждается в результате эффективного излучения (разность восходящего и нисходящего потоков длинноволновой радиации на уровне Земли). Причем минимальная температура воздуха наблюдается перед восходом Солнца, а максимальная отмечается через 2-3 часа после полудня. По экспериментальным данным была установлена зависимость между максимальной температурой почвы, максимальной температурой воздуха и суммарной приходящей солнечной радиацией, наблюдаемыми в течение дня, которую предлагается использовать для расчетов:
Tpmax = l, 2784 • TV,max + 0,4499 • Q + 0,4744,
где Тр тш - максимальная температура на поверхности почвы, Тушх - максимальная температура воздуха. Погрешность вычислений составляет 10,75 %.
В ночное время наблюдается взаимосвязь между минимальной температурой почвы и эффективным излучением. Исходя из данного соображения получено соотношение, которое предлагается использовать для расчетов минимальной температуры почвы:
Тршп = 0, 8024 • Т,шп + 1,4789 • Е + 4,6385, где Т,,тт - минимальная температура на поверхности почвы, Ту Шп - минимальная температура воздуха, Е - эффективное излучение. Для расчета эффективного излучения использовалась двухпотоковая модель длинноволновой радиации [Матвеев, 1984]. Погрешность вычислений составляет 9,76 %.
На рис. 2 приведен пример динамики расчетной и экспериментальной максимальной температуры почвы для метеостанции Барнаул.
Эксперементалычая » Расчетная
Рис. 2. Изменения расчетной и экспериментальной максимальной температуры почвы (°С) в течение вегетационного периода в 2003 г. (Барнаул).
Теплота, поступающая на поверхность почвы, под действием создаваемого градиента температур перераспределяется в толщине почвенного слоя. Для расчета температуры почвы на глубине предлагается использовать уравнение теплопроводности с граничными условиями [Пивоварова, 1977]. Удельная теплоемкость почвы определяется как средневзвешенная удельная теплоемкость воды и почвенного скелета. Нижняя граница помещается, как правило, на глубине, на которой температура либо постоянна, либо зависит от времени известным образом. В качестве верхнего граничного условия записывается соотношение, обеспечивающее «сшивание» решений задачи в почве и в приземном воздухе. Полученные в результате расчетов данные по распределению температуры в почве хорошо согласуются с имеющимися экспериментальными данными и данными из литературных источников. Погрешность составляет 9,8 %.
Основываясь на соотношениях, описанных выше, были проведены корректировка и настройка модели биопродуктивности EPIC для Новосибирской, Тюменской областей и Алтайского края на основе метеоданных и данных о фактической урожайности за 1985-2005 г.
В 2004 г. на конец первой декады июля выполнено прогнозирование урожайности и сроков созревания яровой пшеницы по Тюменской области, а в 2005-2006 гг. по Новосибирской области и Алтайскому краю - по адаптированной модели биопродуктивности EPIC. Прогнозы на 11 июля и 31 августа передавались в Департамент АПК администрации Тюменской области 2004 г. и в Департамент АПК администрации Новосибирской области совместно с ЗапСИБРЦПОД в 2005-2006 гг. Прогнозируемая урожайность для Тюменской области на 1 июля 2004 г. составляет 19 ц/га и на 31 августа 19,4 ц/га, фактическая урожайность по данным статистической отчетности составляет 19,6 ц/га по данным на 31 августа 2004 г. Прогнозируемая биологическая урожайность на 31 августа 2005 г. для Новосибирской области
составила 14,7 ц/га, для Алтайского края 12,6 ц/га, и на 11 июля 2006 г. 8,9 ц/га и 9,6 ц/га соответственно. Динамика прогнозируемой и фактической урожайности яровой пшеницы по Новосибирской области в 2005 году представлена на рис. 3.
13.50
2
3
Л13.00
н
8
X
>5
І1250
12.00
11.50
11.00
1 12.70 12.93
іК^Зб 11.93 11.93
1 11.68 1 ■ I ■
1 Ті
10.50
10.07.2005 20.07.2005 31.07.2005 10.08.2005 20.08.2005 31.08.2005 10.09.2005
■Прогноз
■ Факт
Рис. 3. Динамика прогнозируемой и фактической урожайности яровой пшеницы по Новосибирской области в 2005 г.
Поскольку в регионах Западной Сибири имеется недостаточная густота сети метеостанций, проводилась коррекция листового индекса с использованием зависимости LAI от NDVI по данным сканера TERRA/MODIS. Пример расчетов основных биометрических показателей без коррекции (урожайность 14,06 ц/га) и с коррекцией листового индекса (13,26 ц/га) для Павловского района Алтайского края в 2005 г. представлен на рис. 4.
На рис. 5 представлена динамика прогнозируемой урожайности в зависимости от фактических метеоусловий. Из рисунка видно, что колебания в расчетных значениях биологической урожайности яровой пшеницы становятся минимальными в конце июля. Это подтверждает возможность прогноза урожайности яровой пшеницы на юге Западной Сибири за месяц до начала уборки.
Таким образом, интеграция наземных измерений биометрических параметров растительности, ежедневных метеоданных, значений статистической урожайности и сканерных космических изображений в единую систему позволяет производить экспорт данных в модель биопродуктивности и последующий анализ результатов. Использование методологии моделирования биопродуктивности по модели EPIC совместно с данными дистанционного зондирования повышает точность расчета урожайности основных сельскохозяйственных культур. Выполняются ежедневные расчеты биопродуктивности яровых зерновых культур по модифицированной методике EPIC с корректировкой расчетных значений листового индекса (LAI) на основе спутниковой оценки вегетационного индекса. Периодичность коррекции - не менее двух раз неделю с учетом разреженной облачности на период до начала спада листового индекса, связанного с естественным увяданием листового покрова зерновых культур. На основании расчетов по модели прогнозируются биологическая урожайность культур и оптимальные сроки начала уборки в разрезе районов или субъекта федерации.
На основании метеоданных и космического мониторинга по модифицированной модели рассчитывают основные параметры культуры (биомассу, листовой индекс, вес корня, высоту растений) и прогнозируют сроки биологического созревания зерновых культур и оптимальные сроки начала уборки. Прогноз урожайности для различных зерновых культур выдается за месяц до начала уборки с точностью до 1-2 ц/га.
Применение данной методики перспективно для оценки урожайности и валовых сборов на территории Поволжского, Уральского и Сибирских федеральных округов в интересах орга-
Графики развития Биомассы (ВЮМ), Листового индекса (IAI) и веса кормя (FW) растения
Графики развития Биомассы (БЮМ), Листового индекса (LAI) и беса корня (RW) растения
- LAI
810М
RW
а
б
Рис. 4. Результаты расчетов для Павловского района Алтайского края без коррекции (а) и с коррекцией (б) листового индекса
нов государственной власти, банковских структур, предоставляющих кредиты, и страховых компаний, возмещающих ущерб сельхозпроизводителям в связи с потерями урожая.
список литературы
Евтюшкин А. В., Юшаков В. Н. Разработка элементов технологии мониторинга зерновых культур в Алтайском крае // Тр. междунар. конф. «Современные проблемы информационных технологий и космический мониторинг». Ханты-Мансийск, 2001. С. 154-158.
Евтюшкин А. В., Рычкова Н. В. Изучение динамики вегетационного индекса сельскохозяйственной зоны Западной Сибири по данным NOAA и подспутниковых полигонных измерений // Информационные технологии и космический мониторинг: Материалы III науч.-практ. конф. «Электронная Россия». Ханты-Мансийск, 2004. Екатеринбург: Изд-во «Баско», 2004. С. 54-61.
Матвеев Л. Т. Курс общей метеорологии. Физика атмосферы. Л.: Гидрометеоиздат, 1984.
Пивоварова З. И. Радиационные характеристики климата СССР. Л.: Гидрометеоиздат, 1977.
Полуэктов Р. А. Динамические модели агроэкосистемы. Л.: Гидрометеоиздат, 1991.
Хворова Л. А., Брыксин В. М., Масленников Д. В. Система имитационного моделирования агроэкосистем // Вестник алтайской науки «Проблемы агропромышленного комплекса». Барнаул: АГУ, 2001. Т. 1. С. 277-279.
АЗ 8
"s
х
о 6
X
S
Q. 4 >.
2
О
20. С
\
S 1
1 ■ 1 .11. 1 1 lh
- 4,5 4
3.5 3
2.5 '
2 ,
1.5
0
i.2005 30.06.2005 10.07.2005 20.07.2005 30.07.2005 09.08.2005 19.08.2005 29.08.2005 08.09.2005 Дата прогноза
I ^Суммарные осадки Листовой индекс Прогнозируемая урожайность
Рис. 5. Динамика прогнозируемой урожайности в сравнении с листовым индексом и осадками для метеостанции Здвинск Новосибирской области
Хворова Л. А., Брыксин В. М. Применение математических методов и математического моделирования для оценки агроклиматического потенциала территорий // Изв. Алтайского гос. ун-та. 2002. №1 (23). С. 41-45.
Williams J. R. The Erosion-Productivity Impact Calculator (EPIC). Technical Reference. US Department of Agriculture. 1997.
Материал поступил в редколлегию 12.09.2007
