多时相哨兵2号影像提取撂荒地方法
蒋 涛
(甘肃省地图院, 甘肃 兰州 730000)
摘 要：撂荒地有日渐严重的趋势,快速、精准地掌握撂荒地信息,对保护耕地与粮食安全具有重要意义。本文利用2018—2021年不同季节的哨兵2号卫星影像数据,使用eCognition、ENVI等软件,完成了撂荒地植被指数变化检测法、分类后变化检测法和决策树分类法的提取研究。首先对实验区的影像数据进行多尺度分割、优化分割结果和分类得到耕地范围。接着利用植被指数变化检测法、分类后变化检测法,根据耕地范围的变化和归一化植被指数(NDVI)在夏、秋两个季节的变化提取撂荒地。最后利用决策树分类法对夏季的光谱特征、植被指数与纹理特征进行分析,提取撂荒地。实验表明,植被指数变化检测法提取精度为85%,更适用于撂荒地信息提取。
0 引言
耕地撂荒具有威胁粮食安全、导致土地经济可用性下降、破坏乡村景观、降低农民收入、增加村庄边缘化风险等影响。近几十年来,受社会经济冲击与土地利用政策等影响,耕地撂荒已然成为全球普遍出现的现象,造成了耕地资源的浪费。
最早的撂荒地提取方法有实地调查、当面询问农户耕地状态和填写调查问卷等形式,费时费力,且效果不佳。随着科技进步,发展为多种多样的撂荒地提取方法[1],基本可分为基于图像分类法与变化检测法两大类[2]。
欧盟(European Union,EU)协同定位中心(Co-location Centre,ClC)土地调查,基于Landsat数据应用目视解译分类法进行撂荒地提取,总体精度能够达到80%以上。日本在对东京进行土地调查时采用监督分类方法,利用Landsat 专题绘图仪获取的图像数据(thematic mapper,TM),总体的分类精度达到了78.9%～90%,当采用小范围高分辨率或多时相的影像数据时,提取精度显著提高。我国学者就山东省庆云县多时相遥感数据,利用决策树分类法进行提取,分类精度达到了85%,在提取出撂荒地的基础上,将撂荒地复耕的范围区分出来[3]。马来西亚的一些学者利用Landsat卫星三个季节的观测数据,使用面向对象分类法进行撂荒地提取,精度达到了83.3%～97.4%[4],根据地物的光谱信息和纹理信息进行目标地物提取的同时,还会检测地物的几何信息,从而利用多种信息规则约束目标地类。
美国的WITMER F D W 使用多光谱、多时相Landsat TM影像对波斯尼亚东北部因战争所产生的撂荒地用直接变化检测法进行提取,精度并不理想,一些已知的撂荒地并没有被提取出来[5] 。Kuemmerle利用多个时期的Landsat卫星数据,对喀尔巴千山附近的撂荒地采用分类后变化检测法进行提取,比较分类结果得到变化信息[6] 。采用支持向量机法进行分类后,提取精度能够达到90%,具有极高的分类精度,该方法能够减少环境、传感器、季节以及地形条件等因素对不同时期影像的影响,但是变化检测存在误差累积[7] 。我国学者利用归一化植被指数(normalized difference vegetation index,NVDI)对广东省惠东县撂荒地使用植被指数变化检测法进行提取,提取精度能够达到90%[8] 。我国学者史铁丑2015年利用两年的耕地分布图进行叠加分析提取撂荒地,并与谷歌地图的遥感数据进行对比,提取精度为85%[9] ,其采用2011年前后的影像进行提取,而本文使用最新的遥感影像。
1 研究区及数据
1.1 研究区概况
阜新市清河门区,隶属于辽宁省,位于辽宁西部阜新西南部,南北长12.6 km,东西宽17.25 km,全区面积98.54 km2。实验区西部、西南部荒山荒滩较多,中部和东部较平坦,耕地面积达到5.08×107 m2,主要种植玉米、高粱以及小麦等农作物。
1.2 数据来源
在欧洲空间局官网下载哨兵2号影像数据(https://scihub.copernicus.eu/),周期为2018—2021年,每年的春、夏和秋三个季度的卫星数据与其他数据。在地理空间数据云网站(http://www.gscloud.cn/)下载30 m×30 m数字高程模型(digital elevation model ,DEM)。验证样本数据来源于2021年夏季通过全球定位系统(global positioning system,GPS)实地定位得到的15个撂荒地样本和20个耕地样本。
1.3 影像预处理
下载的哨兵2A L1C 数据只经过了正射校正与几何精度校正,因此利用ESA SNAP软件进行大气校正以及辐射定标,并用Sen2Cor插件将L1C数据转换为L2A数据。将L2A数据导入SNAP软件,同时加载实验区矢量边界的SHP文件,进行B2、B3、B4和B8这4个波段组合,并选“Use Vector as Mask”,依据实验区的SHP文件进行影像裁剪。
2 研究方法
2.1 影像分割
本研究使用历年10月份的卫星数据,进行实验区的耕地提取。设置影像中的光谱权重、分割尺寸、形状因子和紧凑度等参数,采用eCognition进行影像的多尺度分割[10],总体技术路线如图1所示。
图1 研究技术路线图
2.1.1 光谱权重确定
找出B2、B3、B4、B8这四个波段中耕地与非耕地差值最大的波段,从而确定多尺度分割时增加对应波段的权重,进而提高耕地提取的准确性。通过对比有作物耕地与无作物耕地的协方差矩阵和相关系数矩阵,能够得到耕地与非耕地差距最大的波段,将B8的光谱权重设为2,其他B2、B3、B4三个波段的光谱权重设为1。
2.1.2 分割尺度确定
将形状因子设为默认值0.1,紧凑度设为0.5,优选最佳分割尺度的范围。对比分割尺度为40、100、200时的分割效果,如图2所示。通过对比LV数据变化率(rate of change,ROC)变化,结果表明,当分割尺度小于40时,虽然各种地块能够较好且细腻地分开,但是分割结果过于破碎,并不适合进行地物的区分。而当分割尺度大于200时,仅有少数地物被分割,不再出现分割零碎的情况,但是却出现了分割的地块中包含多种地物的情况。为此最终选择95作为分割尺度。
(a)分割尺度为40 (b)分割尺度为100
(c)分割尺度为200
图2 不同分割尺度影像结果
2.1.3 形状因子与紧凑度因子的确定
固定波段权重,进行形状因子的确定,紧凑度设置为0.5,对比形状因子为0.1、0.7和0.9的分割结果,如图3所示,选择适宜的形状因子。可见0.1的形状因子的分类结果进行观察发现地物分类比较完全,耕地的最大面积较为恒定,能够很好地区分耕地与非耕地,因此将形状因子赋值为0.1。
(a)尺度95形状因子0.1 (b)尺度95形状因子0.7
(c)尺度95形状因子0.9
图3 不同形状因子分割效果
影像分类的形状异质性由平滑度与紧凑度构成,二者的权重和也同样为1,平滑度比紧凑度权重高的时候,分割对象的边界较平滑度比紧凑度权重低的分割边界光滑。为了进一步研究紧凑度和平滑度对分割边界的影响,确定紧凑度的赋值,固定另外三个参数,对比紧凑度的值0.1、0.4和0.9的分割结果见图4。耕地的最大面积较为恒定,可见紧凑度与平滑度对耕地提取没有影响,因此将二者的权重设为相同,即都是0.5。
(a)紧凑度0.1 (b)紧凑度0.4
(c)紧凑度0.9
图4 不同紧凑度分割结果
2.2 影像分类及耕地信息提取
光谱特征是利用遥感影像进行面向对象分类的基础特征与直接依据,除了“同谱异物”和“同物异谱”等特殊情况,不同的地物会有不同的光谱特征。纹理特征是遥感影像的灰度分布特征,它能够减少因为“同谱异物”和“同物异谱”所造成的分类误差[11]。为更好地提取耕地还需要设置一些特征,如形状特征,包括外形、轮廓指数等特征,能够提高耕地的提取精度。用eCognition进行面向对象分类,通常基于光谱特征、纹理特征以及一些根据提取地物而添加的专题特征。
2.3 分类后变化检测法提取撂荒地
分类后变化检测法就是通过分析不同年份内同一地区、同一季节的耕地变化来提取出撂荒地。在提取出四年耕地的基础上,将其分类结果导出为设备可视性文件(equipment visibility file,EVF)数据,并利用ENVI Classic将EVF数据转换为SHP数据,导入到ArcGIS中。利用叠加分析对相连的两年进行交集取反,得到的数据便是一年内耕地的变化信息,可以初步视为撂荒地。将四年数据分别以2018—2019年、2019—2020年、2020—2021年三组进行交集取反,得到三年的耕地变化量,最后利用叠加分析中的相交,将三年变化中相交的部分提取出来,其可认定为撂荒地,得到研究区基于分类后变化检测法的撂荒地分布。
2.4 植被指数变化检测法提取撂荒地
归一化植被指数(NDVI)公式为
(1)
式中,N即NIR,代表近红外波段;R表示红外波段。相较于自然生长的植被,耕地因为人为的照顾,比如说施肥、除草等农业活动,植物的长势更好,而当农作物收割后耕地范围内的植被将会因收割而消失,因此耕地的植被变化呈现出周期性,在农作物生长时,植被覆盖率高,相应的NDVI也高[12],但夏秋两季耕地的NDVI变化极为明显。
2.4.1 阈值确定
由不同季节待分类地物NDVI差值统计表1,可以发现夏季与春季耕地的NDVI差值最小值为0.29,夏季与秋季耕地的NVDI差值最小值为0.30,因此能够判断出当前地块的NDVI小于0.3的时候,这块地物便有可能是撂荒地[13]。在此基础上,结合之前提取到的撂荒地,通过不断地调整分割图像,使得到的撂荒地全部落入图像中撂荒地区域,经过测试,发现两个时期NDVI差值中撂荒地分割阈值约为0.38,即NVDI差值小于0.38的区域为撂荒地。
表1 不同季节待分类地物NDVI差值统计
注：NDVI夏-NDVI春表示夏季与春季的NDVI差值;NDVI夏-NDVI秋表示夏季与秋季NDVI的差值。
2.4.2 计算植被指数
将夏秋两季的哨兵数据导入到ENVI中,利用之前提取出的耕地SHP数据,进行掩模处理,再利用波段计算设置算法{[float(b1)-b2]/(b1+b2)}-{[float(b3)-b4]/(b3+b4)},为b1赋值为夏季B8波段,为b2赋值为夏季B4波段,为b3赋值为秋季B8波段,为b4赋值为秋季B4波段,计算夏秋两季的NDVI的差值。
2.4.3 撂荒地提取
使用波段计算设置算法将0.38以上的差值设为0,低于0.38的值设为1,并设置导出数据为TIF格式,最后将TIF导入到ArcGIS中,将0值设为无颜色,1值设为灰色,得到基于植被指数变化检测的撂荒地分布。
2.5 基于决策树分类法的撂荒地提取
决策树模型是一个由多个判断节点组成的树形图[14]。在使用模型进行预测时,会根据输入的参数依次在节点进行判断,直到最后一个节点即为计算结束,得到结果[15]。
2.5.1 分类特征确定
单纯地利用光谱特征进行决策树分类的话,提取出的撂荒地结果可能因为撂荒地生长植物等原因导致精度不高。进而需要其他特征进行进一步的细分,使分类更贴近实际。创建决策树利用的特征通常为光谱特征、植被指数特征、纹理特征[16]。
归一化植被指数则是根据上文,夏季耕地NDVI高,而撂荒地因为人为影响,比如说播撒除草剂等农业行为,植被覆盖率较低,一般小于0.3。同时为了降低岩石、水体、人为建筑等地物的影响,在分类的基础上,再次利用夏秋两季的差值进行区分。
信息熵是一种图像信息量的度量,如果图像没有纹理,信息熵接近于零;如果图像具有良好的纹理,则图像的信息熵最大;如果图像中纹理较少,则图像的信息熵较小。春季信息熵普遍高于夏季和秋季影像中的信息熵,而且春季信息熵耕地和撂荒地的分离性高。因此本文以春季的信息熵为纹理特征来进行撂荒地提取。根据表2,发现3×3窗口时耕地最小值等于耕地的最大值,而5×5窗口时耕地最小值要大于耕地的最大值,并不好区分,因此本次研究使用3×3窗口,并以2.20作为分类依据。
表2 春季哨兵2号信息熵统计
2.5.2 加载算法
将归一化植被指数的节点设为0.30,小于30的为撂荒地,大于30的为其他地物,再将夏秋NDVI差值设为0.38,小于0.38的为撂荒地,大于0.38的为其他地物,最后将信息熵设为2.20,小于2.20的为撂荒地,大于2.20的为其他地物。
算法确定后,利用ENVI的Classification模块中的New Decision Tree建立决策树,第一层节点利用归一化植被指数进行判断,即[float(b1)-b2]/(b1+b2),将b1赋值为近红外波段,b2赋值为红波段,当输出为是的时候得到其他地物,当输出为否的时候,建立下一个节点,利用夏秋NDVI差值进行判断,即{[float(b1)-b2]/(b1+b2)}-{[float(b3)-b4]/(b3+b4)},当输出为否的时候得到其他地物,当输出为是的时候,建立下一个节点,利用信息熵进行判断,当输出为否的时候,是其他地物,输出为是的时候,就能够得到撂荒地。将决策树保存,决策树算法就建立完成了。
2.5.3 提取撂荒地
基于上文设置的决策树,使用Classification模块中的Execute Decision Tree选择之前建立的决策树并运行。将分类结果设为Tif格式,并导入到ArcGIS中,利用ArcGIS将撂荒地提取出来。
3 结果分析
3.1 耕地信息提取精度分析
利用混合矩阵进行精度评定,通过对Kappa 系数来判断分类结果的精度,公式为
(2)
式中,总类别数即总列数,用r来表示;在误差矩阵内第i行、第j列上像元数表示为xij,这即是正确分类的数目;第j列与第i行的总像元数量分别表示为x+j与xi+;总的用于精度评价的数量为N。利用面向对象分类法进行耕地提取的总体精度超过90%,Kappa系数达到0.945 3,结果较好。
3.2 撂荒地提取精度分析
三种提取撂荒地方法的分类精度和Kappa系数分别如下:分类后变化检测法的分类精度为65.28%,Kappa系数为0.63,提取精度较高。利用归一化植被指数以及夏秋两季的遥感影像数据以及信息熵对撂荒地进行提取,得到精度为85%,Kappa系数为0.73的分类结果。最后利用决策树分类法针对NDVI与夏秋NDVI差值提取撂荒地,分类精度达到了65.28%,Kappa值为0.67。结果表明,利用植被指数变化检测法提取撂荒地的精度最高,结果可靠性高。
4 结束语
在面向对象分类法提取耕地范围的基础上,依次使用植被指数变化检测法、分类后变化检测法和决策树分类法对撂荒地进行提取,最后通过Kappa系数对精度进行分析,得到结论如下:
(1)面向对象分类法,当分割尺度设为95、紧凑度因子设为0.5、形状因子设为0.1时分割效果最好。分割后进行分类提取耕地,最终提取出了精度97.03%,Kappa值达到0.945 3的耕地范围。
(2)以提取出的耕地范围为基础,依次使用植被指数变化检测法、分类后变化检测法和决策树分类法,利用光谱特征、植被指数与纹理特征,进行撂荒地提取,提取精度分别为85%、73.82%、65.28%,Kappa系数分别为0.73、0.67、0.63。结果表明,植被指数变化检测法比另外两种提取方法的精度更高,结果更可靠。
实际应用中,高精度耕地现状图和退耕还林图等不易获取,难以大规模应用。为此,可以考虑采用更高分辨率的卫星影像或实地核测进行验证。
参考文献
[1] 陈航,谭永忠,邓欣雨,等.撂荒耕地信息获取方法研究进展与展望[J].农业工程学报,2020,36(23):258-268.
[2] 杨通,郭旭东,岳德鹏,等.基于联合变化检测的耕地撂荒信息提取与驱动因素分析[J].农业机械学报,2019,50(6):201-208.
[3] 肖国峰,朱秀芳,侯陈瑶,等.撂荒耕地的提取与分析:以山东省庆云县和无棣县为例[J].地理学报,2018,73(9):1658-1673.
[4] 杨通,郭旭东,于潇,等.撂荒地监测方法与生态影响述评[J].生态环境学报,2020,29(8):1683-1692.
[5] WITMER F D W. Detecting war-induced abandoned agricultural land in northeast Bosnia using multispectral, multitemporal Landsat TM imagery[J]. International Journal of Remote Sensing, 2008,29(13-14): 3805-3831.
[6] 张路,廖明生,盛辉.基于正交变换的多通道遥感影像变化检测[J].武汉大学学报(信息科学版), 2004,29(5):456-460.
[7] 韩圣其.撂荒地遥感监测与空间分布规律研究[D].北京:北京林业大学,2019.
[8]李升发,李秀彬.耕地撂荒研究进展与展望[J].地理学报,2016,71(3):370-389.
[9]史铁丑.重庆山区耕地撂荒的规模及影响因素研究[D].北京:中国科学院大学(中国科学院地理科学与资源研究所),2015.
[10] 李慧,唐韵玮,刘庆杰,等.一种改进的基于最小生成树的遥感影像多尺度分割方法[J].测绘学报,2015,44(7):791-796.
[11]邓凯,杨灿灿,尹力,等.多源数据协同的城市道路提取[J].测绘通报,2021(10):60-66.
[12]王红岩,汪晓帆,高亮,等.基于季相变化特征的撂荒地遥感提取方法研究[J].遥感技术与应用,2020(3) :596-605.
[13]宋宪强,梁钊雄,周红艺,等.基于决策树与时序NDVI变化检测的耕地撂荒遥感监测——以四川省凉山州普格县为例[J].山地学报,2021,39(6):912-921.
[14]陈利燕,张新长,林鸿,等.跨比例尺新旧居民地目标变化分析与决策树识别[J].测绘学报,2018,47(3):403-412.
[15]万子健,李连营,杨敏,等.车辆轨迹数据提取道路交叉口特征的决策树模型[J].测绘学报,2019,48(11):1391-1403.
引文格式：蒋涛.多时相哨兵2号影像提取撂荒地方法[J].北京测绘,2023,37(4):552-557.
基金项目：国家重点研发计划(2018YFF0215006)
作者简介：蒋涛(1987—),男,甘肃兰州人,大学本科,工程师,从事遥感摄影测量相关工作。
E-mail：297080311@qq.com
转自：“测绘学术资讯”微信公众号
如有侵权，请联系本站删除！
