Long-term Coastline Monitoring in the Thanh Hoa Province (Vietnam) Using Landsat 5 and Landsat 8 Data

全文:

ВВЕДЕНИЕ

Вьетнам – страна с береговой линией протяженностью более 3260 км от Куангниня (Quang Ninh) до Кьензянг (Kien Giang), с множеством устьев крупных рек, на которых развиваются мангровые экосистемы. В последние годы ситуация с береговой эрозией и аккрецией (расширение территории суши в результате переноса и осаждения речных отложений в дельтах рек) в прибрежных районах Вьетнама резко обострилась, вызывая серьезные последствия в состоянии природных экосистем (Cuong и Cu, 2014).

Дистанционное зондирование и космический мониторинг эффективно применяются в решении многих научных вопросов при изучении аэрозольного состава атмосферы, прогнозирования землетрясений (Бондур, 2015; Бондур, Зверев, 2005). Данные оптических спутниковых изображений, таких как Landsat, Aster и др., часто используются для выделения береговых линий, оценки изменений устьев рек. Многие алгоритмы отделения водной поверхности от суши были разработаны для изображений, получаемых с датчиков дистанционного зондирования. Метод оцифровки, основанный на визуальном дешифрировании с использованием изображений в красном, ближнем инфракрасном и коротковолновом инфракрасном диапазонах, хотя и имеет достаточную точность, но трудозатратен и требует много времени (Yang et al., 2015). В исследовании Winasro и др. (Winasro and Budhiman, 2001), Alesheikh и др. (Alesheikh et al., 2007) использовались зеленый, ближний и средний инфракрасный каналы изображений Landsat TM для расчета отношений каналов GREEN/NIR и GREEN/SWIR₁ и дешифрирования границ водной поверхности на основе метода пороговой обработки. Этот метод использовался и в исследованиях (Trung and Khanh, 2016; Tran и Doan, 2017) для определения границы “вода/суша” и оценки изменения береговых и прибрежных территорий во Вьетнаме.

Низкая спектральная отражательная способность воды по сравнению с почвой и растениями, особенно в ближнем инфракрасном и коротковолновом инфракрасном диапазоне (Frazier and Page, 2000), позволяет эффективно использовать водные индексы для повышения точности классификации водных объектов по оптическим спутниковым снимкам. McFeeters (1996) предложил использовать нормализованный разностный водный индекс (Normalized Difference Water Index – NDWI), расчитываемый на основе зеленого и ближнего инфракрасного каналов спутниковых снимков Landsat. Индекс NDWI изменяется в диапазоне от –1 до 1, причем водная поверхность имеет положительные значения и представлена светлыми пикселями, остальные объекты характеризуются отрицательными значениеми и представлены темными пикселями (McFeeters, 1996; Duong et al., 2021). Gao (1996), независимо от McFeeters, также предложил водный индекс NDWI, основанный на использовании спектральной отражательной способности в зеленом и ближнем инфракрасном каналах (Gao, 1996). Модифицированный нормализованный разностный индекс (Modified Normalized Difference Water Index ‒ MNDWI) был разработан Xu (2006) на основе индекса NDWI, который использует зеленый и коротковолновый инфракрасный (SWIR1) каналы для дешифрирования водной поверхности на спутниковых изображениях Landsat (Xu, 2006). Полученные результаты показывают, что индекс MNDWI позволяет классифицировать водные объекты с большей точностью, нежели чем индекс NDWI, особенно на территориях с затопленной растительностью.

Yan и др. использовали расширенный водный индекс (Enhanced Water Index – EWI), чтобы подчеркнуть разницу между пересыхающими руслами рек и окружающей территорией и повысить точность дешифрирования таких специфических водных объектов (Yan et al., 2007). Wang и др. разработали многоканальный водный индекс (Multi-Band Water Index – MBWI), который обеспечивал большую надежность дешифрирования воды в сложных фоновых условиях с высоким уровнем шума (Wang et al., 2018). Rad и др. (2021) предложили индекс ANDWI (Augmented Normalized Difference Water Index), с использованием спектральных каналов в видимом диапазоне, NIR и SWIR1-2, чтобы максимизировать контраст между водной поверхностью и остальной территорией (Rad et al., 2021). Применяются и другие водные индексы, например нормализованный разностный индекс влажности (Normalized Difference Moisture Index – NDMI) (Gao, 1996), New Water Index (Ding et al., 2018), Water Ratio Index (WRI) (Shen and Li, 2010), Index of Water Surfaces (IWS) (Hassani et al., 2015), Automated Water Extraction Index (AWEI) (Feyisa et al., 2014). Отметим, что индекс AWEI имеет две модификации: AWEI_sh ‒ оптимальный для участков с мангровым растительным покровом и AWEI_nsh ‒ для водной поверхности без мангровой растительности. Данный индекс AWEI показал хорошие результаты в исследовании прибрежной зоны, покрытой протяженными мангровыми зарослями, Намдинь-Ниньбинь в северном Вьетнаме (Trinh et al., 2020).

Во многих исследованиях (Acharya et al., 2018a; Acharya et al., 2018b; Mustafa et al., 2017; Sarp and Ozcelik, 2017; Liu et al., 2022; Serban et al., 2022; Laonamsai et al., 2023) использовались различные водные индексы для классификации водных объектов и проводились оценки точности полученной классификации. В дополнение к оптическим изображениям также применялись радиолокационные изображения для оценки площади затопленой территории (Acar et al., 2012; Paz-Delgado et al., 2022) и анализа динамики береговой линии (Trinh et al., 2015; Baghermanesh et al., 2022; Wu et al., 2022).

Отмеченные выше исследования подтверждают эффективность использования дистанционного зондирования для оценки изменчивости береговой линии. Однако прибрежные территории характеризуются различными природными условиями, следовательно, сравнение и выбор соответствующего водного индекса для оценки вариации береговой линии представляет собой важную задачу, имеющую научное и практическое значение. В данной работе многозональные изображения Landsat использовались для расчета водных индексов NDWI, MNDWI, ANDWI, AWEI_sh, AWEI_nsh и BandWet. Метод классификации по порогу для дешифрирования границы водной поверхности на основе индекса ANDWI и последующая векторизация изображений для выделения береговой линии позволили оценить динамику береговой линии за период 1988‒2022 гг. провинции Тханьхоа.

ИСХОДНЫЕ ДАННЫЕ И МЕТОДИКА ОБРАБОТКИ ДАННЫХ

Исследуемая территория и данные дистанционного зондирования

Исследуемая территория. Областью исследования являлась прибрежная зона провинции Тханьхоа (Thanh Hoa), расположенная в Северо-Центральном регионе Вьетнама (рис. 1). Провинция имеет береговую линию протяженностью 102 км и территориальное море площадью 17 000 км², включая глубоководный порт Нги Сон (Nghi Son). Вдоль побережья расположены 5 крупных лиманов.

 

Рис. 1. Местоположение исследуемой территории. Справа вверху – административная карта, внизу – Landsat изображение провинции Тханьхоа.

 

Береговая линия провинции Тханьхоа подвергается динамическим процессам эрозии и расширения за счет речных отложений. Современное просадочное движение земной коры в отсутствии речных отложений приводит к тому, что морская вода вызывает энергичную эрозию берегов. Такие процессы наблюдаются в прибрежной полосе районов Хоанг Хоа (Hoang Hoa), Хау Лок (Hau Loc), Куанг Сюонг (Quang Xuong) и Нги Сон (Nghi Son). Однако для района Нга Соне (Nga Son) опускание территории сопровождается большим количеством поступающих речных взвесей, скорость седиментации частиц высока и в результате континент все больше и больше расширяется в сторону моря (Nguyen и Pham, 2016).

Исходные данные. Для оценки изменения береговой линии провинции Тханьхоа использовались 9 снимков Landsat, в том числе снимки Landsat 5 TM и Landsat 8/9 OLI за период с 1988 по 2022 г. В табл. 1 представлена информация об использованном спутнике дистанционного зондирования и времени съемки исследуемой территории. Изображения Landsat получены в засушливые сезоны разных лет, длящиеся с января по март, в которые погодные условия, как правило, оказываются одинаковыми (рис. 2).

 

Таблица 1. Тип спутника LANDSAT и дата съемки

№	Спутник	Дата съемки
1	Landsat 5 TM	16 марта 1988
2	Landsat 5 TM	20 марта 1992
3	Landsat 5 TM	02 января 1996
4	Landsat 8 OLI_TIRS	01 апреля 2000
5	Landsat 8 OLI_TIRS	20 марта 2004
6	Landsat 8 OLI_TIRS	13 февраля 2009
7	Landsat 8 OLI_TIRS	15 февраля 2014
8	Landsat 8 OLI_TIRS	10 февраля 2018
9	Landsat 8 OLI_TIRS	08 апреля 2022

 

Рис. 2. Мультизональные изображения Landsat исследуемой территории, полученные за период 1988–2022 гг.

 

Методика обработки данных

Данные дистанционного зондирования собирались и обрабатывались непосредственно на платформе облачных вычислений Google Earth Engine (GEE). GEE ‒ облачная платформа для геопространственного анализа, которая позволяет пользователям напрямую обрабатывать данные с помощью инструмента Code Editor без загрузки на персональный компьютер. В базе данных GEE находятся изображения, собранные за последние 40 лет со спутников дистанционного зондирования Landsat, MODIS, Sentinel. Работа с данными дистанционного зондирования на платформе GEE весьма удобна и позволяет сэкономить время и усилия при обработке изображений.

Для классификации водной поверхности использовались водные индексы NDWI, MNDWI, ANDWI, AWEI_sh, AWEI_nsh и BandWet. Индекс NDWI рассчитывался в соответствии со следующей формулой (Gao, 1996):

$NDWI=\frac{{\rho }_{GREEN}-{\rho }_{SWIR1}}{{\rho }_{GREEN}+{\rho }_{SWIR1}}$, (1)

где ρ_GREEN, ρ_SWIR1 ‒ значения коэффициента отражения в зеленом и коротковолновом инфракрасном (SWIR₁) каналах мультиспектрального изображения Landsat. Для изображения Landsat TM эти спектральные каналы соответствуют каналу 2 (зеленый) и каналу 5 (SWIR₁), а для изображения Landsat 8/9 эти каналы соответствуют каналу 3 (зеленый) и каналу 6 (SWIR₁).

Модифицированный индекс NDWI (MNDWI) вычислялся по формуле (2), в которой использовалась спектральная отражательная способность в зеленом (GREEN) и коротковолновом инфракрасном (SWIR₁) каналах. Индекс MNDWI определяется по следующей формуле (Xu, 2006):

$MNDWI=\frac{{\rho }_{GREEN}-{\rho }_{SWIR1}}{{\rho }_{GREEN}+{\rho }_{SWIR1}}$. (2)

Индексы ANDWI (Rad, 2021), AWEI_sh, AWEI_nsh (Feyisa, 2014) и BandWet (Zou, 2018) построены на основе использования спектральной отражательной способности в видимых (синий, зеленый, красный) и инфракрасных каналах (NIR, SWIR₁, SWIR₂) для увеличения контраста между водными и другими объектами. Индексы ANDWI, AWEI_sh, AWEI_nsh и BandWet вычисляются для изображений Landsat по следующим формулам:

Augmented Normalized Difference Water Index (ANDWI):

$ANDWI=\frac{{\rho }_{BLUE}+{\rho }_{GREEN}+{\rho }_{RED}-{\rho }_{NIR}-{\rho }_{SWIR1}-{\rho }_{SWIR2}}{{\rho }_{BLUE}+{\rho }_{GREEN}+{\rho }_{RED}+{\rho }_{NIR}+{\rho }_{SWIR1}+{\rho }_{SWIR2}}$, (3)

Automated Water Extraction Index (AWEI):

$AWE{I}_{sh}={\rho }_{BLUE}+2.5{\rho }_{GREEN}-1.5({\rho }_{NIR}+{\rho }_{SWIR1})-0.25{\rho }_{SWIR2}$, (4)

$AWE{I}_{nsh}=4({\rho }_{BLUE}-{\rho }_{SWIR1})-0.25{\rho }_{NIR}+2.75{\rho }_{SWIR2}$, (5)

BandWet:

$BandWet=\left\{\begin{array}{l}1 если MNDWI>NDWI или MNDWI>EVI и EVI<0.1\\ 0 если наоборот\end{array}\right\}$, (6)

где индекс EVI (Enhanced Vegetation Index) рассчитывается по формуле (7) (Huete et al., 2002):

$EVI=\frac{2.5\times ({\rho }_{NIR}-{\rho }_{RED})}{{\rho }_{NIR}+6\times {\rho }_{RED}-7.5\times {\rho }_{BLUE}+1}$. (7)

На основе полученных водных индексов выполнялось обнаружение береговой кромки методом автоматического определения порога Otsu. Метод Otsu разделяет пикселы полутонового изображения на два класса – “водная поверхность” и “суша”, таким образом, чтобы внутриклассовая дисперсия была минимальной (Otsu, 1979).

Блок-схема методики оценки изменения береговой линии по Landsat снимкам представлена на рис. 3.

 

Рис. 3. Блок-схема методики оценки изменений береговой линии по Landsat снимкам.

 

РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЕ

На основе спутниковых снимков Landsat 1988-2022 гг. были рассчитаны водные индексы NDWI, MNDWI, ANDWI, AWEI_sh, AWEI_nsh и BandWet, результаты представлены на рис. 4–9. На изображениях с индексам NDWI, MNDWI, ANDWI, AWEI_sh и AWEI_nsh водная поверхность представлена яркими белыми пикселями, территории суши представлена темным цветом. Наоборот, для водного индекса BandWet пиксели водной поверхности представлены темным цветом, объекты суши ‒ белым цветом.

 

Рис. 4. Пространственное распределение индекса NDWI на исследуемой территории.

 

Рис. 5. Пространственное распределение индекса MNDWI на исследуемой территории.

 

Рис. 6. Пространственное распределение индекса ANDWI на исследуемой территории.

 

Рис. 7. Пространственное распределение индекса AWEInsh на исследуемой территории.

 

Рис. 8. Пространственное распределение индекса AWEIsh на исследуемой территории.

 

Рис. 9. Пространственное распределение индекса BandWet на исследуемой территории.

 

В табл. 2 показаны минимальные, максимальные и средневзвешенные значения водных индексов NDWI, MNDWI, ANDWI, AWEI_sh, AWEI_nsh и BandWet прибрежной зоны изучаемой территории. В табл. 3 показаны значения порога Otsu при дешифрировании водных и неводных объектов на основе разных индексов.

 

Таблица 2. Максимальные, минимальные и средневзвешенные значения водных индексов

Значение	Водный индекс
	NDWI	MNDWI	ANDWI	AWEIsh	AWEInsh	BandWet
Минимальное	–0.57	–0.42	–0.57	–1.14	–0.67	–0.42
Максимальное	0.70	0.81	0.70	0.67	1.39	6.46
Средневзвешенное	0.08	0.27	0.08	–0.07	0.18	0.28

 

Таблица 3. Пороговые значения Otsu для разных водных индексов

Дата съемки	Водный индекс
	NDWI	MNDWI	ANDWI	AWEIsh	AWEInsh	BandWet
16 марта 1988	–0.02	0.24	0.18	–0.02	0.20	0.18
20 марта 1992	–0.03	0.23	0.17	–0.01	0.22	0.18
02 января 1996	–0.01	0.26	0.20	0.01	0.22	0.16
01 апреля 2000	–0.01	0.23	0.18	–0.03	0.19	0.19
20 марта 2004	–0.05	0.21	0.17	0.00	0.24	0.19
13 февраля 2009	–0.04	0.21	0.17	0.00	0.24	0.20
15 февраля 2014	–0.07	0.16	0.13	–0.03	0.21	0.20
10 февраля 2018	–0.05	0.17	0.13	–0.05	0.18	0.23
08 апреля 2022	–0.08	0.15	0.11	–0.07	0.18	0.24

 

Оценка точности дешифрирования границы водной поверхности по индексам NDWI, MNDWI, ANDWI, AWEI_sh, AWEI_nsh и BandWet выполнялась на основе статистического индекса Каппа и общей точности классификации, результаты представлены в табл. 4. Для этого 217 пикселей водной поверхности и 190 пикселей суши были выбраны случайным образом. Общая точность классификации и индекс Каппа рассчитывались с помощью программного обеспечения платформы облачных вычислений GEE. Общая точность классификации рассчитывалась как отношение суммарного количества правильно отклассифицированных объектов (сумма диагональных элементов матрицы ошибок) к общему числу объектов выборки. Индекс ANDWI показал наибольшую точность при классификации водной поверхности, выраженную как в максимальном значении коэффициента Каппа (0.961), так и в общей точности классификации (0.921). Поэтому в дальнейшей обработке при картографировани изменения береговой линии использовался только водный индекс ANDWI.

 

Таблица 4. Точность классификации водной поверхности по разным водным индексам

Водный индекс	NDWI	MNDWI	ANDWI	AWEIsh	AWEInsh	BandWet
Коэффициент Каппа	0.943	0.951	0.961	0.946	0.944	0.946
Общая точность	0.886	0.901	0.921	0.891	0.887	0.892

 

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

 

Рис. 10. Обзорная карта динамики береговой линии провинции Тханьхоа.

 

Детальная динамика береговой линии северной части провинции Тханьхоа

Прибрежная часть территории, граничащая с провинциями Ниньбинь и Намдинь, претерпела существенные изменения за рассматриваемый период. Общая протяженность поля аккреции составила 5.24 км, что соответствует площади 32.34 км². Наибольшая скорость прироста суши отмечена за период с 1996 г. по 2000 г. и составляла 457.5 м/год. Средняя скорость аккреции составляла 150 м/год. Начиная с 2000 г. и по настоящее время береговая линия района Нгиа Хунг подвергается незначительной эрозии. Однако, в целом, на северной части побережья провинции Тханьхоа продолжается значительное накопление речных отложений, хотя скорость аккреции стала ниже, по сравнению с периодом 1988–2000 гг. (рис. 11). Значительная динамика береговой линии этой территории связана с тем, что здесь находится устье крупнейшей реки северного Вьетнама – реки Дэй.

 

Рис. 11. Детальная карта изменения береговой линии северной части провинции Тханьхоа.

 

Детальная динамика береговой линии центральной части провинции Тханьхоа

Общая протяженность аккреционного побережья в районе Хоангхоа (Hoang Hoa) составляла 0.56 км, что дает среднюю за период с 1988 по 2022 г. скорость аккреции 16.5 м/год. В периоды 1992–1996 гг. и 2014–2018 гг. максимальная скорость аккреции достигала 120 м/год, средняя скорость аккреции равна 30 м/год (рис. 12).

 

Рис. 12. Детальная карта изменения береговой линии центральной части провинции Тханьхоа.

 

Береговая линия в пределах черты города Сам Сон (Sam Son) в течение 1988–1996 гг. менялась медленно, в основном за счет приращения территории, а затем в период 1996–2022 гг. преобладала эрозия побережья со средней скоростью около 10 м/год.

Детальная динамика береговой линии южной части провинции Тханьхоа

В районе Куанг Сюонг (Quang Xuong) береговая линия увеличивалась в течение периода 1988–2004 гг. со скоростью 21.3 м/год и достигла общей протяженности 0.34 км. Далее с 2004 г. по 2009 г. прибрежная зона подвергалась медленной эрозии. Затем, с 2009 по 2022 г., снова береговая линия в этом районе начала захватывать морскую территорию со средней скоростью 7.7 м/год (рис. 13).

 

Рис. 13. Детальная карта изменения береговой линии южной части провинции Тханьхоа.

 

В сельскохозяйственной коммуне Куанг Хай процесс эрозии береговой линии сменялись процессом аккреции и снова эрозии береговой кромки. Так, в период с 1988 г. по 1992 г., береговая линия размывалась со средней скоростью 7.5 м/год. Начиная с 1992 г. по 1995 г. имела место аккреция и протяженность береговой линии увеличилась до 110 м (что в среднем составляет 27.5 м/год). Далее, в течение периода 1996–2000 гг. вновь происходила эрозия со средней скоростью 20 м/год. За период с 2000 г. по 2009 г. протяженность поля аккреции достигла 0.32 км, что соотвествует средней скорости более чем 35 м/год. Однако начиная с 2009 г. по настоящее время в этом районе наблюдается максимальная эрозия со средней скоростью более 30 м/год.

В целом на южном побережье провинции Тханьхоа за исследуемый период преобладала эрозия береговой линии с максимальной скоростью 38 м/год и средней скоростью около 10 м/год.

ВЫВОДЫ

Получена динамика изменения береговой линии провинции Тханьхоа (северо-центральный Вьетнам) за более чем тридцатилетний период (1988 – 2022 гг.) с помощью спутниковых снимков Landsat 5 TM и Landsat 8 OLI_TIRS. Шесть водных индексов, включая NDWI, MNDWI, ANDWI, AWEI_sh, AWEI_nsh и BandWet, сравнивались для оценки возможности дешифрирования береговой линии на основе метода порога Otsu. Точность классификации водной поверхности по индексам оценивалась с помощью коэффициента Каппа и общей точности классификации, что позволило выбрать индекс ANDWI, как имеющий наибольшую точность.

Выполнена векторизация береговых линий на бинаризованных по порогу Otsu изображениях ANDWI. Границы береговых линий накладывались друг на друга, в результате получены обзорная и детальные карты динамики береговой линии провинции Тханьхоа. Анализ полученных результатов показал, что береговая линия претерпела сильное изменение за период 1988-2022 гг., в ходе развития таких процессов, как аккреция и эрозия. Аккреция береговой линии активно происходила в северной прибрежной зоне провинции, особенно в тех участках, где река Дэй впадает в Восточное море. Максимальная скорость аккреции на этом участке и достигала 457 м/год. С другой стороны, южный прибрежный район провинции подвергается эрозии с максимальной скоростью 38 м/год и средней скоростью около 10 м/год.

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

ИСТОЧНИК ФИНАНСИРОВАНИЯ

Работа выполнена в рамках научного проекта: “Применение дистанционного зондирования и ГИС-технологий в изучении изменения береговой линии провинции Тханьхоа за период 2000–2023 гг.”, код: T2023-03-04TD. Авторы благодарят Вьетнамскую сельскохозяйственную академию (VNUA – Vietnam National University of Agriculture) за финансирование проекта.

作者简介

Thi Giang Le

Vietnam National University of Agricuture

Email: trinhlehung@lqdtu.edu.vn
越南, Hanoi

Le Hung Trinh

Le Quy Don Technical University

编辑信件的主要联系方式.
Email: trinhlehung@lqdtu.edu.vn
越南, Hanoi

V. Zablotskii

Moscow State University of Geodesy and Cartography

Email: trinhlehung@lqdtu.edu.vn
俄罗斯联邦, Moscow

Quoc Vinh Tran

Vietnam National University of Agricuture

Email: trinhlehung@lqdtu.edu.vn
越南, Hanoi

Xuan Bien Tran

Hanoi University of Natural Resources and Environment

Email: trinhlehung@lqdtu.edu.vn
越南, Hanoi

Thi Phuong To

Thanh Dong University

Email: trinhlehung@lqdtu.edu.vn
越南, Hai Duong province

Van Phu Le

Le Quy Don Technical University

Email: trinhlehung@lqdtu.edu.vn
越南, Hanoi

Van Tho Le

Thai Nguyen University of Agriculture and Forestry

Email: trinhlehung@lqdtu.edu.vn
越南, Thai Nguyen

参考

补充文件