Индикация многолетних изменений в растительном покрове залежных земель лесостепи на основе рядов вегетационного индекса NDVI

Определение тенденций в развитии растительного покрова является одной из ключевых задач в исследованиях состояния природной среды. Участки пашни, выведенные из сельскохозяйственного оборота, или залежные земли, являются одним из немногих типов объектов в районах интенсивного аграрного использования, на основе которых эти тенденции могут быть изучены.

На оставленных аграрных угодьях после прекращения сельскохозяйственной деятельности начинаются процессы восстановления почвенно-растительного покрова, протекание которых во многом определяется текущими природно-климатическими условиями в регионе [1, 2].

Изменения растительного покрова, наблюдающиеся на залежах, находят отражение в их спектральном отклике [3, 4]. В связи с этим становится актуальным количественный анализ спектрально-отражательных свойств залежных земель для пространственновременного моделирования изменений в их растительности [5, 6]. Учитывая, что развитие растительного покрова носит непрерывный характер на протяжении вегетационного сезона каждого года, одним из условий его объективного анализа является изучение многолетних рядов спектрально-отражательных характеристик [7–9]. Решение данной задачи во многом связано с использованием разновременных данных дистанционного зондирования, которые благодаря периодичности получения и одновременному охвату значительных территорий часто выступают практически единственным источником информации для реконструкции изменений в состоянии растительности [10, 11].

Информационные продукты MOD13Q1[12] являются одним из современных типов спутниковых данных. Они создаются на основе регулярно получаемых космических снимков MODIS [13] и на текущий момент являются единственным видом спутниковой информации, позволяющим построить и проанализировать для конкретных угодий непрерывные многолетние ряды вегетационных индексов, с начала 2000-х гг. по настоящее время. Данные MOD13Q1 содержат 16дневные композитные изображения вегетационного индекса NDVI, рассчитанные на основе атмосферно и радиометрически откорректированных спутниковых снимков, что является ключевым условием объективного сравнения спектральных признаков, измеренных в разные временные срезы. Основным ограничением применения данных MOD13Q1 является их относительно невысокое пространственное разрешение, со- ставляющее 250 м / пиксель, вследствие которого необходимо учитывать минимальную площадь анализируемых объектов.

Вегетационный индекс NDVI, основанный на коэффициентах спектральной яркости в красной и ближней инфракрасной областях спектра [14], является показателем, достаточно широко используемым для анализа состояния растительности [15, 16]. Тем не менее, вопросы его применения для моделирования сукцессионных процессов, т.е. последовательных смен одних растительных сообществ другими, на оставленных аграрных угодьях во многом остаются открытыми. Это обусловлено региональными особенностями протекания сукцессий на залежах вследствие определенного сочетания природноклиматических условий и недостаточной изученностью их спектрального отклика. Особый интерес представляет анализ изменений растительности на основе рядов многолетних значений вегетационного индекса. Это связано с тем, что относительно длинные ряды NDVI, охватывающие интервал порядка двух десятилетий, накоплены только к настоящему времени, т.к. сенсор MODIS поставляет данные с 2000 года.

Центральное Черноземье, расположенное преимущественно в зоне лесостепи, относится к регионам активного аграрного использования. Оценка изменений растительности на залежных землях, типичных для его территории, является актуальной задачей для получения объективных представлений о региональных трендах развития растительного покрова и обусловленных ими тенденциях в динамике спектрально-отражательных признаков.

Предыдущими исследованиями [10] установлено, что одной из ключевых особенностей современного развития растительного покрова региона является формирование древесных насаждений на травянистых пространствах. Оценка скорости данного процесса на залежах и его пространственных особенностей является задачей, решение которой возможно исключительно на основе разновременной спутниковой информации.

Цель исследования, таким образом, заключалась в количественном анализе спектрально-отражательных свойств оставленных аграрных угодий лесостепи на основе информационных продуктов MOD13Q1 и оценке их возможностей для определения пространственно-временных тенденций динамики растительности залежных земель.

Задачами исследования являлись:

• 1)    оценка информационных возможностей изображений MOD13Q1 для анализа степени покрытия залежей древесными сообществами;

• 2)    расчет и анализ многолетних рядов вегетационного индекса NDVI для залежей с различными типами формирующихся на них лесных насаждений;

• 3)    пространственно-временное моделирование изменений в растительности залежей на основе

многолетних рядов спектрального вегетационного индекса.

1.    Материалы и методы

Исследование проведено на основе экспериментальных данных, собранных на территории Центрального Черноземья (ЦЧР), включая шесть областей: Белгородскую, Курскую, Воронежскую, Липецкую, Тамбовскую и Орловскую.

Работа была выполнена в несколько этапов:

• -    формирование и предварительная обработка многолетней серии изображений MOD13Q1 в период 2000–2018 гг., покрывающих территорию исследования;

• -    формирование экспериментальной выборки объектов и создание их векторной основы;

• -    вычисление многолетнего ряда NDVI для каждой отобранной залежи;

• -    математический анализ зависимости между фактической величиной покрытия залежей древесными сообществами, формирующимися на них, и значениями вегетационного индекса NDVI, измеренными на основе данных MOD13Q1;

• -    анализ параметров многолетних рядов NDVI для залежей Центрального Черноземья с различными типами формирующихся на них древесных сообществ: лиственных, хвойных, смешанных;

• -    геоинформационное моделирование пространственно-временных особенностей в значениях вегетационного индекса залежей в пределах лесостепи ЦЧР в период 2000–2018 гг.

Сформированная многолетняя серия композитных изображений MOD13Q1 (источник данных: охватывала период съемки с начала функционирования сенсора MODIS: 2000–2018 гг. На каждый год отбирали данные на 13 временных срезов, полностью охватывающих вегетационный период с начала апреля по конец октября с интервалом в 16 дней. В общей сложности было обработано 247 изображений, каждое из которых включало территорию лесостепи Центрального Черноземья.

Предварительная обработка изображений MOD13Q1 заключалась в вырезании из них фрагментов, покрывающих территорию исследования, и перевод в проекцию UTM. После этого из них были извлечены растры значений NDVI на каждую исследуемую дату.

Для анализа была сформирована экспериментальная выборка объектов, включающая несколько сотен оставленных аграрных угодий, расположенных на изучаемой территории. Объекты отбирали с использованием снимков сверхвысокого пространственного разрешения (1 м / пиксель) на основе ряда критериев:

• -    угодье должно быть залежью на конечную анализируемую дату (2018 г.) и использоваться под пашню в конце 1990 – начале 2000 гг.;

• -    минимальная площадь угодья - 20 га, определена с учетом пространственного разрешения

спутниковых данных MOD13Q1, по которым считали спектральные признаки. При этом средняя площадь составила в итоге более 63 га;

• -    в выборке объектов должны быть представлены залежи со всем типами древесных сообществ, формирующихся на них: лиственными, хвойными, смешанными.

Для анализа стремились отбирать залежи, имеющие геометрически максимально правильную форму, чтобы большинство пикселей располагалось внутри их контуров. При этом отбор таких залежей производился с учетом визуального анализа наложения их контуров на спутниковые данные.

В общей сложности было отобрано и изучено 486 оставленных аграрных угодий (табл. 1). С использованием геоинформационных средств методом ручной оцифровки для них подготовили векторную основу, к которой на следующем этапе была присоединена информация о состоянии растительного покрова: типе древесных сообществ, величине покрытия залежей древесной растительностью в 2018 г. и в другие годы исследуемого периода, на которые удалось провести ее оценку.

Табл. 1. Характеристики проанализированных залежных земель

Тип лесных насаждений на залежах	Площадь средняя, га	Число	Площадь общая, га
Лиственный	65,1	297	19344,6
Хвойный	47,5	70	3327,7
Смешанный	68,8	119	8181,4
Всего	63,5	486	30853,6

Для оценки состояния залежей на этом этапе были использованы следующие типы спутниковых данных: снимки сверхвысокого пространственного разрешения из Google Earth, многозональные снимки Landsat TM, ETM+ начала 2000 гг., современные снимки Landsat OLI.

Состояние залежей на конечную анализируемую дату, 2018 г., определялось на основе спутниковых данных сверхвысокого пространственного разрешения (1 м / пиксель), полученных из интернет-сервисов открытого доступа, преимущественно Google Earth.

Соответствующее разрешение снимков по ряду дешифровочных признаков позволяет достоверно выявлять залежные земли среди пашни и оценивать необходимые для исследования параметры их растительного покрова. Состояние угодий на рубеже 1990–2000 гг. определялось на основе многозональных снимков Landsat TM, ETM+, полученных в период вегетации.

Тип лесных насаждений идентифицировался на основе снимков Landsat OLI 2018 г. (пространственное разрешение 30 м / пиксель). Его определение производили в комбинациях инфракрасных (2,10– 2,30 мкм, 1,56 – 1,66 мкм) и красного (0,64 – 0,68 мкм)

каналов, в которых лиственные и хвойные породы могут быть достоверно разделены. По снимкам сверхвысокого пространственного разрешения 2018 года для каждой залежи определили величину проективного покрытия кронами древесных сообществ. При наличии детальных снимков на более ранние даты для них также определяли величину проективного покрытия крон. Его оценка осуществлялась методом анализа снимков на каждый исследуемый объект (примеры на рис. 1).

а)

покрытия кронами древесных сообществ:

0,5 (а); 0,07 (б)

б)

Рис. 1. Примеры залежных земель с различной величиной

Следующий этап состоял в вычислении для залежей значений вегетационного индекса на основе изображений MOD13Q1. Для всех оставленных аграрных угодий на их основе методом зональной статистики рассчитали многолетний ряд NDVI, включающий значения индекса за вегетационный период (начало апреля – конец октября) каждого года с 2000 по 2018 гг.

Таким образом, многолетний ряд индекса NDVI представлял последовательность значений, каждое из которых являлось его величиной, усредненной за вегетационный сезон:

(Уп) = У1, У 2, У3..........Уп где yn – среднее значение NDVI за вегетационный период отдельных лет, n – число анализируемых лет, равное в нашем случае 19, с 2000 по 2018 гг.

Метод зональной статистики позволяет вычислять статистические параметры для конкретных признаков (в нашем случае для NDVI), рассчитанные в пределах определенных границ. Для оставленных аграрных угодий были рассчитаны значения NDVI, усредненные в границах их контуров.

Таким образом, была сформирована необходимая выборка, включающая данные о параметрах растительного покрова и соответствующих им значениях вегетационного индекса. Наличие значений NDVI, измеренных в разные сроки вегетационного периода, позволило оценить информативность сезонных значений индекса для анализа величин проективного покрытия залежей кронами деревьев.

Форма и параметры зависимости NDVI от проективного покрытия крон были оценены с использованием регрессионного и графического анализа в программе STATISTICA. Для анализа были использованы величины проективного покрытия крон и значения вегетационного индекса, измеренные в 2018 г.

Для залежей с каждым анализируемым типом лесных насаждений (лиственные, хвойные, смешанные) вычислили усредненный со всех угодий многолетний ряд индекса по значениям, измеренным за вегетационный сезон. Оценка наличия и значимости трендовой составляющей была выполнена в программе R с использованием непараметрического теста Манна–Кендалла.

В исследовании был предложен новый подход для оценки скорости зарастания залежей древесными сообществами. Он основан на пространственно-временном анализе коэффициента наклона линии тренда NDVI, который может выступать индикатором скорости формирования древесных сообществ на залежах.

С использованием геоинформационных пространственно-статистических инструментов в программе ArcGIS вычислили картосхему, характеризующую территориальное изменение коэффициента наклона трендовой линии NDVI в пределах Центрального Черноземья. Картосхему рассчитали и проанализировали для залежей с лиственными древесными сообществами. Это наиболее распространенный тип залежей, репрезентативно представленный во всех его частях, что являлось необходимым условием корректного пространственного анализа. На основе полученной картосхемы удалось установить ряд особенностей динамики растительного покрова оставленных аграрных угодий в период 2000–2018 гг.

2.    Результаты

Оценка информационных возможностей изображений MOD13Q1 показала (табл. 2, 3), что значения вегетационного индекса NDVI, вычисленные на их основе, достаточно тесно коррелируют с величиной проективного покрытия залежей кронами древесных сообществ.

Табл. 2. Параметры корреляции между проективным покрытием крон и значениями NDVI для залежей с лиственными породами

Коэффициент корреляции	NDVI
	13.06– 26.06	27.06– 12.07	13.07 – 28.07	29.07 – 13.08	13.08 – 28.08
Пирсона	0,57	0,59	0,39	0,67	0,69
Спирмена	0,59	0,63	0,46	0,71	0,74

При этом статистически значимая корреляция отмечена во все сроки вегетационного периода. Для залежей с хвойными насаждениями характерны более высокие коэффициенты корреляции, чем с лиственными породами.

Табл. 3. Параметры тесноты связи между проективным покрытием крон и значениями NDVI для залежей с хвойными породами

Коэффициент корреляции	NDVI
	13.06 – 26.06	27.06 – 12.07	13.07 – 28.07	29.07– 13.08	13.08 – 28.08
Пирсона	0,73	0,80	0,55	0,80	0,82
Спирмена	0,70	0,79	0,56	0,79	0,82

Другой важной особенностью является то, что для залежей с насаждениями лиственного породного состава коэффициент корреляции Спирмена выше, чем значения коэффициента корреляции Пирсона, что предварительно указывает на нелинейный характер связи между величиной покрытия залежей древесными сообществами и значениями вегетационного индекса.

Из табл. 2, 3 также следует, что наиболее высокими коэффициентами корреляции с величиной покрытия залежей кронами древесных сообществ характеризуются значения вегетационного индекса периода середины – второй половины августа (13–28 августа). Т.е. значения NDVI этих сроков наиболее информативны для пространственно-временного анализа доли древесных насаждений, присутствующих на залежах лесостепи Центрального Черноземья.

Регрессионный анализ подтвердил нелинейную форму зависимости между значениями NDVI и величиной покрытия крон для залежей с формирующимися лиственными насаждениями (рис. 2).

Рис. 2. Изменение значений NDVI в зависимости от величины проективного покрытия крон для залежей с лиственными древесными сообществами

Она может быть описана статистически значимой логарифмической зависимостью (уровень значимости меньше 0,05).

Этот же вид кривой может быть использован и для описания связи между исследуемыми параметрами для оставленных аграрных угодий с хвойными лесными насаждениями (рис. 3). Для описания зависимости в обоих случаях были использованы значения вегетационного индекса, измеренные в период середины августа.

Рис. 3. Изменение значений NDVI в зависимости от величины проективного покрытия крон для залежей с хвойными древесными сообществами

В то же время при более детальном анализе графиков NDVI, а также учитывая, что значения проективного покрытия имеют предельную величину, зависимость между ним и вегетационным индексом должна приближаться к S-образной, сигмоидной кривой.

Высокая теснота связи между величиной покрытия залежей древесными сообществами и значениями вегетационного индекса выступила теоретической основой для оценки тенденций изменения растительного покрова залежей в период 2000–2018 гг. на основе значений NDVI.

Оценка многолетних рядов индекса, усредненных на основе данных со всех залежей одного и того же типа, показала, что для оставленных аграрных угодий лесостепи Центрального Черноземья в период 2000– 2018 гг. наблюдалась положительная тенденция динамики индекса, статистически значимая на уровне меньше 0,05 (рис. 4–6, табл. 4).

Годы

Рис. 4. Многолетняя динамика вегетационного индекса для залежей с формирующимися насаждениями лиственного породного состава

Она характерна для залежей с лиственными, хвойными и смешанными лесными насаждениями.

Значимость трендовой составляющей, оцененная на основе критерия Манна–Кендалла, показала положительные значения статистики tau для всех типов оставленных аграрных угодий (табл. 4). Статистически значимая величина показателя tau указывает на наличие положительного тренда, отрицательная – является индикатором тенденции снижения.

Учитывая, что полученные значения вычислены на основе данных, собранных в разных частях региона, можно сделать выводы о количественных параметрах тенденций динамики растительности, наблюдающихся на залежах лесостепи с начала XXI века по настоящее время.

Табл. 4. Параметры статистики tau для разных типов залежей лесостепи Центрального Черноземья

Тип древесных сообществ на залежах	Значение tau	Уровень значимости
Лиственный	0,09	0,03
Хвойный	0,25	0,00
Смешанный	0,14	0,00

Из рис. 4–6 видно, что значения индекса возрастали на протяжении анализируемого периода, но с разной интенсивностью.

Рис. 5. Многолетняя динамика вегетационного индекса для залежей с формирующимися насаждениями хвойного породного состава

Рис. 6. Многолетняя динамика вегетационного индекса для залежей с формирующимися насаждениями смешанного породного состава

Снижение значений вегетационного индекса в отдельные годы (2002, 2009, 2010) можно объяснить климатическими факторами. В частности, 2010 год в регионе характеризовался наиболее жарким и сухим вегетационным периодом. На промежуточной стадии формирования древесных сообществ на залежах, когда древесной растительностью покрыта только их часть, а остальная занята травянистым пространством, фактор сухой травы мог являться причиной снижения NDVI.

Учитывая, что вегетационный индекс достаточно тесно коррелирует с долей древесной растительности, присутствующей на залежах, многолетняя динамика NDVI позволила оценить интенсивность формирования древесных сообществ в период 2000–2018 гг. Для ее анализа предложено использовать величину коэффициента наклона трендовой линии, вычисленной на основе многолетнего ряда вегетационного индекса. Вследствие прямой зависимости между NDVI и проективным покрытием залежей кронами, более высокие значения коэффициента трендовой линии должны соответствовать и более высокой интенсивности зарастания залежей лесными насаждениями.

Исходя из примерно одинакового возраста изученных залежей и графиков динамики NDVI (рис. 4–

6), можно сделать вывод, что залежи с хвойными породами характеризуются наиболее высокой скоростью формирования древесных сообществ.

Высокая корреляция вегетационного индекса NDVI с величиной покрытия залежей кронами лесных насаждений, формирующихся на них, позволила провести пространственно-временную оценку скорости зарастания залежей региона лесной растительностью. Она была осуществлена для наиболее крупной и типичной для региона группы залежей с лиственными древесными сообществами.

Пространственно-временной анализ скорости формирования лесных насаждений на залежах предложено выполнить на основе величины коэффициента наклона трендовой линии в период 2000–2018 гг. Этот показатель вычислили для каждой из 297 изученных залежей лесостепи ЦЧР.

Интерполяция коэффициента наклона трендовой линии, проведенная в геоинформационной системе ArcGIS с использованием радиальной базисной функции, позволила получить пространственновременную картину данного показателя (рис. 7), которая характеризует скорость зарастания оставленных аграрных угодий древесными сообществами.

Рис. 7. Пространственное распределение коэффициентов наклона трендовой линии NDVI (умножен на 10000), вычисленных для залежей лесостепи Центрального Черноземья в период 2000–2018 гг.

На основе полученной картосхемы выявлен ряд особенностей, отражающих динамику растительного покрова залежей в начале XXI века.

Во-первых, практически на всей территории лесостепи Центрального Черноземья для залежей наблюдается положительная динамика вегетационного индекса, что указывает на аналогичную тенденцию зарастания залежей древесными сообществами. Ее отсутствие характерно только для южной части ЦЧР, относящейся преимущественно к подзоне южной лесостепи и северной степи.

Во-вторых, наблюдается выраженная территориальная тенденция в изменении коэффициента наклона трендовой линии NDVI с северо-запада на юго-восток и северо-востока на юго-запад. Соответствующие изменения наблюдаются при переходе от условий северной лесостепи к степной зоне [17]. Т.е. значения NDVI позволили установить, что скорость зарастания залежных земель лесной растительностью имеет выраженный внутризональный аспект.

Основной причиной в различной скорости роста древесных сообществ на залежах могут выступать природно-климатические условия, в большей степени на севере и меньшей степени на юге способствующие формированию лесных насаждений. В то же время этот вопрос требует отдельного исследования.

В-третьих, вариация в величине коэффициента наклона трендовой линии в пределах региона достаточно существенная, что указывает на значительные внутрирегиональные различия в скорости формирования древесных сообществ на залежах.

Необходимо отметить, что в силу отсутствия других данных, например, фактических сведений об изменении доли древесных сообществ на изученных залежах, подобная пространственно-временная оценка изменений растительного покрова была возможна исключительно на основе рядов NDVI, измеренных на основе данных MOD13Q1.

Таким образом, анализ многолетних рядов вегетационного индекса позволил выявить и изучить ряд особенностей формирования лесных сообществ на залежах лесостепи Центрального Черноземья, наблюдающихся в начале XXI века.

Выводы

Изучены возможности значений NDVI, вычисленных на основе информационных продуктов MOD13Q1, для количественной оценки величины проективного покрытия залежей кронами древесных сообществ, формирующихся на них, и реконструкции его многолетних изменений.

Из сезонных величин NDVI, измеренных в период вегетации, значения середины августа наиболее информативны для оценки покрытия залежей лесными насаждениями, формирующимися на них. Высокая корреляция отмечена для залежей с лиственными и хвойными породами.

Для залежей с формирующимися хвойными насаждениями теснота связи между величиной крон и значениями NDVI несколько выше, чем для залежей с лиственными породами.

Характер связи между проективным покрытием крон и значениями вегетационного индекса может быть описан логарифмической моделью, но приближается к обратной сигмоидной кривой.

Многолетние ряды NDVI для залежей лесостепи с разными типами лесных насаждений, формирующимися на них (лиственными, хвойными, смешанными) в период 2000–2018 гг., характеризуются положительной статистически значимой тенденцией к увеличению. Наибольший рост межгодовых значений вегетационного индекса характерен для оставленных агарных угодий с хвойными породами.

Пространственное распределение коэффициента наклона трендовой линии NDVI залежных земель в пределах Центрального Черноземья может быть использовано для картографирования современных географических внутризональных различий на территории лесостепи.