Алгоритм оценки долговременных вариаций MODIS NDVI Текст научной статьи по специальности «Науки о Земле и смежные экологические науки»

Journal of Siberian Federal University. Engineering & Technologies, 2018, 11(1), 61-68

УДК 004.021:528.852.8

Algorithm of Assessment

of the MODIS NDVI Long-Term Variations

Bator V. Sodnomov*, Alexander A. Ayurzhanaev, Bair Z. Tsydypov and Endon Zh. Garmaev

Baikal Institute of Nature Management SB RAS 6 Sakhyanovoy Str., Ulan-Ude, 670047, Russia

Received 21.02.2017, received in revised form 16.05.2017, accepted 30.07.2017

An algorithm was developed for assessment of long-term variations of vegetation characteristics. The algorithm use NDVI data from spectroradiometer MODIS of Terra and Aqua satellites. The algorithm includes pre-processing procedures - the restoration of missing values, smoothing using Savitsky-Golay filter. To analyze the areas with snow cover the minimum length of the growing season is determined for the full length of the time series. Obtained after removing the seasonal component NDVI values are used to construct a linear regression and determine the trend. As a result of applying the algorithm the map of the spatial distribution of NDVI linear trends was created from 2000 to 2016 for the West Transbaikalia. Examples of visual verification of vegetation cover changes using satellite images of ultra-high spatial resolution are presented.

Keywords: algorithm, MODIS, NDVI, trend, Transbaikalia, vegetation.

Citation: Sodnomov B.V., Ayurzhanaev A.A., Tsydypov B.Z., Garmaev E.Zh. Algorithm of assessment of the MODIS NDVI long-term variations, J. Sib. Fed. Univ. Eng. technol., 2018, 11(1), 61-68. DOI: 10.17516/1999-494X-0009.

Алгоритм оценки долговременных вариаций MODIS NDVI

Б.В. Содномов, А.А. Аюржанаев, Б.З. Цыдыпов, Е.Ж. Гармаев

Байкальский институт природопользования СО РАН Россия, 670047, Улан-Удэ, ул. Сахьяновой, 6

Разработан алгоритм оценки долговременных вариаций характеристик растительного покрова, использующий значения вегетационного индекса NDVI спектрорадиометра MODIS спутников Terra и Aqua. Алгоритм содержит процедуры предварительной обработки данных -восстановление отсутствующих значений, сглаживание фильтром Савицкого-Голея. Для анализа территорий со снежным покровом определяется минимальная длина вегетационного

© Siberian Federal University. All rights reserved Corresponding author E-mail address: [email protected]

периода для всей длины временного ряда. Полученные после удаления сезонной компоненты значения NDVI используются для построения линейной регрессии и определения тренда. В результате применения алгоритма создана карта пространственного распределения линейных трендов NDVI с 2000 по 2016 гг. для территории Западного Забайкалья. Представлены примеры визуальной верификации изменения растительного покрова с использованием спутниковых изображений сверхвысокого пространственного разрешения.

Ключевые слова: алгоритм, MODIS, NDVI, тренд, Забайкалье, растительность.

Важнейшей информацией, которая позволяет исследовать происходящие на поверхности Земли процессы, являются данные дистанционного зондирования Земли [1]. С их помощью решается обширный класс задач, связанных с оценкой изменения климата, мониторингом стихийных бедствий, изучением динамики геосистем и т.д. [2, 3]. При этом исследование вариаций характеристик растительного покрова приобретает особое значение на фоне происходящих климатических изменений [4]. Для отслеживания динамики растительного покрова широко используются вегетационные индексы, среди которых наиболее распространен нормализованный относительный вегетационный индекс NDVI (Normalized Difference Vegetation Index), основывающийся на отражении и поглощении растительностью электромагнитного излучения в красном и ближнем инфракрасном диапазонах соответственно.

В настоящее время широкое распространение получили съемочные системы TM, ETM+ и OLI (спутников серии Landsat), AVHHR (NOAA), MODIS (Terra и Aqua), данные которых находятся в свободном доступе. Снимки со спутников Landsat обладают более высоким пространственным разрешением по сравнению с AVHHR и MODIS, но имеют частоту съемки один раз в 16 сут и часто оказываются неприменимы для анализа из-за наличия облачности. Для задач непрерывного мониторинга растительного покрова оптимальным представляется использование данных MODIS высокого временного разрешения (частота съемки может достигать до четырех раз в сут) и более высоким пространственным разрешением по сравнению с AVHHR c широким охватом съемки. При этом необходима разработка алгоритмов для обработки большого объема данных спутниковой информации.

Одним из наиболее используемых средств анализа временных рядов вегетационных индексов является программа TIMESAT [5]. Данная программа позволяет аппроксимировать временные ряды с целью определения фенологических параметров растительности с использованием асимметричной гауссовой функции, двойной логистической функции, фильтра Савицкого-Голея. TIMESAT обладает развитым инструментарием для удаления пиков, сглаживания данных, использования информации о типах ландшафтов. Недостатки TIMESAT: необходимость конвертации данных в бинарный формат с использованием стороннего программного обеспечения, отсутствие возможности полиномиальной аппроксимации трендовой компоненты временного ряда и анализа вклада различных фаз вегетационного периода в общий тренд.

Цель данной работы является создание алгоритма оценки многолетних изменений состояния растительности с использованием временных рядов вегетационного индекса MODIS NDVI. В результате разработан алгоритм обработки рядов NDVI высокого временного разрешения, проведены расчеты трендов NDVI для обширной территории Западного Забайкалья, представлены примеры изменения растительного покрова с использованием снимков сверхвысокого пространственного разрешения.

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

В работе использованы данные спектрорадиометра MODIS, установленного на спутниках Terra и Aqua. Данные доступны для загрузки из портала NASA EOSDIS. Исходными данными послужили временные ряды NDVI тематического продукта MOD13Q1, который помимо нормализованного разностного вегетационного индекса NDVI содержит усовершенствованный вегетационный индекс EVI [3]. Продукт представлен растровыми 16-дневными композитами с пространственным разрешением 250 м в синусоидальной проекции. Для создания композита используется метод максимального значения MVC (Maximum Value Compositing). Композит формируется из максимальных значений NDVI за 16 дней, таким образом исключая влияние атмосферной составляющей. Кроме того, для сокращения оптической дистанции выбирают данные со спутника, находящегося ближе надиру, что также уменьшает атмосферное влияние.

На рис. 1. представлена блок-схема разработанного алгоритма.

Значения NDVI из тематического продукта MOD13Q1 записывают в куб данных размерностью 4800*4800*N, где первые две координаты являются пространственными, третья - временной, где N - количество композитов. Данные содержатся в иерархическом hdf-файле и пред-ставленыввидеМОУМООО 0.

/ Запись полученных / значений в tif-формате

Рис. 1. Блок-схема алгоритма обработки данных Fig. 1. Flowchart of the data processing algorithm

Проводится предварительная обработка данных, заключающаяся в восстановлении отсутствующих значений и сглаживании временного ряда NDVI. Отсутствующие значения во временном ряду пикселя, обозначенные специальным флагом, заменяются средним от предыдущего и следующего значений вегетационного индекса. Если отсутствуют подряд два значения (т.е. нет данных о NDVI в двух композитах), то этот пиксель отмечается флагом «Baddata» и в дальнейшей обработке не участвует. Предусмотрена опция замены отсутствующего значения средним по всему временному ряду.

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

Чтобы исключить влияние на тренд значений вегетационного индекса открытой почвы и снега, временной ряд ограничивается значениями NDVI для периода вегетации [6]. За начало (SOS) иконец (EOS)вегетационногопериодапринято значение, когда NVDI становится соответственно выше и нижопнрооового знечения. Вработах [4, O] в пачестое пгооговыхзначений вегетационного периода пиеноеы зоочгния NDVH, рвеные ОДО еЯД. Поскольку дляназных сенсорое псиогтаое зннченпы NDVS мтжссотличаться, то в олгоритме оно задается пользователем. Далее овы е дечастсяы лиеельцгяеь ве геяецивышаге псриода(Ъ08),пос ттянноп дыя пхждого пикселе.За ^оипко перивда вяьвооциивбоЯиpaeтcянaиOвкee пеорнее зненение,аыкеы нец - наиболее раннее, т.е. длительность периода выбирается минимальной (рис. 2).

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

На рис. 3 представлен пример реализации алгоритма для одного временного ряда NDVI.

Для каждого временного ряда проводится оценка статистической значимости трендов с использованием критерия Фишера. В качестве нулевой гипотезы принимаем, что тренд равен 0. Значимым считаем тренд с 95%-ным уровнем (p < 0,05). Пользователь имеет возможность задавать уровень значимости.

Рис. 2. Схема определения вегетационного периода Fig. 2. Scheme ofthedefinition ofgrowing season

8000

-Raw data

-Savitzky-Golay filter

-Moving average filter

-Linear trend

2000

7000

6000

4000

5000

3000

2000 2002 2004 2006 2008 2010 2012 2014 2016

Year

Рис. 3. Пример обработки временного ряда NDVI. Черным цветом показаны исходные данные, красным -отфильтрованный ряд, синим-нелинейныйтревд, коричневым - линейныйтренд

Fig. 3. Ехатр1нргесesaing of NDVHimeseriHS.Black calor showsfreiHitial deta, eed а filtsred series, blue -non-linear trend, brown - linear trend

Реализованавозможность оценкивкладафазвеге тацнонногопериода в долговременную динамику NDVI. Так, определяются тренды максимума вегетации POS, суммы значений в фазе роста от SOS до POS и в фазе окончания вегетации от POS до EOS.

Выходные данные ремлизар ии рлгорнамспредераврнюнаяврасмровйм азобмажении tif-формаас.

Разработанныйалгоритм апробировотдля территории ЗападногоЗабайкалья. Входными даннымипостужили KoanooHTbiNDVI зт пертодо ^(Р0ОпоОШТ г. В резултаттеоопиксельной обработки поаучено болел Лбмлн ркгрогсионвыхулгтманийи тценньыих коэффициенты наклона (тренд).

Построена карта пространственного распределения трендов NDVI в Западном Забайкалье (рис. 4). Зоны с отрицательным трендом NDVI, характеризующим устойчивое снижение фото-синтетически активной биомассы, приурочены к степным растительным сообществам, которые в основном расположены в межгорных котловинах. Территории с наибольшей площадью отрицательных трендов NDVI располагаются в долинах рр. Селенга, Уда, Хилок и Баргузин. Зоны с положительным трендом вегетационного индекса NDVI приурочены к лесным ландшафтам.

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

Результаты и обсуждение

Рис. 4. Карта лшгейныхтрендов NDVIb Забаажвлье сеОООпо 201Л гг. 3 еленым икрасшлм щетасн обозначены статистичесни аначтмыо (р «с 0,05) положителхнып и отлицателкных пкетды с оотхетспкендо

Fig. 4. MapTflinear IrooifebfNDV1 in ТеansbaiCalic Caam 2B00 to 20c6.GTeeTandred coloraindicatea^nificant (p < 0.05)feaitiveapd neTftive CfeTdi

ь)

Рис. 5. Изменение ландшафтов в долине р. Топка (Бичурский район, Республика Бурятия): а - 2005 г.; б - 2010 г.; в - 2015 г.; 1 - участок вырубки леса; 2 - участок лесовосстановления; 3 - разработка Окино-Ключевского угольного карьера

Fig. 5. Changesof landscapes in the valley of the river Topka: a - 2005; b- 2010; c - 2015; 1 - the site of deforestation; 2 - the site of reforestation; 3 - development of the Okino-Klyuchevskiy coal mine

б)

а)

в)

Рис. 6. Тренды NDVI с 2000 по 2015 гг.: а - участок вырубки леса; б - участок лесовосстановления; в -разработка Окино-Ключевского угольного карьера

Fig. 6. NDVI Trends from 2000 to 2015: a - a site of deforestation; b - a site of reforestation; c - development of the Okino-Klyuchevskiycoalmine

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

В результатз яаеоманизиртвннно0 ябрдсотко времгнной серии ксмпясиаов N1^1 обнаружены а8у учясткл о гястремальиа толожртелсными о^ндамоМОУ1 (облясти тзваеосстановления) общой поощодью 1ль,° ям2, 47р участлов с экстремалыюневятевными трзкдами охватывгава .^^км2. Рбсрюбоятптлюмеяадикуоценми лендшрфвныxизмс66но2меяоxoдимо уточнить нз^чслаолр^еь^с^ноь^^г^ь^?^ физико-геягрефыажскил oьoлсинoсясй для oбшиpвoн лорри-ториаЗа&таалья.

На рис. 6 представлены графики NDVI, характерные для различных случаев изменения ландшафтов. Вырубка леса (рис. 6а) проявляется в резком снижении NDVI с 0,7 (хвойные леса) до 0,4 (степная растительность). Лесовосстановление характеризуется постепенным ростом NDVI (рис. 66). Разработка угольного карьера проявляется в снижении NDVI до значений 0,15-0,2, характерных для открытых почв (рис. 6в).

Выводы

Разработан алгоритм оценки долговременной динамики растительности на основе анализа временных рядов NDVI с широким пространственным охватом. Применение алгоритма позволило установить разнонаправленный характер трендов NDVI для Западного Забайкалья. Положительные тренды NDVI приурочены к лесным массивам, в степях межгорных котловин наблюдаются негативные тренды. В результате работы алгоритма выявлены участки антропогенного воздействия на растительность Забайкалья (вырубки леса, лесовосстановление).

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

Исследование выполнено при финансовой поддержке РФФИ в рамках проекта № 17-05-01059 и в рамках государственного задания БИП СО РАН.

Список литературы

[1] Шовенгердт Р. А. Дистанционное зондирование. Методы и модели обработки изображений. М.: Техносфера, 2010. 582 с. [Shovengerdt R.A. Remote sensing. Methods and models of image processing. Moscow, Tekhnosfera, 2010, 582 p. (in Russian)]

[2] Шукилович А.Ю., Федотова Е.В., Маглинец Ю.А Применение сенсора MODIS для оперативного мониторинга земель сельскохозяйственного назначения. Журнал СФУ. Техника и технологии, 2016, 9(7), 1035-1044 [Shukilovich A.Yu., Fedotova E.V., Maglinets Yu.A. Using MODIS data for agricultural areas monitoring, J. Sib. Fed. Univ. Eng. technol., 2016, 9(7), 1035-1044 (in Russian)]

[3] Huete A., Justice C., Van Leeuwen W. Modis vegetation index (MOD13): Algorithm theoretical basis document, Version 3. USGS Land Process Distributed Active Archive Center, 1999. 129 p.

[4] Медведева М.А., Савин И.Ю., Барталев С.А., Лупян Е.А. Использование данных NOAA-AVHRR для выявления многолетней динамики растительности Северной Евразии. Исследование Земли из космоса, 2011, 4, 55-62. [Medvedeva M.A., Savin I.Yu., Bartalev S.A., Lupyan E.A. Using NOAA-AVHRR data to identify long-term dynamics of vegetation of the Northern Eurasia, Earth observation from space, 2011, 4, 55-62 (in Russian)]

[5] Jonsson P., Eklundh L. TIMESAT - A program for analyzing time-series of satellite sensor data, Computers & Geosciences, 2004, 30, 833-845.

[6] Fensholt R., Proud S.R. Evaluation of Earth observation based global long term vegetation trends - Comparing GIMMS and MODIS global NDVI time series Remote Sensing of Environment, 2011, 119, 131-147. http://dx.doi.org/10.1016Zj.rse.2011.12.015

[7] Pang G. et al. Using the NDVI to identify variations in, and responses of, vegetation to climate change on the Tibetan Plateau from 1982 to 2012 Quaternary International, 2016, http://dx.doi. org/10.1016/j.quaint.2016.08.038