Extraction of Sugarcane Plantation Distribution in Western Guangdong by Combining Remote Sensing and Statistic Data and Its Spatiotemporal Analysis

  • YANG Yingpin , 1, 2 ,
  • WU Zhifeng , 1, * ,
  • HUANG Qiting 3 ,
  • LUO Jiancheng 4, 5 ,
  • WU Tianjun 6 ,
  • DONG Wen 4 ,
  • HU Xiaodong 7 ,
  • XIAO Wenju 1
Expand
  • 1. School of Geography and Remote Sensing, Guangzhou University, Guangzhou 510006, China
  • 2. Key Laboratory of Natural Resources Monitoring in Tropical and Subtropical Area of South China, Ministry of Natural Resources, Guangzhou 510670, China
  • 3. Agricultural Science and Technology Information Research Institute, Guangxi Academy of Agricultural Sciences, Nanning 530007, China
  • 4. Aerospace Information Research Institute, Chinese Academy of Sciences, Beijing 100094, China
  • 5. University of Chinese Academy of Sciences, Beijing 100049, China
  • 6. School of Science, Chang'an University, Xi'an 710064, China
  • 7. Image Sky International Co., Ltd., Suzhou 215163, China
*WU Zhifeng, E-mail:

Received date: 2022-08-12

  Revised date: 2022-10-16

  Online published: 2023-04-27

Supported by

China Postdoctoral Science Foundation(2021M690769)

National Natural Science Foundation of China(42201413)

National Natural Science Foundation of China(42071316)

National Natural Science Foundation of China - Guangdong Joint Fund(U1901219)

Key Laboratory of Natural Resources Monitoring in Tropical and Subtropical Area of South China, Ministry of Natural Resources(2022NRM0004)

Abstract

Western Guangdong is one of the three major sugarcane producing areas in China. Sugarcane in western Guangdong is mainly distributed in Zhanjiang, with more than 2 million sugarcane farmers. In 2020, the sugarcane planting area in Zhanjiang reached 130 030 hectares. Mapping sugarcane plantation and analyzing its spatiotemporal characteristics in western Guangdong are of great values for making polices in sugarcane industry, optimizing the distribution of sugarcane plantation, and promoting production efficiency. Remote sensing technology provides an efficient way to acquire land cover information. In this study, the sugarcane plantation distribution information in 2000, 2008, and 2020 was acquired based on Landsat remote sensing data and statistics data in sugarcane planting areas. Following steps were implemented: preprocessing of Normalized Difference Vegetation Index (NDVI) time series, construction of reference NDVI of sugarcane, extraction of the amplitude and the maximum of NDVI time series, and identification of sugarcane using the Time Weighted Dynamic Time Warping (TWDTW) method. The TWDTW method calculated the distance between NDVI time series of unknown pixels and sugarcane pixels, and a distance threshold was set via the statistics data to acquire the sugarcane plantation distribution. Based on the extracted distribution of sugarcane plantation, the kernel density of sugarcane distribution was calculated to analyze the spatial clustering characteristics of sugarcane planting areas. Landscape pattern indexes such as the percentage of landscape, average patch area, patch density, and aggregation index were calculated to analyze the spatial distribution characteristics of sugarcane planting patches. The topographic characteristics of sugarcane planting areas were also analyzed based on DEM data. The results showed that: ① the TWDTW model could realize sugarcane identification with high accuracy by combing remote sensing time series data and statistics data. In 2000, 2008 and 2020, the average accuracy of sugarcane mapping reached 87.62%; ② Sugarcane was mainly distributed in Suixi, Leizhou, and Xuwen in western Guangdong. The distribution of sugarcane planting in Suixi and Leizhou presented a pattern of high-density aggregation in multi centers; ③ From 2000 to 2020, the average area of sugarcane patches increased, the patch density decreased, and the aggregation index increased in Suixi and Leizhou, which indicated that the layout of sugarcane plantation had been significantly adjusted in these areas, and sugarcane production showed a trend of intensive production; ④ In Suixi and Leizhou, most of sugarcane was planted in flat areas, showing a great potential to develop mechanized production.

Cite this article

YANG Yingpin , WU Zhifeng , HUANG Qiting , LUO Jiancheng , WU Tianjun , DONG Wen , HU Xiaodong , XIAO Wenju . Extraction of Sugarcane Plantation Distribution in Western Guangdong by Combining Remote Sensing and Statistic Data and Its Spatiotemporal Analysis[J]. Journal of Geo-information Science, 2023 , 25(5) : 1012 -1026 . DOI: 10.12082/dqxxkx.2023.220585

1 引言

甘蔗是世界上最重要的糖料作物,贡献了全球80%左右的食糖产量,也是栽培作物中可发酵量最高的能源作物[1],随着全球人口的日益增长,食糖和能源的需求与日俱增。我国甘蔗产量位居世界第三,年产量约1亿t[2],广东省西部(简称粤西)是我国甘蔗三大主产区和优势生产基地之一,粤西甘蔗主要分布在湛江,蔗农200万余人,2020年湛江市甘蔗种植面积达1.30×105 hm2,约占全省甘蔗种植总面积的81.82%[3-4]。掌握粤西甘蔗种植空间分布信息及其时空特征,对制定甘蔗种植管理政策、促进甘蔗集约化生产、合理规划生产布局具有重要价值,对促进我国甘蔗产业可持续发展具有深远意义。
传统田间调查方式需耗费大量时间和人力,难以追溯历史种植信息。遥感具有范围大、探测周期短、成本低的优势,在作物类型识别和空间分布监测中应用广泛[5-7]。遥感影像空间特征、光谱特征和物候特征是作物识别的三大理论基础[8-12],通过获取作物类型样本并计算多维特征,利用监督分类模型可实现作物类型识别。不同作物因生长季相节律的差异,在多时相遥感影像上呈现出相异的物候特征,随着作物生长发育过程,NDVI时序曲线上升、到达饱和、下降的过程体现了作物萌芽生长、生长旺盛、成熟衰老的全生命周期过程。基于物候特征的作物识别方法具有较强的机理可解释性,近年来备受关注[12-17]
在甘蔗种植分布遥感提取方法方面,大体可以分为基于机器学习的甘蔗识别方法和基于遥感时序匹配的甘蔗识别方法。基于机器学习的甘蔗识别研究,通过多源遥感数据提取光谱、物候、空间等多维特征,利用支持向量机、随机森林等监督分类模型,进行甘蔗识别。Vieira等[18]和Zhou等[19]通过图像分割获取田块对象,并提取田块的光谱、纹理和空间特征,利用决策树分类方法实现对象级甘蔗识别。Wang等[20]基于Landsat7/8和Sentinel-2数据构建耕地NDVI时间序列曲线,提取开始生长时间、生长速率、生长季长度、生长季结束时间等关键物候特征,通过决策树分类算法实现广西甘蔗分布提取。杨颖频等[17]利用Landsat8和Sentinel-2数据构建作物NDVI时序曲线,提取物候特征和多时相多光谱特征,利用随机森林模型进行分类特征优选并实现甘蔗识别。Jiang等[21]利用Sentinel-1的VV、VH特征对广东湛江的甘蔗进行早期识别,并评估时序长度对早期识别精度的影响。谢鑫昌等[22]基于Landsat8近红外、短波红外、蓝波段反射率和NDVI指数,并将数字高程模型(DEM)数据作为辅助特征,通过随机森林分类方法,实现广西2014—2018年的甘蔗识别。基于机器学习的甘蔗识别方法,大多针对分类特征的提取、筛选和优化进行深入研究,该类方法属于统计学方法,分类模型的训练过度依赖样本和影像数据,训练的分类模型可迁移性较弱。基于遥感时序匹配的甘蔗识别方法,通过甘蔗样本数据构建甘蔗NDVI参考时序曲线,利用时间序列相似度计算模型,衡量未知地物NDVI时序与甘蔗NDVI时序的相似性。Zheng等[23]基于时间加权的动态时间规整(TWDTW)模型,逐像元计算NDVI时序与甘蔗参考NDVI时序的距离,协同遥感与统计数据提取了2016—2019年的巴西甘蔗分布。TWDTW模型是在动态时间规整(DTW)方法基础上发展起来的针对植被时序曲线匹配的模型,通过对2个长度不同的时间序列进行延伸和压缩,比较两条曲线波形的相似性,该方法具有抗物候差异、仅需少量样本的特点,在大范围作物识别中具有明显优势[24]。基于遥感时序匹配的甘蔗识别方法,充分利用了甘蔗的生长物候特征,具有较强的作物生理学含义。
近年来,遥感与统计数据协同的方法在作物种植分布提取中得到发展和应用,与仅利用遥感数据进行作物分布提取不同,该方法先利用遥感数据找到像元点之间的相似性,再以统计数据为基准,对总体数量进行空间化,基于这套思路,研究人员在水稻、玉米、小麦、甘蔗等作物分布提取方面开展了应用,并取得了较高的提取精度[23,25-27]
尽管已有研究利用遥感对粤西地区甘蔗分布进行提取[19,21],但多以区县为研究对象,且尚无研究针对粤西甘蔗的历史分布情况进行提取,也缺乏对粤西甘蔗空间分布特征的分析。本研究基于甘蔗生长物候特征,利用遥感时序匹配模型方法,协同遥感数据与统计数据,选取2000、2008、2020年3个年开展粤西地区甘蔗种植分布提取,基于核密度估算方法分析甘蔗种植区空间聚集特征,利用景观格局指数计算方法[28],定量分析甘蔗种植斑块的景观格局特征,结合DEM数据分析蔗区地形特征。本研究旨在为该地区甘蔗种植过程生产资源调度、甘蔗种植布局调整、糖料供给能力评估提供数据基础和科学指导。

2 研究区概况和数据源

2.1 研究区概况

本研究选取粤西甘蔗主产区湛江市为研究区(图1)。湛江地处祖国大陆最南端、广东省西南部(109°40′E—110°58′E,20°13′N—21°57′N),包括整个雷州半岛及半岛北部的一部分,辖区总面积132 63 km2。湛江东濒南海,南隔琼州海峡与海南省相望,西临北部湾,西北与广西壮族自治区的合浦、博白等县毗邻。湛江的地势大致中轴高、东西侧低,南北高而中间低,起伏和缓,多为平原和台地。湛江地处北回归线以南的低纬地区,属于热带北缘季风气候,终年受海洋气候的调节,冬无严寒,夏无酷暑,年平均气温在22.7 ~23.5 ℃。湛江最主要的农作物为甘蔗和水稻,甘蔗通常3月开始萌芽生长,6—9月为伸长期,11月到达成熟期,砍收期从11月持续至次年2月。水稻可种植两季,3月底4月初播种早稻,7—8月早稻收割、晚稻播种,10月底晚稻收割。
图1 研究区地理位置概况

Fig. 1 Geographical location of the study area

2.2 数据源

本研究采用的数据源包括Landsat系列遥感数据、作物类型样本数据、统计年鉴数据、DEM数据。

2.2.1 遥感时间序列数据

甘蔗种植分布提取使用的遥感数据源来自Landsat系列卫星,包括Landsat-5和Landsat-8。Landsat影像空间分辨率为30 m,时间分辨率为16 d。在Google Earth Engine云平台检索整体云量小于80%的Landsat地表反射率数据,计算NDVI特征,利用QA_Pixel质量控制波段筛选云污染像元,对NDVI数据进行掩膜处理,批量下载NDVI时间序列数据。本研究以10年为间隔,对2000—2020年的甘蔗分布进行提取并开展时空分析,由于2010年和2009年的时序影像数据质量较差,利用2008年替代2010年开展研究,最终选取2000、2008、2020年3个年的遥感影像进行甘蔗分布提取。考虑到甘蔗的生长季长度,监测时间窗口设定为当年1月至次年2月,影像列表如表1所示。
表1 Landsat影像数据

Tab. 1 Landsat images

序号 2000年Landsat5成像时间 2008年Landsat5成像时间 2020年Landsat8成像时间
1 2000-01-24 2008-02-15 2020-01-15
2 2000-03-28 2008-03-02 2020-01-31
3 2000-04-29 2008-04-19 2020-03-03
4 2000-05-15 2008-05-05 2020-04-20
5 2000-05-31 2008-06-22 2020-05-06
6 2000-07-02 2008-07-24 2020-06-07
7 2000-08-19 2008-08-25 2020-06-23
8 2000-09-04 2008-09-10 2020-07-09
9 2000-10-06 2008-09-26 2020-07-25
10 2000-11-07 2008-10-28 2020-08-10
11 2000-11-23 2008-11-13 2020-09-11
12 2000-12-09 2008-11-29 2020-09-27
13 2000-12-25 2008-12-15 2020-11-30
14 2001-02-11 2009-01-16 2021-01-01
15 2009-02-17 2021-01-17
16 2021-02-02
17 2021-02-18

2.2.2 作物类型样本数据

作物类型样本数据来源于野外地面调查以及Google Earth高分辨率影像目视解译。2020年7月研究团队在湛江开展实地调查,获取了2008个作物类型采样数据,包括甘蔗样本860个,非甘蔗样本1148个,其中香蕉样本124个、水稻样本156个、火龙果样本20个、蔬菜样本36个、撂荒地样本64个。基于地面调查样本数据,观察Google Earth高分辨率影像上甘蔗的空间特征,通过目视解译并结合NDVI时间序列曲线形态特征,补充样本数量并构建作物类型样本数据集。如图2所示,生长旺盛期时的甘蔗由于叶片高度覆盖和叶片间的细微孔隙,影像上呈现出较饱和的绿色调和粗颗粒状纹理。构建的样本集中共包含2020年1039个甘蔗样本、1553个非甘蔗样本,2008年1059个甘蔗样本、1490个非甘蔗样本; 2000年1008个甘蔗样本、1522个非甘蔗样本。
图2 甘蔗样地照片及Google Earth高分辨率影像

Fig. 2 Photo of sugarcane field and the high-spatial-resolution Google Earth image

2.2.3 统计年鉴数据

2000、2008、2020年甘蔗种植面积统计数据来源于《广东农村统计年鉴》[29],统计数据的空间尺度为区县尺度(图3)。本研究以甘蔗种植面积统计数据为约束条件,辅助甘蔗种植分布遥感提取。
图3 湛江各区县甘蔗种植面积统计数据

Fig. 3 Statistics of the sugarcane planting areas in counties of Zhanjiang

2.2.4 DEM数据

下载SRTM DEM数据(https://lpdaac.usgs.gov/products/srtmgl1v003/[30],空间分辨率为30 m,用于分析甘蔗种植区的地形特征。湛江最高海拔为384 m,徐闻县北部、雷州市南部、廉江市北部地势相对较高,其余多为低丘陵区、缓坡台地区和平原区(图4)。
图4 湛江高程分布

Fig. 4 Elevation map of Zhanjiang

3 研究方法

本研究的基本框架如图5,甘蔗种植分布信息的提取流程包括NDVI时序数据预处理、甘蔗参考NDVI时序曲线构建、NDVI时序特征提取、基于TWDTW模型的甘蔗识别。在甘蔗种植分布提取结果的基础上,进一步分析甘蔗种植区空间聚集特征、甘蔗种植斑块景观格局特征和甘蔗种植区地形特征。
图5 研究基本框架

Fig. 5 Framework of this study

3.1 协同遥感与统计数据的甘蔗种植分布提取

3.1.1 NDVI时序数据预处理

利用多时相NDVI构建NDVI时序曲线,为了进一步消除异常观测值的影响,对NDVI时序曲线中的异常值进行剔除。异常值剔除方案如下:对每个NDVI观测点(Pi)建立大小为3的滑动窗口,窗口内包含Pi的前(Pi-1)、后(Pi+1)观测点。对Pi-1Pi+1观测点进行线性内插,记作Pi',比较观测点Pi和插值点Pi'的NDVI取值大小,当PiPi'的NDVI之差超过0.2并且Pi距离上一次观测的天数超过14天,或者PiPi'的NDVI之差超过0.5并且Pi距离上一次观测的天数超过21天时,将Pi判定为异常观测(图6)。以表1中Landsat数据的观测时相为采样间隔,对NDVI时序曲线进行线性内插,填补缺失值。
图6 NDVI时序曲线异常值剔除示意

Fig. 6 Schematic diagram of eliminating outlier of NDVI temporal curve

3.1.2 甘蔗参考NDVI时序曲线构建

在研究区随机选取了300余个甘蔗样本,基于预处理后的NDVI时序曲线,获取各个观测时相上所有样本NDVI的四分位数,提取其中处于25%~75%位置之间的NDVI值,建立NDVI取值区间,计算各个时相上甘蔗NDVI平均值,构建甘蔗NDVI参考时序曲线,2020年甘蔗NDVI参考时序曲线如图7所示。图8绘制了研究区双季稻、菠萝、花生等其他典型作物类型的NDVI时序曲线,与甘蔗NDVI时序曲线形态均体现出较大差异。
图7 甘蔗NDVI参考时序曲线

Fig. 7 Reference NDVI time series of sugarcane

图8 研究区除甘蔗外其他典型作物类型的NDVI时序曲线

Fig. 8 NDVI time series of typical crop types other than sugarcane in the study area

3.1.3 NDVI时序特征提取

基于NDVI时序曲线提取耕地复种指数、NDVI振幅和NDVI全年最大值。复种指数的提取方式如下,识别NDVI时序曲线极值点,当极值点NDVI小于0.4、2个相邻极值点之间的时间间隔大于60 d、振幅大于0.35时,认为这2个极值点之间存在一次完整的作物生长过程,通过以上方法获取2000年、2008年、2020年的耕地复种指数空间分布。
NDVI时序曲线的振幅反映了作物在生长过程中绿叶覆盖的变化情况,NDVI最大值则体现了作物旺盛期生长态势。基于上述预处理后的NDVI时序曲线,提取时序曲线的振幅和NDVI最大值。利用上述样本统计甘蔗NDVI时序振幅和NDVI最大值,绘制直方图。如图9所示,甘蔗NDVI时序曲线的振幅大于0.4,NDVI最大值大于0.6。
图9 甘蔗NDVI时序曲线振幅和最大值直方图

Fig. 9 Histogram of the amplitude and maximum of sugarcane NDVI time series

3.1.4 基于TWDTW模型的甘蔗识别

将研究区甘蔗总体种植面积统计数据除以900(Landsat像元大小为30 m×30 m),获取研究区应有的甘蔗像元总体个数N。筛选出复种指数等于1、NDVI时序振幅大于等于0.4、NDVI时序最大值大于等于0.6的像元,作为可能为甘蔗的备选像元,排除其他不符合该特征的点,例如桉树等常绿林振幅不超过0.4,被预先排除掉。以1月至次年2月为监测时间窗口,利用TWDTW模型逐像元计算NDVI时序曲线与甘蔗NDVI参考时序曲线的“距离”,记作Dist图10),对所有像元的Dist进行排序,按照Dist从小到大的顺序,动态设定Dist阈值,统计Dist取值小于等于该阈值的像元个数N’,当N’恰好等于或刚刚大于总体个数N时,将该阈值设定为甘蔗的Dist阈值。将Dist小于等于该阈值的像元判定为甘蔗类型。
图10 NDVI时序相似性计算结果示例

Fig. 10 Calculation of the similarity of NDVI time series

3.2 甘蔗种植区空间聚集特征分析

本研究利用核密度分析方法分析甘蔗种植区的空间聚集程度。核密度分析方法借助一个移动窗口,对点格局的密度进行估计,其原理[31]如下:以每个栅格样点xi为中心,搜索半径范围h内符合特定属性的点;通过核函数计算出每个点对该栅格的密度贡献值,距离越近,密度值越大;对该半径范围内所有点的密度贡献值累加;输出每个栅格的密度值(式(1))。
f x = 1 n h i = 1 n K x - x i h
式中:f(x)为空间位置x处的核密度计算函数;h为搜索半径范围;n为分析范围内的点数;K为默认的权重核函数;x-xix到点xi之间的距离。

3.3 甘蔗种植斑块景观格局分析

本研究将景观生态学中的景观格局分析方法用于分析甘蔗的空间格局特征,研究2000—2020年甘蔗种植斑块的景观格局变化情况。本研究选取的景观格局指数包括景观比例、平均斑块面积、斑块密度、聚集度指数,景观格局指数通过Fragstat软件计算,各指数含义如下:
(1)景观比例(Percentage of Landscape,PLAND):甘蔗景观占整个景观的面积相对比例;
(2)平均斑块面积(Average Patch area, AP):甘蔗景观面积与甘蔗斑块数量的比值;
(3)斑块密度(Patch Density, PD):单位面积上的甘蔗斑块数,反映斑块破碎化程度和景观空间异质性程度;
(4)聚集度指数(Aggregation Index,AI):基于同类型斑块像元间公共边的长度来计算,数值越大,表示斑块之间距离越近,反映不同斑块个体空间分布的聚集程度。

3.4 甘蔗种植区地形特征分析

地形是区域土壤、气候、水文等自然因素的主导因素,可通过高程、坡度、坡向来表达地形的分异规律。本研究在甘蔗种植分布提取结果的基础上,统计分析甘蔗种植区的地形特征,为甘蔗的耕种管理、甘蔗种植区的选择规划提供科学指导。按照中国农业区划委员会发布的《土地利用现状调查技术规程》[32],将甘蔗种植区的坡度划分为0°~2°平坡区、2°~6°缓坡区和大于6°斜坡区这3个等级,分析甘蔗种植区的坡度等级。

4 结果及分析

4.1 甘蔗种植区空间分布提取及精度验证结果

基于上述研究方法提取粤西地区甘蔗种植空间分布结果如图11所示。
图11 2000—2020年粤西甘蔗种植区空间分布提取结果

Fig. 11 Extracted sugarcane plantation distribution of western Guangdong in 2000-2020

利用样本数据验证2000、2008、2020年的甘蔗提取精度,选取甘蔗种植面积排名前三的雷州市、徐闻县、遂溪县进行精度验证。精度验证结果如表2所示,3年平均甘蔗制图精度为87.62%,总体精度为90.90%。
表2 甘蔗识别精度混淆矩阵

Tab. 2 Confusion matrix of the accuracy of sugarcane identification

年份 区县 真实类型 识别类型 合计 制图精度/% 用户精度/% 总体精度/%
甘蔗 非甘蔗
2000年 遂溪县 甘蔗 190 29 219 86.75 90.48 92.41
非甘蔗 20 432 454 95.15 93.71
雷州市 甘蔗 254 36 290 87.59 83.28 91.08
非甘蔗 51 634 685 92.55 94.63
徐闻县 甘蔗 140 22 162 86.42 83.83 91.01
非甘蔗 27 356 383 92.95 94.18
2008年 遂溪县 甘蔗 191 29 220 86.82 83.77 90.69
非甘蔗 37 452 489 92.43 93.97
雷州市 甘蔗 287 38 325 88.31 84.16 90.85
非甘蔗 54 626 680 92.06 94.28
徐闻县 甘蔗 140 23 163 85.89 85.37 90.29
非甘蔗 24 297 321 92.52 92.81
2020年 遂溪县 甘蔗 202 29 231 87.45 83.13 90.50
非甘蔗 41 465 506 91.90 94.13
雷州市 甘蔗 303 31 334 90.72 82.34 91.04
非甘蔗 65 673 738 91.19 95.60
徐闻县 甘蔗 117 15 132 88.64 81.25 90.25
非甘蔗 27 281 309 90.94 94.93
利用区县尺度的统计年鉴数据对甘蔗种植面积遥感估算结果进行精度验证,验证结果如图12所示。总体精度RMSE为4.8×103 hm2,甘蔗面积统计数据与遥感估算结果的相关系数为0.93。
图12 甘蔗遥感估算面积与统计数据对比

Fig. 12 Comparison between the extracted sugarcane area by remote sensing and statistics data

甘蔗识别的误差来源主要包含以下方面:部分园地、荒地等其他地类与甘蔗NDVI时序曲线形态存在相似性,导致基于TWDTW模型计算的曲线距离较小,从而被误判为甘蔗;此外,云雨天气使得光学时间序列数据稀疏或关键时相点数据缺失,导致观测的作物NDVI时序曲线异常,例如双季稻7—8月处于早稻收割、晚稻播种的关键期,理论上,在该时段内NDVI时序曲线呈现一个低谷特征,但云雨天气导致未观测到该低谷特征,部分水稻被误判为甘蔗。

4.2 甘蔗种植区空间聚集特征分析结果

基于核密度计算的甘蔗种植区空间聚集分布如图13所示。2000—2020年遂溪县甘蔗种植区呈现多中心扩散的趋势,并且逐渐形成了多个高密度种植区。2000年甘蔗种植区主要分布在乌塘镇和岭北镇,2008年甘蔗种植区扩散至城月镇、洋青镇、岭北镇、港门镇,2020年甘蔗种植区进一步扩散至杨柑镇、沙古镇、界炮镇,多个乡镇呈现高密度的甘蔗种植区分布格局。2000—2020年雷州市甘蔗种植区在空间分布上呈现出扩散趋势,甘蔗种植高密度区的空间范围呈现出先扩张后收缩的趋势,在位置上也表现出一定变化。2000年甘蔗种植区主要分布在纪家镇、唐家镇和客路镇,2008年调风镇、雷高镇也出现了甘蔗种植高密度区,2020年甘蔗种植高密度区迁移至英利镇、龙门镇、北和镇。2000—2020年徐闻县的高密度甘蔗种植区从徐闻县西部迁移至北部,甘蔗种植区的空间范围呈现出先扩张后收缩的趋势。2000年甘蔗种植区主要集中在南山镇和城北乡,2008年甘蔗高密度种植区迁移至徐闻县北部的下桥镇和曲界镇,随后2020年甘蔗种植面积下降,以下桥镇局部区域为主,曲界镇零星分布。2000—2020年麻章区的甘蔗种植区均聚集在麻章镇,空间分布上变化较小。廉江市的甘蔗在空间上分布零散,未呈现高密度聚集的甘蔗种植区。
图13 2000—2020年粤西甘蔗种植区核密度空间分布

Fig. 13 Spatial distribution of the kernel density of sugarcane plantation area of western Guangdong from 2000 to 2020

以上甘蔗种植区空间聚集分布变化表明,遂溪县甘蔗产业的发展力度最大,呈现多中心发展和高密度聚集的空间格局。雷州市出现了较明显的种植空间分布调整,甘蔗种植区从北部聚集发展至南部、北部2个聚集中心。徐闻县也发生了显著的甘蔗种植布局调整,种植区从西部迁移至北部,种植面积也有所下降。该甘蔗种植区分布密度的分析结果对于糖厂规划选址具有重要的参考价值。

4.3 甘蔗种植斑块景观格局分析结果

湛江各区县2000、2008、2020年甘蔗景观格局变化情况如表3所示。横向比较各个区县的景观格局特征,遂溪县的甘蔗景观比例最大、平均斑块面积最大、聚集度指数最高,说明甘蔗为当地优势作物类型,呈现聚集化、规模化的种植格局,有利于发展机械化、集约化生产。
表3 湛江2000、2008、2020年甘蔗景观格局指数计算结果

Tab. 3 Calculation results of sugarcane landscape pattern index in Zhanjiang in 2000, 2008 and 2020

区县 年份 景观格局指数
景观比例/% 平均斑块面积/hm2 斑块密度/(个/km2 聚集度指数
遂溪县 2000 16.33 1.14 13.36 61.45
2008 21.78 1.69 12.03 67.13
2020 23.53 2.31 9.51 70.99
雷州市 2000 9.47 0.56 15.79 48.85
2008 12.74 0.82 14.60 62.15
2020 15.05 1.08 13.04 64.98
徐闻县 2000 9.44 0.51 17.47 47.19
2008 10.71 1.00 10.09 64.09
2020 7.61 0.67 10.68 59.35
麻章区 2000 7.49 0.78 9.00 57.31
2008 7.50 0.77 9.11 58.17
2020 6.31 0.47 12.65 47.49
廉江市 2000 1.13 0.20 5.23 26.76
2008 2.07 0.35 5.47 42.15
2020 1.77 0.24 6.90 31.26
吴川市 2000 0.43 0.17 2.36 21.64
2008 0.39 0.18 2.09 23.00
2020 0.58 0.15 3.64 17.06
霞山区 2000 1.32 0.27 4.54 33.18
2008 1.14 0.30 3.61 35.50
2020 0.22 0.13 1.60 12.93
坡头区 2000 0.82 0.21 3.67 27.83
2008 0.71 0.20 3.36 25.23
2020 0.52 0.17 2.86 22.14
赤坎区 2000 0.42 0.16 2.37 21.07
2008 0.35 0.18 1.84 22.49
2020 0.22 0.14 1.57 15.16
虽然雷州和徐闻也为甘蔗生产大县,但甘蔗景观平均斑块面积相较遂溪县偏小,斑块密度较高、聚集度指数较小,说明甘蔗景观斑块破碎程度相对较高,蔗田分布相对零散,有待进一步调整甘蔗种植布局,以便于规模化种植。麻章区、廉江区、吴川市、霞山区、赤坎区、坡头区,甘蔗种植面积较小,平均斑块面积、斑块密度和聚集度指数均处于较低水平,蔗田分布十分零散,不利于发展机械化生产。
从时间维度纵向比较,2000—2020年,遂溪县和雷州市的甘蔗景观平均斑块面积增长、斑块密度降低、聚集指数增长,2020年遂溪县的甘蔗平均斑块面积达2.31 hm2,雷州市的甘蔗平均斑块面积达1.08 hm2,均较2000年增长约1倍,表明遂溪县和雷州市发生了较显著的甘蔗种植布局调整和优化,甘蔗种植向空间分布聚集化的方向发展。在这期间,徐闻县经历了甘蔗种植面积先增长后减少的过程,平均斑块面积先增长后下降、斑块密度先减小后增长、聚集程度先升高后降低。2000—2020年麻章区、廉江市、吴川市、霞山区、坡头区、赤坎区的甘蔗景观格局变化相对较小。

4.4 甘蔗种植区地形特征分析结果

选取甘蔗生产大县遂溪县、雷州市、徐闻县,重点分析甘蔗种植区的地形特征(高程、坡度),探究了2000—2020年甘蔗种植区的空间分布变化与地形的关系(图14图15)。
图14 遂溪县、雷州市、徐闻县甘蔗种植区高程统计分布

Fig. 14 Histogram of the elevation of sugarcane planting areas in Suixi, Leizhou and Xuwen

图15 遂溪县、雷州市、徐闻县甘蔗种植区坡度统计分布

Fig. 15 Histogram of the slope of sugarcane planting areas in Suixi, Leizhou and Xuwen

遂溪县甘蔗大多种植于高程在50 m以下、坡度2°以内的平坡区,一般无水土流失现象,在农业机械化种植方面具有较大优势,有利于降低人工成本。2000—2020年新增的甘蔗种植区大多分布于低高程、地势平坦的地带(图14(a)图15(a))。雷州市甘蔗种植区高程大多小于100 m,超过一半的甘蔗种植于坡度2°以内的平坡区,将近40%的甘蔗种植于2°~6°的缓坡区,发生轻度土壤侵蚀和水土流失的概率较大,受坡度制约,开展甘蔗生产机械化作业的潜力受限,加大了甘蔗种植的人工成本。2000—2020年新增甘蔗种植区大多分布于低高程、地势平坦的地带(图14(b)图15(b))。
徐闻属低丘台地地形,北部地势较高,东、西、南三面沿海倾斜。甘蔗种植区的地形特征相对多样,2000年甘蔗种植区主要分布在高程80 m以内、坡度2°以内的平坡区,2000—2020年发生了较大的甘蔗种植区空间布局调整,2008年、2020年甘蔗种植区大多分布于高程100 m以上、坡度2°~6°的缓坡丘陵区(图14(c)图15(c)),给甘蔗生产机械化作业增加了一定的难度。
经统计,各区县甘蔗在东、西、南、北各个坡向的种植分布比例相差不大。由于南坡太阳辐射强度高于北坡,更加有利于甘蔗生长,建议各区县的甘蔗种植可以适当地向南坡调整,对提高甘蔗生产力具有积极作用。

5 讨论

本文协同遥感与统计数据开展粤西地区甘蔗种植分布提取及时空分析研究。通过构建甘蔗参考NDVI时序曲线并计算遥感时序相似性距离,开展统计数据的空间化,与当前大多数甘蔗种植分布提取的方法不同,该方法协同了遥感与统计数据两类数据。下面对本文研究方法的不确定性和未来的研究方向讨论如下:
(1)基于物候特征的作物识别方法是当前研究热点,构建作物NDVI时序曲线是其中的关键环节,然而云雨天气容易导致观测数据缺失,给甘蔗种植分布提取带来不确定性,如何结合多源光学遥感数据,以及结合光学和SAR遥感数据重建NDVI时序数据,是日后需要解决的关键技术问题。
(2)本研究以NDVI时序曲线作为判别甘蔗的特征,但是部分园地、荒地等其他地类与甘蔗NDVI时序曲线形态存在相似性,导致误判为甘蔗类型,后续研究将开展植被指数的筛选和构建,选择敏感特征,丰富特征维度,有助于提高识别精度。
(3)考虑到数据的可获得性,本研究采用了 30 m分辨率的Landsat数据开展甘蔗历史分布信息的提取,然而30 m的空间尺度对于细小、破碎的地块,不可避免地存在混合像元问题,如何对混合像元进行解混从而减小甘蔗识别的不确定性,有待进一步探究。
(4)本研究以统计数据为约束开展甘蔗种植分布提取,该方法对统计数据具有依赖性。在最新统计数据公布之前,及时提取甘蔗种植分布,对于提高信息时效性具有重要价值。后续可以探究的研究思路如下:先基于本文方法提取上一年的甘蔗种植分布,以甘蔗轮作规律(甘蔗为多年生作物,种植后可连续收获2—3年)为先验知识,结合当前最新年份的少量样本开展变化检测,对甘蔗种植分布数据进行更新。
(5)在甘蔗种植分布时空分析方面,本文采用2000、2008、2020年3个年的甘蔗种植分布提取结果来表征2000—2020年的甘蔗种植时空分布情况,时间间隔相对较长,后续研究将提高信息提取的频次,以年为单位,开展连续年份的甘蔗种植分布变化提取。
(6)本文从空间聚集程度、景观格局特征、地形特征3个方面对甘蔗种植区分布开展时空分析,后续研究将进一步从自然因素、人文因素两个方面,分析影响甘蔗种植空间分布的驱动力,自然因素包括气候变化、土壤条件、水热条件等,人文因素包括人口、经济政策、市场需求等。

6 结论

粤西是我国甘蔗三大优势产区之一,甘蔗种植空间分布及其变化信息可以为生产资源调度、甘蔗种植布局优化、甘蔗产业政策制定等提供数据支撑,对促进我国甘蔗产业可持续发展具有重要价值。遥感为获取大范围、长时序的甘蔗种植分布信息提供了有效手段。本研究基于Landsat NDVI时间序列数据,利用甘蔗生长物候特征,通过TWDTW模型匹配NDVI时序曲线的相似性,结合甘蔗种植面积统计年鉴数据,提取2000、2008、2020年的甘蔗种植空间分布信息。分析甘蔗种植区空间聚集特征,分析2000—2020年甘蔗种植斑块景观格局特征,基于DEM数据进一步分析蔗区地形特征。主要结论如下:
(1)基于Landsat NDVI时序数据,利用TWDTW模型匹配NDVI时序相似性,可以实现较高精度的甘蔗空间分布信息提取,2000、2008、2020年3年平均甘蔗制图精度达87.62%。
(2)粤西甘蔗主要分布在遂溪县、雷州市和徐闻县,2000—2020年遂溪县、雷州市的甘蔗种植分布呈现扩散趋势,徐闻县甘蔗种植分布呈现收缩趋势。2020年粤西甘蔗种植主要聚集在遂溪县、雷州市北部、雷州市南部、徐闻县北部。
(3)甘蔗种植区的空间聚集特征分析结果表明,2000—2020年遂溪县甘蔗种植区呈现多中心发展和高密度聚集的空间格局,雷州市甘蔗种植区从北部聚集发展至南部、北部2个聚集中心,徐闻县高密度甘蔗种植区从西部迁移至北部。
(4) 2000—2020年粤西甘蔗种植面积最大的 2个区县(遂溪县和雷州市)甘蔗景观平均斑块面积增长、斑块密度降低、聚集指数增长,说明这些地区发生了较显著的甘蔗种植布局调整和优化,甘蔗种植正在向空间分布聚集化的方向发展。
(5)比较各个区县的甘蔗种植区地形特征,遂溪县和雷州市大部分甘蔗种植区地势平坦,2000—2020年新增的甘蔗种植区大多分布于低高程、小坡度的地带,有利于发展机械化作业的潜力,徐闻县大部分甘蔗种植在缓坡区,机械化作业条件较差。
[1]
谭宗琨, 欧钊荣, 何燕. 全球蔗糖主产国甘蔗产量与气象条件关系的初步研究[J]. 中国农业气象, 2007, 28(1):71-75.

[ Tan Z K, Ou Z R, He Y. A preliminary study on relationship between sugarcane yields and meteorological conditions in main sucrose production countries worldwide[J]. Chinese Journal of Agrometeorology, 2007, 28(1):71-75. ]

[2]
国家统计局农村社会经济调查司. 中国农村统计年鉴[M]. 北京: 中国统计出版社, 2005.

[Department of Rural Surveys, National Bureau of Statistics. China rural statistical yearbook[M]. Beijing: China Statistics Press, 2005. ]

[3]
马改艳, 徐学荣. 对当前我国甘蔗产业发展形势的分析与思考[J]. 云南农业大学学报(社会科学版), 2013, 7(6):29-35.

[ Ma G Y, Xu X R. Analysis and deliberation on the current development situation of China's sugar industry[J]. Journal of Yunnan Agricultural University (Social Science), 2013, 7(6):29-35. ]

[4]
严艳荣. 粤西地区甘蔗产业对外贸易的发展现状及对策[J]. 桉树科技, 2019, 36(2):57-62.

[ Yan Y R. Development status and countermeasures of foreign trade of sugarcane industry in western Guangdong[J]. Eucalypt Science & Technology, 2019, 36(2):57-62. ] DOI:10.13987/j.cnki.askj.2019.02.011

DOI

[5]
Wardlow B D, Egbert S L, Kastens J H. Analysis of time-series MODIS 250 m vegetation index data for crop classification in the US Central Great Plains[J]. Remote Sensing of Environment, 2007, 108(3):290-310. DOI:10.1016/j.rse.2006.11.021

DOI

[6]
Löw F, Michel U, Dech S, et al. Impact of feature selection on the accuracy and spatial uncertainty of per-field crop classification using Support Vector Machines[J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2013, 85:102-119. DOI:10.1016/j.isprsjprs.2013.08.007

DOI

[7]
Zhong L H, Hu L N, Zhou H. Deep learning based multi-temporal crop classification[J]. Remote Sensing of Environment, 2019, 221:430-443. DOI:10.1016/j.rse.2018.11.032

DOI

[8]
赵英时. 遥感应用分析原理与方法[M]. 北京: 科学出版社, 2003.

[ Zhao Y S. Principle and method of remote sensing analysis[M]. Beijing: Science Press, 2003. ]

[9]
宋茜, 周清波, 吴文斌, 等. 农作物遥感识别中的多源数据融合研究进展[J]. 中国农业科学, 2015, 48(6):1122-1135.

DOI

[ Song Q, Zhou Q B, Wu W B, et al. Recent progresses in research of integrating multi-source remote sensing data for crop mapping[J]. Scientia Agricultura Sinica, 2015, 48(6):1122-1135. ]

DOI

[10]
胡琼, 吴文斌, 宋茜, 等. 农作物种植结构遥感提取研究进展[J]. 中国农业科学, 2015, 48(10):1900-1914.

DOI

[ Hu Q, Wu W B, Song Q, et al. Recent progresses in research of crop patterns mapping by using remote sensing[J]. Scientia Agricultura Sinica, 2015, 48(10):1900-1914. ]

DOI

[11]
唐华俊, 吴文斌, 杨鹏, 等. 农作物空间格局遥感监测研究进展[J]. 中国农业科学, 2010, 43(14):2879-2888.

[ Tang H J, Wu W B, Yang P, et al. Recent progresses in monitoring crop spatial patterns by using remote sensing technologies[J]. Scientia Agricultura Sinica, 2010, 43(14):2879-2888. ]

[12]
Dong J W, Xiao X M, Menarguez M A, et al. Mapping paddy rice planting area in northeastern Asia with Landsat 8 images, phenology-based algorithm and Google Earth Engine[J]. Remote Sensing of Environment, 2016, 185:142-154. DOI:10.1016/j.rse.2016.02.016

DOI PMID

[13]
Han J C, Zhang Z, Luo Y C, et al. The RapeseedMap10 database: Annual maps of rapeseed at a spatial resolution of 10 m based on multi-source data[J]. Earth System Science Data, 2021, 13(6):2857-2874. DOI:10.5194/essd-13-2857-2021

DOI

[14]
Reed B C, Brown J F, VanderZee D, et al. Measuring phenological variability from satellite imagery[J]. Journal of Vegetation Science, 1994, 5(5):703-714. 10.2307/3235884

DOI

[15]
Zhang X Y, Friedl M, Schaaf C, et al. Climate controls on vegetation phenological patterns in northern mid- and high latitudes inferred from MODIS data[J]. Global Change Biology, 2004, 10(7):1133-1145. DOI:10.1111/j.1529-8817.2003.00784.x

DOI

[16]
Gao F, Zhang X Y. Mapping crop phenology in near real-time using satellite remote sensing: Challenges and opportunities[J]. Journal of Remote Sensing, 2021, 2021:1-14. DOI:10.34133/2021/8379391

DOI

[17]
杨颖频, 吴志峰, 骆剑承, 等. 时空协同的地块尺度作物分布遥感提取[J]. 农业工程学报, 2021, 37(7):166-174.

[ Yang Y P, Wu Z F, Luo J C, et al. Parcel-based crop distribution extraction using the spatiotemporal collaboration of remote sensing data[J]. Transactions of the Chinese Society of Agricultural Engineering, 2021, 37(7):166-174. ]

[18]
Vieira M A, Formaggio A R, Rennó C D, et al. Object Based Image Analysis and Data Mining applied to a remotely sensed Landsat time-series to map sugarcane over large areas[J]. Remote Sensing of Environment, 2012, 123:553-562. DOI:10.1016/j.rse.2012.04.011

DOI

[19]
Zhou Z, Huang J F, Wang J, et al. Object-oriented classification of sugarcane using time-series middle-resolution remote sensing data based on AdaBoost[J]. PLoS One, 2015, 10(11):e0142069s. DOI:10.1371/journal.pone.0142069

DOI

[20]
Wang J, Xiao X M, Liu L, et al. Mapping sugarcane plantation dynamics in Guangxi, China, by time series Sentinel-1, Sentinel-2 and Landsat images[J]. Remote Sensing of Environment, 2020, 247:111951. DOI:10.1016/j.rse.2020.111951

DOI

[21]
Jiang H, Li D, Jing W L, et al. Early season mapping of sugarcane by applying machine learning algorithms to sentinel-1A/2 time series data: A case study in Zhanjiang City, China[J]. Remote Sensing, 2019, 11(7):861. DOI:10.3390/rs11070861

DOI

[22]
Xie X C, Yang Y C, Tian Y, et al. Sugarcane planting area and growth monitoring based on remote sensing in Guangxi[J]. Chinese Journal of Eco-Agriculture, 2021, 29(2):410-422. DOI:10.13930/j.cnki.cjea.200419

DOI

[23]
Zheng Y, dos Santos Luciano A C, Dong J, et al. High-resolution map of sugarcane cultivation in Brazil using a phenology-based method[J]. Earth System Science Data, 2022, 14(4):2065-2080. DOI:10.5194/essd-14-2065-2022

DOI

[24]
Maus V, Câmara G, Cartaxo R, et al. A time-weighted dynamic time warping method for land-use and land-cover mapping[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2016, 9(8):3729-3739. DOI:10.1109/JSTARS.2016.2517118

DOI

[25]
Hu Q, Yin H, Friedl M A, et al. Integrating coarse-resolution images and agricultural statistics to generate sub-pixel crop type maps and reconciled area estimates[J]. Remote Sensing of Environment, 2021, 258:112365. DOI:10.1016/j.rse.2021.112365

DOI

[26]
Kluger D M, Wang S, Lobell D B. Two shifts for crop mapping: Leveraging aggregate crop statistics to improve satellite-based maps in new regions[J]. Remote Sensing of Environment, 2021, 262:112488. DOI:10.1016/j.rse.2021.112488

DOI

[27]
Wang Y, Zhang Z X, Zuo L J, et al. Mapping crop distribution patterns and changes in China from 2000 to 2015 by fusing remote-sensing, statistics, and knowledge-based crop phenology[J]. Remote Sensing, 2022, 14(8):1800. DOI:10.3390/rs14081800

DOI

[28]
黄青, 唐华俊, 吴文斌, 等. 农作物分布格局动态变化的遥感监测——以东北三省为例[J]. 中国农业科学, 2013, 46(13):2668-2676.

DOI

[ Huang Q, Tang H J, Wu W B, et al. Remote sensing based dynamic changes analysis of crop distribution pattern—Taking northeast China as an example[J]. Scientia Agricultura Sinica, 2013, 46(13):2668-2676. ]

[29]
《广东农村统计年鉴》编辑委员会. 广东农村统计年鉴[M]. 北京: 中国统计出版社, 2021.

[ Editorial committee of Guangdong rural statistical yearbook. Guangdong rural statistical yearbook[M]. Beijing: China Statistics Press, 2021. ]

[30]
NASA JPL. NASA Shuttle Radar Topography Mission Global 1 arc second. NASA EOSDIS Land Processes DAAC, 2013. DOI:10.5067/MEaSUREs/SRTM/SRTMGL1.003

DOI

[31]
王远飞, 何洪林. 空间数据分析方法[M]. 北京: 科学出版社, 2007.

[ Wang Y F, He H L. Method of spatial data analysis[M]. Beijing: Science Press, 2007. ]

[32]
全国农业区划委员会. 土地利用现状调查技术规程[S]. 1984:9.

[ National Agricultural Zoning Committee. Technical Regulations for Land Use Status Survey[S]. 1984:9. ]

Outlines

/