Метод фенологической нормализации временных рядов данных ДЗЗ

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

В качестве исходных данных используются восстановленные многолетние ежедневные бесснежные композитные изображения AQUA/TERRA MODIS пространственного разрешения 230 м для красного и ближнего инфракрасного спектральных каналов. На рисунке 1 приведены этапы формирования используемого набора данных. Предварительная обработка включает в себя формирование масок облачного покрова на основе данных в красном, ближнем инфракрасном, голубом и коротковолновом спектральных диапазонах. Информация о геометрии наблюдения позволяет фильтровать измерения, полученные при некорректных условиях съёмки, и формировать области расположения теней от облаков. Пространственный гистограммный анализ позволяет уточнять границы классов искажённых объектов в масках облачного покрова. Маски содержат набор классов, позволяющих фильтровать снежный покров, различные типы облачности, атмосферную дымку и некорректные наблюдения.

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

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

Рисунок 1 ‒ Логическая схема формирования восстановленных многолетних
ежедневных бесснежных композитных изображений AQUA/TERRA MODIS

Рисунок 2 — Восстановленные ежедневные временные ряды наблюдений за разные годы	Рисунок 3 — Фенологическая нормализация временных рядов

На рисунке 2 приведены восстановленные временные ряды наблюдений за разные годы для одной точки, принадлежащей лесному типу растительного покрова. На примере видно, что в разные годы даты начала наблюдения варьируют в пределах от 99 дня года в 2016 г., до 165 дня года в 2008, даты окончания наблюдения от 231 дня года в 2008 до 286 дня в 2020, а период наблюдения от 65 дней в 2008 г. до 165 в 2020 г. При этом значения спектральных характеристик и вегетационных индексов, полученных в одной фенологической фазе развития растения, также будут сильно отличаться.На рисунке 3 изображено фенологическое совмещение временных рядов, состоящее из 17 разновременных срезов, нормированных на длину вегетационного периода. При этом первый фенологический срез соответствует началу вегетации SOS (Start Of Season), последний — окончанию вегетационного сезона EOS (End Of Season), а середина временного ряда соответствует максимуму вегетации MAX. Логическая схема алгоритма создания разновременных фенологических срезов представлена на рисунке 4.

Рисунок 4 — Логическая схема создания разновременных  фенологических срезов временных рядов

Метод фенологической нормализации временных рядов направлен на создание разновременных композитных изображений, каждое из которых содержит значения, соответствующее определённой фенологической фазе. При этом влияние межгодовых вариаций сводится к минимуму. Для этого на первом этапе ежегодно определяются значения максимума вегетации maxNDVI, и весенние и осенние минимумы LminNDVI и RminNDVI. Для нерастительных объектов с нетипичной годовой динамикой значений КСЯ так же определяются весенние и осенние максимумы LmaxNDVI и RmaxNDVI.

Рисунок 5 — Восстановленные ежедневные временные ряды наблюдений за разные годы (смена класса растительного покрова)	Рисунок 6 — Фенологическая нормализация временных рядов (смена класса растительного покрова)

Многолетняя статистика значений экстремумов может отражать как естественную изменчивость наступления фенологических фаз и условий начала и окончания периода наблюдений в году (см. рисунок 3), так и ситуации смены класса растительного покрова вследствие таких явлений как пожары, вырубки, заболачивание и т.д. (рисунок 5).На рисунке 5 видно, что в промежуток между 2008 и 2012 г. вследствие пожара участок, принадлежащий изображённому лесному покрову, оказался повреждён, и в 2012 г. характерной сезонной динамики для растительного покрова не наблюдается. Максимум вегетации в два раза ниже, чем в 2004 и 2008 г., а значения NDVI ни в один день 2012 года не достигает значений LminNDVI и RminNDVI для 2004 и 2008 лет. После происходит постепенное восстановление повреждённой лесной территории и, возможно, изменение класса растительного покрова. В такой ситуации фенологическая нормализация временных рядов может привести к потере информации за 2012 г., и размытию различий между разными растительными сообществами. Поэтому второй этап (см. рисунок 4) метода фенологического совмещения временных рядов включает в себя дифференциацию случаев межгодовой изменчивости фенологических характеристик и ситуаций смены класса растительного покрова.

Сезонные экстремумы за период наблюдения с 2001 по 2022 г. позволяют оценить стабильность класса растительности: среднее значение максимума NDVI за 22 года ¯maxNDVI и величину стандартного отклонения D(maxNDVI). Для стабильных классов определяются медианные значение весеннего и осеннего многолетнего минимума medianLmin и medianRmin. Если величина стандартного отклонения D(maxNDVI) превышает 5 % от ¯maxNDVI, предполагается смена класса и наличие двух многолетних весенних и осенних экстремумов, соответствующих разным классам растительного покрова. Минимумы для объектов с высокой межгодовой вариацией downLmin, upLmin, downRmin, upRmin определяются методом “k-средних”.

Многолетние экстремумы используются как ключевые точки для масштабирования временных рядов и нормирования их на одну величину:

                                                                               (1)

где Day_i — порядковый номер дня года во временном ряду для фенологической фазы i; StepL — параметр масштабирования слева от максимума вегетации; StepR — параметр масштабирования справа от максимума вегетации; N — параметр дробления временного ряда; (2N+1) — число фенологических срезов; n — порядковый номер композита в нормализованном наборе данных.

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

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

а	б
в	г

Рисунок 7 — Пространственное распределение дат наступления фенологических фаз для 2021 г. на территории России: a — 1/9 фенологический срез, соответствующий SOS; б — 3/9 фенологический срез; в — 5/9 фенологический срез, соответствующий MAX; г — 9/9 фенологический срез, соответствующий EOS

В таблице 1 приведены средние даты и величины стандартного распределения, соответствующие некоторым фенологическим композитам. Для первого из девяти построенных композитов, условно обозначенного в нормализованном наборе данных как 1/9 ширина распределения даты наступления фазы SOS наибольшая (см. таблицу 1), так как на нее сильнее всего влияют географо-климатические условия и наличие облачного покрова, так как от облачности зависит наличие данных и начало периода наблюдений в конкретном вегетационном сезоне (см. рисунок 7а). На рисунке 7б приведён третий композит, соответствующий достижению половины амплитуды значений NDVI на этапе роста растения. К середине сезона различия сглаживаются. Рисунок 7в соответствует максимуму вегетации, и разброс дат для этой фазы ниже, чем для SOS и EOS. На даты наступления фазы EOS (рисунок 7г) в большой степени влияет установление постоянного снежного покрова.

Таблица 2.2.1.1 — Даты, соответствующие разновременным композитным изображениям

Номер композита	Фенологическая фаза	Средняя дата	Стандартное отклонение
1/9	SOS	127	32
3/9	1/2 maxNDVI	155	22
5/9	maxMAX	198	19
7/9	1/2 maxNDVI	244	20
9/9	EOS	276	22

 

Гистограмма дат для 3/9 композита (рисунок 8a) приведена на рисунке 8б. Видно, что основная часть наблюдений распределена равномерно между 135-м и 165-м днями. На рисунке 8в представлен композит для 135-го дня года, а на рисунке 8г для 165-го дня года в RGB синтезе RED-NIR-RED. Эти даты в равной степени участвуют в создании фенологического композита, однако можно видеть, насколько сильно они различаются между собой. На 135-й день на северо-восточной части России ещё лежит снег, когда как на 165-й день, растительность на многих территорий уже давно миновали фазу SOS, и достигла максимума вегетации (жёлтый цвет на рисунке 7в).

а	б
в	г

Рисунок 8 — Сравнение фенологического среза и ежедневных данных. RGB синтез: RED-NIR-RED: a — 3/9 фенологический срез; б — гистограмма распределения дат для 3/9 фенологического среза; в — однодневный композит для 135-го дня года; г — однодневный композит для 165-го дня года

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

Значения NDVI, полученные на основе фенологического среза 3/9 (рисунок 9а) демонстрируют более равномерное распределение значений, чем для однодневных композитов (рисунки 9б и в), и отображают не географо-климатические вариации, а вариации значений NDVI для разных типов растительного покрова. На рисунке 9г приведено сравнение гистограмм для 135-го дня года, и соответствующих точек фенологического среза. Распределение значений для фенологического среза заметно уже, чем для 135-го дня. При этом видно, что часть территории на 135-го день только вышла из-под снега. Аналогичное сравнение для 165-го дня приведено на рисунке 165. Гистограмма однодневного композита шире, чем для фенологического композита, и смещена в область высоких значений NDVI. Любой однодневный композит, или построенный за определённый фиксированный временной промежуток, на большой территории будет содержать территории, находящиеся на разных этапах фенологического развития. Это можно увидеть на рисунке 10.

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

			б
a			в
г		д
Рисунок 9 — Сравнение фенологического среза и ежедневных данных, NDVI: a — 3/9 фенологический срез; б — 135-й день года; в — 165-й день года; г — гистограмма значений NDVI для фенологического среза и 135-го дня года; д — гистограмма значений NDVI для фенологического среза и 165-го дня года
Рисунок 10 — Сравнение гистограмм распределения значений NDVI для набора фенологических срезов и ежедневных композитов, соответствующих датам с наибольшей частотой значений разновременных композитных изображений: a — фенологические срезы 1/9, 3/9, 5/9; б — фенологические срезы 7/9 и 9/9; в — однодневные композиты за 135-й, 165-й, 195-й дни года; г — однодневные композиты за 255-й и 285-й дни года

На рисунках в и г представлены гистограммы для однодневных изображений, соответствующих средним значениям дат разновременных композитов (см. таблицу 1). Видно, что «движение» распределений для однодневных срезов смазанное, распределения более широкие и содержат несколько пиков.

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