2000–2020年珠三角30m空间分辨率植被覆盖度时空变化数据集

其他数据论文 II 区论文（已发表） • 版本 ZH2 Vol 8 (4) 2023

A dataset of spatial-temporal changes of vegetation coverage with a spatial resolution of 30m in the Pearl River Delta from 2000 to 2020

沈明潭，谭炳香，侯瑞霞，黄逸飞

SHEN Mingtan, TAN Bingxiang, HOU Ruixia, HUANG Yifei

DOI: 10.11922/11-6035.csd.2023.0129.zh
CSTR: 32001.14.11-6035.csd.2023.0129.zh
: 10.57760/sciencedb.j00001.00877

： 2023 - 07 - 13

： 2023 - 10 - 09

： 2023 - 08 - 14

： 2023 - 12 - 26

map

5204 38 0

摘要&关键词

摘要：伴随着珠三角社会经济和人口密度的快速增长及城镇不断扩张，导致局部地区植被生境不断缩减，因此准确估计植被覆盖度、动态监测珠三角植被的分布格局和时空变化特征，对珠三角的生态安全至关重要。本研究在Google Earth Engine调用Landsat数据集，采用影像中值合成的方式获取研究区高质量无云影像，分别合成珠三角2000、2005、2010、2015和2020年5个时期的影像，保证植被覆盖度反演基础数据源的科学合理性；选择基于NDVI的像元二分法分别估测出珠三角5个时期的植被覆盖度，利用变异系数分析植被覆盖度变化的波动性，选用Slope趋势分析来探究植被覆盖度变化的时空变化特征，采用Hurst指数预测研究区植被覆盖度未来的发展趋势，系统地反映出珠三角地区植被覆盖度时空变化情况。本数据集通过高分辨率影像进行精度验证，精度为0.90，可为珠三角植被覆盖度分布格局和时空变化等方面的研究提供基础数据支撑。

关键词：植被覆盖度；时空变化；珠三角

Abstract & Keywords

Abstract: With the rapid growth of social economy and the expansion of cities and towns, coupled with increasing population density in the Pearl River Delta, vegetation habitats in some areas have been continuously reduced. Therefore, it is crucial to accurately estimate vegetation coverage and dynamically monitor the distribution patterns and temporal-spatial changes in vegetation in the Pearl River Delta for the sake of its ecological security. In this study, we utilized the Landsat dataset through Google Earth Engine to acquire high-quality cloud-free images in the study area. By means of median image synthesis, we synthesized the images of the Pearl River Delta in 2000, 2005, 2010, 2015 and 2020, respectively. To ensure the scientific rationality of the basic data source of FVC inversion, we selected the pixel dichotomy method based on NDVI to estimate the FVC in the five periods of the Pearl River Delta, and used the variation coefficient to analyze the fluctuation of FVC change in the five periods of the Pearl River Delta. Moreover, and we adopted the Slope trend analysis to explore the spatial-temporal characteristics of FVC changes. The Hurst index was employed to predict the future development trends of vegetation coverage in the study area, offering a comprehensive representation of the temporal and spatial changes of vegetation coverage in the Pearl River Delta region. The precision of this dataset was validated against high-resolution images with an accuracy of 0.90, which can provide fundamental data support for the research on the distribution patterns and spatio-temporal changes of vegetation coverage in the Pearl River Delta.

Keywords: fraction vegetation cover; space-time change; Pearl River Delta

数据库（集）基本信息简介

数据库（集）名称	2000–2020年珠三角30m空间分辨率植被覆盖度时空变化数据集
数据作者	沈明潭，谭炳香，侯瑞霞，黄逸飞
数据通信作者	谭炳香（tan@ifrit.ac.cn）
数据时间范围	2000–2020年
地理区域	地理范围为21°30'-24°40'N，111°21'-114°53'E。地理区域范围为肇庆市、江门市、佛山市、广州市、中山市、珠海市、东莞市、深圳市、惠州市。
空间分辨率	30 m
数据量	823 MB
数据格式	*.tif
数据服务系统网址	https://doi.org/10.57760/sciencedb.j00001.00877
基金项目	广东省林业局“珠三角森林城市群质量提升关键技术研究”（2021KJCX009）
数据库（集）组成	本数据集共8个tif文件， PRD_FVC_2000.tif、PRD_FVC_2005.tif、PRD_FVC_2010.tif、PRD_FVC_2015.tif、PRD_FVC_2020.tif分别为2000年、2005年、2010年、2015年和2020年5个时期的年度植被覆盖度影像； PRD_FVC_CV.tif为植被覆盖度变异系数影像；PRD_FVC_Slope_F.tif为植被覆盖度变化趋势影像； PRD_FVC_Hurst.tif为植被覆盖度未来变化趋势影像。

Dataset Profile

Title	A dataset of spatial-temporal changes of vegetation coverage with a spatial resolution of 30m in the Pearl River Delta from 2000 to 2020
Data corresponding author	TAN Bingxiang (tan@ifrit.ac.cn)
Data authors	SHEN Mingtan, TAN Bingxiang, HOU Ruixia, HUANG Yifei
Time range	2000 – 2020
Geographical scope	Zhaoqing City, Jiangmen City, Foshan City, Guangzhou City, Zhongshan City, Zhuhai City, Dongguan City, Shenzhen City, and Huizhou City (21°30'-24°40'N, 111°21'-114°53'E)
Spatial resolution	30 m
Data volume	823 MB
Data format	*.tif
Data service system	<https://doi.org/10.57760/sciencedb.j00001.00877>
Sources of funding	“Research on Key Technologies for Improving the Quality of Forest Urban Agglomeration in the Pearl River Delta”, supported by Guangdong Forestry Bureau, China (2021KJCX009).
Dataset composition	This dataset consists of 8 tif files: tif, PRD_FVC_2005.tif, PRD_FVC_2010.tif, PRD_FVC_2015.tif, PRd_fvc_2015. tif, PRd_fvc_2015. tif and PRD_FVC_2020.tif corresponding to five years, namely 2000, 2005, 2010, 2015 and 2020, respectively. PRD_FVC_CV.tif contains the images of vegetation coverage variation coefficient. PRD_FVC_Slope_F.tif contains the trend images of vegetation coverage changes. PRD_FVC_Hurst.tif contains the imaged of the future change trends of vegetation coverage.

map

引言

珠三角在国家整体发展中具有重要的战略地位，自改革开放以来经历了快速的经济增长、人口增长和城市化扩张，导致当地植被急剧退化和破坏，从而严重影响了生态系统的服务和功能，因此动态监测珠三角植被的时空变化和分布格局对当地的生态建设具有重要的参考意义。在已公开发布包含珠三角地区的植被覆盖度数据集产品中，高吉喜等基于空间分辨率为250 m的MODI3Q1 NDVI（Normalized differential vegetation index）数据，采用NDVI最大合成法的方式进行中国区域的植被覆盖反演^[1]。与该数据集相比，本研究选用空间分辨率更高的30 m Landsat5 TM和Landsat8 OLI数据，利用中值合成法获取研究区的影像数据，更能精确估测出珠三角地区的植被覆盖度，更科学有效地反映出珠三角地区的植被覆盖状况。

植被覆盖度（Fraction Vegetation Cover, FVC）是指地面植被（包括树叶、茎和树枝）的垂直投影面积占统计区域总面积的百分比，可以反映地表植被的生长状况，是监测生态环境的重要基础数据，在评估地表覆盖变化、土地利用方式、景观和生态系统功能等方面有重要指示作用^[2]，被广泛用于监测植被变化、区域环境质量评估、气候预测、土壤侵蚀评估等。遥感技术因其获取速度快、覆盖范围广、能够提取大量信息而广泛应用于研究植被动态监测、生态环境评估和土地覆盖分类等领域^[3]，是监测区域和全球FVC最高效的方法^[4]，可以动态监测植被覆盖变化及其影响因素。通过长时间序列和大规模的FVC研究来评估特定区域的生态环境条件和质量，在揭示区域环境条件的演变过程和预测未来变化趋势等方面具有重要意义。

1 数据来源和处理方法

本研究采用以下主要步骤对珠三角植被覆盖度进行分析：（1）在GEE（Google Earth Engine）平台调用Landsat数据，基于NDVI像元二分法估算珠三角2000、2005、2010、2015和2020年植被覆盖度；（2）利用变异系数分析珠三角5个时期植被覆盖度变化的波动性；（3）利用Slope趋势分析来探究植被覆盖度变化的时空变化特征；（4）采用Hurst指数预测研究区植被覆盖度未来的发展趋势。

1.1 数据来源

Landsat5 TM和Landsat8 OLI数据具有较高的空间分辨率（30 m），但由于研究区地处亚热带区域，受到长时间大面积高云量的影响，难以获取研究区完整清晰的影像。因此本研究基于GEE平台，利用一年中云量低于20%的数据，按照成像时间进行排列，选取位于中间时段的像元，合成出质量较高的研究区影像（见表1）。调用“LANDSAT/LT05/C01/T1_SR”数据集，分别利用2000-01-01至2000-12-31期间的141景Landsat5 TM影像合成2000年珠三角年度中值影像，2005-01-01至2005-12-31期间的119景Landsat5 TM影像合成2005年珠三角年度中值影像，2010年的数据量较少，需要扩大时间段才能得出质量较好的影像数据，因此利用2009-11-01至2010-12-31期间的101景Landsat5 TM影像合成2010年珠三角年度中值影像。调用“LANDSAT/LC08/C01/T1_SR”数据集，分别利用2015-01-01至2015-12-31期间的162景Landsat8 OLI影像合成2015年珠三角年度中值影像，2020-01-01至2020-12-31期间的152景Landsat8 OLI影像合成2020年珠三角年度中值影像。由于按照时间段内云量低于20%的数据进行筛选，每个时间段内符合的景数不同，且能完整合成出研究区清晰无云的影像数据需要的景数也不同，因而每个期间所调用的影像景数会有差别。所用的Landsat数据集都已经过辐射定标和大气校正，分别裁剪出珠三角5个时期的合成影像（见图1）。

表1 影像合成信息

合成影像 Resultant image	卫星数据集 Satellite data set	时间范围 Time range	景数 Number of scenes
2000年影像	LANDSAT/LT05/C01/T1_SR	2000-01-01—2000-12-31-	141
2005年影像	LANDSAT/LT05/C01/T1_SR	2005-01-01-—2005-12-31-	119
2010年影像	LANDSAT/LT05/C01/T1_SR	2009-11-01-—2010-12-31-	101
2015年影像	LANDSAT/LC08/C01/T1_SR	2015-01-01-—2015-12-31-	162
2020年影像	LANDSAT/LC08/C01/T1_SR	2020-01-01-—2020-12-31-	152

图1 珠三角Landsat8 标准假彩色影像

1.2 数据处理方法

1.2.1 植被覆盖度估测

本研究基于GEE平台生成珠三角年度中值植被覆盖度影像^[5-6]。主要的计算过程为：（1）输入珠三角范围；（2）构建Landsat遥感影像的时间序列数据集。分别调用LANDSAT/LT05/C01/T1_SR、LANDSAT/LC08/C01/T1_SR数据集，根据时间段进行筛选数据，然后进行去云处理；（3）生成珠三角年度中值遥感影像数据集。对构建后的影像数据集按不同的年份进行分组，逐像元计算生成每个时期的年度中值影像；（4）计算NDVI。调用normalized Difference函数计算NDVI（公式1）；（5）去除水体。为了避免水体对植被覆盖度估测的影响，调用mask函数将水体进行掩膜处理；（6）计算珠三角FVC。选用基于NDVI的像元二分法（公式2），分别计算每个时期的植被覆盖度；（7）导出计算结果。

1.2.2 归一化差值植被指数

Kriegler等人^[7]提出了一种简单的波段变换：用近红外光反射率减去红光反射率，再除以近红外光反射率加红光反射率，称为归一化差异植被指数（NDVI）。NDVI对植被的生长状态和空间分布特征非常敏感^[8]，可以在很大程度上消除地形和大气的干扰，并与植被覆盖度有显著的线性相关^[9]，在大规模植被动态监测等方面具有明显优势，NDVI值越高，植被覆盖度越高，是目前应用最广泛的植被指数^[10]。计算公式为：

\(NDVI=\frac{NIR-R}{NIR+R} \) （1）

式中，NDVI是归一化差异植被指数，NIR和R分别为传感器在近红外光和红光区域记录的反射率。

1.2.3 像元二分模型

植被覆盖度与NDVI之间具有显著的线性相关^[11]，利用植被指数与像元分解模型相结合的方法，具有一定的物理意义、反演算法简单^[12]，而且在很大程度上弱化了大气、土壤背景等因素的影响，被广泛应用于区域和全球范围内植被覆盖度提取^[13-14]。因此本研究利用基于NDVI的像元二分模型^[15]对研究区植被覆盖度进行估算：

\(\begin{array}{c}FVC=\frac{NDVI–{NDVI}_{soil}}{{NDVI}_{veg}–{NDVI}_{soil}}\end{array}\) （2）

式中FVC为植被覆盖度，本研究选取置信区间为[5％, 95％]的NDVI值作为NDVI_soil 和NDVI_veg^[16] 。NDVI_soil 表示无植被覆盖像元NDVI值，NDVI_veg 表示纯植被覆盖像元NDVI值。

1.2.4 变异系数

对研究区5个时期的植被覆盖度数据逐像元计算变异系数，探究研究区植被覆盖度的稳定性^[17]。CV值越大，植被变化的波动性越大，将波动性划分为4种类型，见表2。计算公式如下所示：

\(CV=\frac{\sqrt{\frac{1}{n}\sum _{i=1}^{n}{\left({FVC}_{i}-\stackrel{-}{FVC}\right)}^{2}}}{\stackrel{-}{FVC}}\) （3）

式中CV植被覆盖度的变异系数，n是研究期数，FVC_i 第i年的植被覆盖度值，\(\stackrel{-}{FVC}\)是研究期内植被覆盖度均值。

表2 植被覆盖度波动变化类型划分表

CV 类型 Type	划分依据 Rules
低波动	CV ≤0.4
中波动	0.4<CV≤0.8
较高波动	0.8<CV≤1.6
高波动	CV >1.6

1.2.5 Slope趋势分析

采用一元线性回归分析法中的斜率公式计算要素（y）的变化趋势（Slope）并进行显著性检验，根据Slope值和P值对植被覆盖度的变化趋势进行划分，将变化趋势分为7种类型^[18]，见表3。其中Slope的计算公式如下：

\(Slope=\frac{n×\sum _{i=1}^{n}\left(i×{y}_{i}\right)-\sum _{i=1}^{n}i\sum _{n=1}^{n}{y}_{i}}{n×\sum _{i=1}^{n}{i}^{2}-{\left(\sum _{i=1}^{n}i\right)}^{2}}\) （4）

式中，Slope表示要素的年际变化趋势；i为时期变量；n是研究期数，本研究为5个时期；y_i为第i年的要素。Slope>0和Slope<0分别表示要素序列随着时间增加和减少，Slope的绝对值越大表明要素变化越大。

采用F检验对变化趋势进行显著性检验，当检验结果P<0.05为显著变化，0.05<P<0.1为弱显著变化，P>0.1为不显著变化。统计量计算公式如下：

\(F=U×\frac{n-2}{Q}\) （5）

\(U= {\textstyle \sum_{i=1}^{n}} \left ( \hat{y_{i} } -\bar{y} \right ) ^{2} \) （6）

\(Q= {\textstyle \sum_{i=1}^{n}} \left ( y_{i} - \hat{y_{i} } \right ) ^{2} \) （7）

式中，U为回归平方和，Q为残差平方和，y_i为第i年的要素，\(\hat{y_{i} } \)为要素的回归值，\(\stackrel{-}{y}\)为要素的平均值，n表示期数。

表3 植被覆盖度变化趋势类型划分表

类型 Type	划分依据 Rules
类型 Type	FVC变化趋势 Slope	显著性 P
显著退化	<0	P<0.05
弱显著退化	<0	0.05≤P<0.1
退化不显著	<0	P≥0.1
无变化	0
显著改善	>0	P<0.05
弱显著改善	>0	0.05≤P<0.1
改善不显著	>0	P≥0.1

1.2.6 Hurst指数

基于R/S分析法的Hurst指数^[19]，能够用于定量描述时间序列信息的长期相关性，目前已广泛运用于分析区域植被覆盖未来变化趋势，Hurst指数取值范围存在以下3种情况：（1）若0.5<H<1，表明植被覆盖度时间序列具有长期相关特征，未来和研究期间内的变化趋势一致；（2）若H=0.5，表明植被覆盖度时间序列为相互独立的随机序列，未来和研究期间内的变化趋势无关；（3）若0<H<0.5，表明植被覆盖度时间序列具有反持续性，未来和研究期间内的变化趋势相反^[20]。

时间序列{FVC(t)}, t=1, 2, 3,…, n，对于任意一个正整数τ=1，定义均值系列：

\({FVC}_{\tau }=\frac{1}{\tau }\sum _{i=1}^{\tau }FVC\left(t\right) \tau =1,2,3,\dots ,n\) (8)

累计离差：

\(X\left(t,\tau \right)=\sum _{\mu =1}^{i}\left(FVC\left(\mu \right)-{FVC}_{\pi }\right) 1\le t\le \tau \) (9)

极差：

\(R\left(\tau \right)=\underset{1\le t\le \tau }{\mathit{max}}X\left(t,\tau \right)-\underset{1\le t\le \tau }{\mathit{min}}\left(t,\tau \right) \tau =1,2,3,\dots ,n\) (10)

标准差：

\(S\left(\tau \right)={\left[\frac{1}{\tau }\sum _{t=1}^{\tau }{\left(FVC\left(t\right)-{FVC}_{t}\right)}^{2}\right]}^{\frac{1}{2}} \tau =1,2,3,\dots ,n\) (11)

考虑比值\(R\left(\mathrm{\tau }\right)\)/\(S\left(\tau \right)\)≅R/S，若存在R/S∝τ^H，则说明时间序列{φ(t)}, t=1, 2, 3,…, n存在Hurst现象，其中H值称为Hurst指数。参考刘明霞等^[21]将植被覆盖度未来变化趋势分为7类，见表4。

表4 植被覆盖度未来发展趋势类型划分表

FVC预测 Forecast	划分依据 Rules
FVC预测 Forecast	FVC观测值趋势	Hurst指数
由退化变改善	<0	<0.5
退化趋势不确定	<0	=0.5
持续性退化	<0	>0.5
不变	=0
由改善变退化	>0	<0.5
改善趋势不确定	>0	=0.5
持续性改善	>0	>0.5

2 数据样本描述

2.1 数据组成

本数据集包括：（1）PRD_FVC_2000.tif、PRD_FVC_2005.tif、PRD_FVC_2010.tif、PRD_FVC_2015.tif、PRD_FVC_2020.tif，分别为2000年、2005年、2010年、2015年和2020年5个时期的年度植被覆盖度数据；（2）PRD_FVC_CV.tif为植被覆盖度变异系数数据；（3）PRD_FVC_Slope_F.tif为植被覆盖度变化趋势数据；（4）PRD_FVC_Hurst.tif为植被覆盖度未来变化趋势数据。上述数据空间分辨率为30 m，地理坐标系为WGS1984，投影坐标系为WGS_1984_UTM_Zone_49N，数据保存为GeoTIFF格式，整个数据集解压后大小约为2.75 G。

2.2 样本描述

基于GEE平台，利用公式1和公式2进行珠三角植被覆盖度反演，得到2000年、2005年、2010年、2015年和2020年5个时期的年度植被覆盖度影像（见图2）。基于珠三角5个时期的植被覆盖度反演结果，利用公式3对珠三角植被覆盖度变化的波动性进行分析，根据表2将变异系数划分为4种等级，结果见图2(a)；利用公式4计算出珠三角5个时期的植被覆盖度的时空变化趋势，采用F检验（公式5、6和7）对变化趋势进行显著性检验，根据表3的将变化趋势划分为7个等级，结果见图2(b)；利用公式4计算出珠三角植被覆盖度的变化趋势，再利用Hurst指数（公式8、9、10和11）来反映植被覆盖度变化趋势的持续性，从而揭示珠三角植被覆盖度未来的发展趋势，根据表4将珠三角发展趋势划分为7种类型，结果见图2(c)。

图2 珠三角5期植被覆盖度（2000、2005、2010、2015、2020）、变异系数(a)、变化趋势(b)、Hurst指数(c）空间分布

3 数据质量控制和评估

本数据集参考植被覆盖度的相关研究成果，在GEE平台利用中值合成法得到研究区清晰无云的影像，提高反演数据源的质量和精度。植被覆盖度的实际采样点观测值与像元二分法的估测值有较好的一致性，利用像元二分法对植被覆盖度进行反演是科学可靠的方法^[22-23]。因此选择基于NDVI的像元二分法进行珠三角地区植被覆盖度反演，根据研究区的气候条件和植被覆盖状况，在5%的置信水平下，采用累计频率在5%~95%的NDVI值作为上下阈值，在一定程度上减少遥感影像的噪点影响^[24]。基于NDVI的二分像元模型对图像辐射校正的影响不敏感，具有较高的精度和易操作性，不需要实际的植被覆盖度数据建模，仅需简单计算，易于解释^[18]，目前广泛应用于区域、省、流域、国家或全球范围等不同尺度的植被覆盖变化研究。参考Song等学者的研究^[25]，利用相关系数R²和均方根误差RMSE进行产品的精度验证和一致性检验。

3.1 精度验证

本研究在珠三角范围内随机生成100个采样点，通过world imagery 2020-07-22高分辨率影像进行精度验证，将相同坐标位置的2020年植被覆盖度遥感反演值和高分辨率影像目测值进行回归分析，两者相关系数R²=0.90，均方根误差RMSE=0.03，表明本数据集的精度基本满足区域性植被覆盖度研究的需要。

3.2 一致性检验

关于产品的一致性检验，目前公开发布的包含珠三角地区的植被覆盖度产品都为粗分辨率的，比如GEOV2的2000 m植被覆盖度产品，GLASS的500 m植被覆盖度产品，高吉喜团队的中国区域250 m植被覆盖度数据集（2000–2022）。本数据集以高吉喜团队2020年植被覆盖度数据为对比产品，在珠三角范围内随机生成1000个点，提取本研究2020年植被覆盖度数据和对比产品相同位置的值进行一致性检验，相关系数R²=0.69，均方根误差RMSE=0.09，表明两个产品整体上有较好的一致性，但根据world imagery 2020-07-22影像可以看出，本研究的FVC数据相较于对比产品表现出了更多的细节（见图3）。一致性不高的原因有：（1）虽然都采用了像元二分法，但NDVI的取值方式不一样，对比数据集采用的是NDVI最大值合成法，而本研究利用中值法合成影像，然后再计算NDVI；（2）两个数据集的空间分辨率不一样；（3）数据源不一样。

图3 world image 影像（a）、本研究2020年FVC（b）、对比产品2020年7月FVC（c）

4 数据价值

本数据集提供了珠三角2000、2005、2010、2015和2020年5个时期植被覆盖度的反演结果、变异系数、变化趋势和Hurst指数，可用于定量研究珠三角植被覆盖度的时空变化特征和预测未来变化趋势，探究珠三角植被覆盖度变化的机理，更好地了解珠三角植被变化的动态过程，对植被资源的管理及生态环境的保护具有重要的指导意义。在珠三角快速的城市化进程中，本数据集具有重要的实际应用价值。在城市规划方面，可以评估城市化对植被生态系统的影响，有助于指导城市的规划和发展。在生态保护建设方面，植被覆盖度研究有助于科学地认识生态系统及生态过程，可以指导生态环境的保护和恢复工作。在气候变化研究方面，监测植被变化可以深入了解气候变化对生态系统的影响，为可持续发展提供关键信息。

5 数据使用方法和建议

本数据集以GeoTIFF格式进行存储，在常用的GIS软件（ArcGIS, QGIS等）及其他数据分析软件（R, MATLAB等）都可进行查看和再分析，可以实现2000、2005、2010、2015和2020年珠三角植被覆盖度的定量评价和时空变化特征分析，可以应用于珠三角植被覆盖度动态变化与影响因子之间的相关性探究，以及定量分析自然因子和人为因子对珠三角植被覆盖度变化的相对影响程度等方面的研究。植被覆盖度的取值范围为[0,1]，为了监测地表植被覆盖细微的时空变化，同时参考相关的植被覆盖度数据集^[1]，本研究FVC值保留6位有效数字，用户在使用过程中可根据实际需要自行决定保留几位有效数字。

[1]

高吉喜, 史园莉, 张宏伟, 等. 中国区域250米植被覆盖度数据集（2000-2022）[DS/OL]. 国家青藏高原科学数据中心, 2022. DOI: 10.11888/Terre.tpdc.300330. [GAO J X, SHI Y L, ZHANG H W, et al. China regional 250m fractional vegetation cover data set (2000-2022) [DS/OL]. National Tibetan Plateau Data Center, 2022. DOI: 10.11888/Terre.tpdc.300330.]