CN111598045B - 一种基于对象图谱和混合光谱的遥感耕地变化检测方法 - Google Patents

Abstract

本发明提供一种基于对象图谱和混合光谱的遥感耕地变化检测方法，包括：步骤1，获取至少五年同期遥感影像数据；步骤2，根据所述遥感影像数据构建影像的对象图谱；构建的过程包括数据预处理、影像分割及特征提取和影像分类；步骤3，对对象图谱进行混合光谱分析；步骤4，通过叠置分析实现耕地变化进行检测。本发明结合耕地对象图谱的差异特征和混合光谱计算的分异特性，并根据耕地变化类型的演化特征建立起相应的知识规则结构，对变化过程中耕地与其他地物类型的转变进行描述，提高了耕地变化细节的检测精度，进而提高整个遥感耕地变化检测方法的精度，使得新检测方法适应性更好。

Description

一种基于对象图谱和混合光谱的遥感耕地变化检测方法

技术领域

本发明属于遥感影像检测技术领域，具体涉及一种基于对象图谱和混合光谱的遥感耕地变化检测方法。

背景技术

随着经济的快速发展、污染加剧、人口增长等众多原因，全球面临着耕地减少、土壤质量下降等粮食安全问题，因此，保护耕地是推进实现粮食安全过程中的一项重要任务；那么，能够高效且精确地检测、把握耕地状况变化，会对土地资源类型变化的精准把控、土地资源发展规划等都具有重要的指导意义。

传统的耕地变化检测是依赖于耕地监测，这种方式存在着工作量大、成本高、周期长、时效性差等问题，已无法满足耕地动态变化检测对工作效率和精度的需求，而结合遥感技术和地理信息系统技术的耕地变化检测方法正是一种快速、准确、实时有效的耕地变化检测方法。

目前利用遥感技术进行耕地变化检测的方法主要有三种：光谱直接比较法、分类后变化检测法和时间序列变化检测法。

光谱直接比较法是直接对同一区域不同时相遥感影像的光谱信息进行处理并比较来确定变化位置、范围与类型的一种方法。目前常用的直接比较法有影像差值法、对象变化矢量分析法和基于像斑光谱向量相似度法等。通过这些方法能够直接确定变化的位置,避免大范围分类,提高了检测效率,但是只能得到变化和未变化的结果，不能直接获得变化的类型。

分类后变化检测方法是一种对所观测的地理范围内的所有数据进行单相分类，比较分析各阶段之间的差异，得到最后的变化检测结果的方法。该方法显著的特点是，先分类后检测。分类后变化检测法是一种常用的变化检测方法。例如，Abd El-Kawy等使用遥感数据视觉解释与遥感影像的监督分类相结合的方法，利用最大似然法对四期Landsat图像进行分类，为埃及尼罗河西部地区提供了近期和历史的土地利用与覆盖(LULC)地图。随着科技的发展，分类后变化检测方法与机器学习相结合，出现了一些新的算法，如支持向量机(Support Vector Machine,SVM)分类变化检测法、神经网络分类变化检测法、决策树分类变化检测法等。这种方法可以回避不同时相和不同传感器带来的归一化问题，但是其变化检测的结果依赖于分类的精度，而且每个时间点提取的特征是相互隔离的，没有考虑到地表覆盖的时间上的依赖性。

时间序列变化检测方法是一种主要针对长时间序列遥感影像的时间趋势分析并进行变化检测的方法。即对同一地理位置上的像素点构造出相应的时间序列，对特征点的时间维变化进行分析和提取，估计并确认特征点，对子序列进行分割和分类，最终得到时间序列变化检测结果。时间序列变化检测法也有着广泛的用途。例如，Alcantara开发一种使用MODIS归一化植被指数(NDVI)时间序列绘制活动土地和休耕土地的方法，并提供全欧洲范围内第一幅农田废弃和耕种的分布图。这种方法主要用于处理中低空间分辨率图像，能够充分利用和处理地物特征的时间相关性，可以有效地提高时间序列分类的精度，但由于时间序列分类算法的能力有限，不擅长对于细微变化的检测和变化细节的提供，不适合不规则变化的检测。

综上所述，目前关于耕地变化检测技术主要面临以下问题：

(1)耕地变化类型的识别。耕地发生变化的位置能够直接得到确定，对于变化类型的划分和变化类型的识别都需要进一步的申请，在申请和应用中，尤其是对三种耕地变化类型存在较大的“辨”识难度，即，新增耕地、撂荒耕地和退耕耕地。

(2)地表覆被的互动变化。目前关于耕地变化检测技术主要着重于土地利用与覆盖变化层次上的申请，忽略了耕地变化过程中各地物要素和覆被类型的互动变化，缺少对不同耕地变化类型联合起来的讨论。

(3)耕地变化细节的检测。对于实际发生耕地变化的区域地表覆被复杂，样本可分离度低，容易与裸地、草地和林地等地类混分，对于变化细节的检测很大程度上要依赖与影像分类精度的提高。

因此，亟需提出一种新的关于耕地变化的检测方法，以解决现有技术中存在的技术缺陷。

发明内容

本发明的目的是提供一种基于对象图谱和混合光谱的遥感耕地变化检测方法，以克服现有检测精度低的问题。

为了实现上述目的，本发明提供如下技术方案：

一种基于对象图谱和混合光谱的遥感耕地变化检测方法，包括如下步骤：

步骤1，获取至少五年同期遥感影像数据；

步骤2，根据所述遥感影像数据构建影像的对象图谱；构建的过程包括数据预处理、影像分割及特征提取和影像分类；

步骤3，对所述对象图谱进行混合光谱分析；

步骤4，通过叠置分析实现耕地变化进行检测。

进一步的，所述步骤2中的数据预处理为进行辐射定标、大气校正、影像融合和几何校正的处理。

进一步的，所述步骤2中的影像分割采用多尺度分割方法，所述多尺度分割是基于异质性最小算法，通过特征选择获取分割影像，所述特征包括光谱特征、几何特征和纹理特征。

进一步的，在进行多尺度分割处理时，分割参数为50，形状因子为0.3，紧致度因子为0.5。

进一步的，所述步骤2中的影像分类为多层次分类，包括支持向量机、最近邻分类和随机森林。

进一步的，所述步骤2的构建的过程还包括对影像分类结果进行精度评价，所述精度评价方式为混淆矩阵法，能够获得用户精度、生产者精度、总体精度和Kappa系数。

进一步的，通过所述精度评价确定采用随机森林分类得到的分割影像。

进一步的，所述步骤3中将对象图谱中的混合像元分解成多个端元进行建模，得到线性光谱混合模型。

进一步的，所述步骤3中混合像元分解得到的端元为林地、草地和裸地；通过对多个所述端元的丰度计算实现对混合光谱的分析；

线性光谱混合模型为：

式中：R_iλ为第λ波段第i混合像元的光谱反射率；0＜λ≤m，m为波段的数量；C_kλ为第k个端元在第λ波段的光谱反射率；

为端元的数目，f_ki为对应于第i混合像元的第k个端元所占的分量值，又称丰度；ε_iλ为残余误差值。

进一步的，所述步骤4是基于建立的知识结构规则利用叠置分析进行耕地变化检测的判断。

与现有技术相比，本申请的技术方案能够实现的有益技术效果：

本发明结合耕地对象图谱的差异特征和混合光谱计算的分异特性，利用混合光谱指数完成对易混地类的分离，并根据耕地变化类型的演化特征建立起相应的知识规则结构，对变化过程中耕地与其他地物类型的转变进行描述，提出了一种基于对象分类和混合像元分解的耕地变化检测的新方法，充分挖掘耕地变化的要素特征，提高了耕地变化细节的检测精度，进而提高整个遥感耕地变化检测的精度，使得本发明的检测方法适应性更好。

附图说明

构成本申请的一部分的说明书附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。其中：

图1为本发明耕地变化检测方法流程图；

图2为本发明实施例中研究区地理位置及其遥感影像图；

图3为本发明实施例中研究区特征分割过程的分析对比图；

图4为本发明实施例中研究区土地利用分类结果图；

图5为本发明实施例中三种端元波谱响应曲线图；

图6为本发明实施例中在2018年三种端元的丰度图；

图7为本发明实施例中2019年三种端元的丰度图；

图8为本发明实施例中林草光谱解混结果图；

图9为本发明实施例中耕地变化检测流程图；

图10为本发明实施例中五年间耕地变化检测结果图。

具体实施方式

下面将对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明的一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员所获得的所有其他实施例，都属于本发明保护的范围。

下面将参考附图并结合实施例来详细说明本发明。需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。

在本发明中，以我国现有土地覆盖/土地利用分类系统为基础，参考中国科学院系统下《中国1：100万土地利用图》分类体系，将土地利用类型分为耕地、林地、草地、水域、建设用地和未利用土地6个一级类型。

将涉及到耕地变化的区域分为三种类型：新增耕地、撂荒耕地和退耕耕地。新增耕地是指新增加的种植农作物的土地，包括其他土地利用类型转换为耕地类型的区域；撂荒耕地是2年或2年以上没有耕种农作物的区域；退耕耕地指受国家政策影响，耕地转换为草地、林地和湿地的区域。

新增耕地的出现在时间上属于相对快速的过程，具体体现在其他类型土地转换为耕地。撂荒耕地的出现在时间上属于相对缓慢的过程，随着时间变迁，撂荒耕地的状态为：耕地到裸地到草地到林地。由于无人管理，撂荒地体现出形状不规则、较为破碎的特征，且林草混杂，在影像中混合像元较多。退耕耕地的出现在时间上属于可能快速也可能缓慢的过程，具体表现为耕地转换为了林地、草地或湿地。退耕还林的区域包括一些不适合耕作的坡耕地，退耕还草相对于撂荒地中耕地变为草地的过程时间上更加快速。退耕耕地相对于撂荒耕地，退耕耕地有人为管理，在影像中形状会较为规整且像元较为纯净。

由于耕地被撂荒后完全恢复到自然状态可能需要几十年，而在本申请的检测方法中分析的年份跨度相对较小，所提取到的撂荒耕地都是五年内新产生的，可认为本申请的检测方法申请中被撂荒的耕地都处于自然演替的初期，即在土地类型中属于裸地或者是生物量较低的林草混杂地块。

本发明提出了一种基于对象图谱和混合光谱的遥感耕地变化检测方法，能够在利用基于对象分类方法优势的同时，通过混合像元分解后的各端元丰度深入到像元内部。耕地变化检测方法如图1所示，具体包括三部分：对象图谱构建、混合光谱分析和耕地变化检测。

步骤一、获取遥感影像

针对要分析的区域，从数据中获取多年的遥感影像，在本发明的实施例中，选取北京市密云区作为研究区域如图2所示。密云区地理位置上西起东经116°39′33"，东至117°30′25"，东西长69公里；南起北纬40°13′7"，北至北纬40°47′57"，位于北京市东北部，属燕山山地与华北平原交接地，地势自北向西南倾斜,呈簸箕形，属暖温带半湿润大陆性季风气候区。

密云区的耕地主要分布在西南部平原地区，北部山区较少。其中，水浇地和菜地主要分布在西南部平原地区，西部河谷地区也有少量分布，旱地则分布较分散，全区平原、丘陵、山地区均有分布。

密云区是北京市重要地表饮用水源地和生态涵养区，为了守护首都“生命之水”，近些年有关部门采取了一系列退耕禁养等保水措施，另一方面随着城市化进程的推进，其农业产业结构也在不断调整，当地耕地的数量和分布都在发生深刻的变化。

步骤二、数据预处理

遥感影像的数据是北京市密云县的高分一号WFV影像和高分六号WFV影像，来源于中国资源卫星应用中心网站。查询华北农作历，7月和8月是农作物的主要生长阶段，确定影像时间分别为2015年8月15日(GF-1)、2016年8月26日(GF-1)、2017年7月8日(GF-1)和2018年8月20(GF-6)日和2019年7月31日(GF-6)。

GF-6WFV传感器可在八个光谱通道中观察到地球反射的太阳辐射，数据的空间分辨率为16m.与GF-1卫星相比，新增了紫波段、黄波段、植物红边波段1和植物红边波段2，这些新增波段可用于提高植被特别是农作物的检测精度。

由于中国资源卫星应用中心提供的数据为未经过任何处理的原始影像1级产品，需要对其进行辐射定标、大气校正、影像融合和几何校正的处理。

步骤三、构建对象图谱

为了构建耕地对象图谱，需要对影像中“图”信息进行对象化提取和多尺度表达和对与耕地相关的“谱”信息进行反演，从而能够建立起遥感像元波谱与目标几何结构的相互转换关系。

面向对象的遥感图像分析能够基于对象开展影像特征分析和分类处理，实现“对象级”的信息提取，主要由两个阶段组成：影像分割及特征提取和分类。

3.1影像分割及特征提取

影像分割是将影像分成均匀区域的方法，这些区域因其特定的属性如纹理、颜色、形状、大小和灰度等级而不同。本发明采用的是应用最为广泛的eCognition中的多尺度分割方法(Multiresolution Segmentation)。多尺度分割是基于最小异质性算法的自下而上区域合并。从一个种子像素开始，相邻区域合并为对象，代表两个相邻区域之间差异的异质性在这个过程中起着重要作用。图像区域的异质性通常分为两类：形状异质性和光谱异质性，形状异质性又包含光滑度和紧凑度。对应在算法中的参数分别是颜色因子、平滑度因子和紧致度因子。通过多次调整颜色因子和紧致度因子进行分割，然后主要依赖目视判别就能得到分割较好的结果，为了得到最佳的分割效果需要确定合适的参数因子，实现影像对象整体异质性最小。

在本实施例中，利用多尺度分割工具对数据预处理后的2015年至2019年五年的密云县高分影像进行多次分割实验得出结论，将分割参数设置为50、形状因子设置为0.3、颜色因子设置为0.7、平滑度因子设置为0.5、紧致度因子设置为0.5时分类效果更好。

提取地物光谱特征、几何特征和纹理特征等指标，每个对象得到30个特征和指数，对此分类特征进行特征优选，结果显示特征空间5维(5个特征)时类别与类别之间的可分性最大，考虑各个类别的样本和影像不同，因此，不同年份所选的对应的5个特征分别如表1所示。

表1每年影像对象优选特征及描述

特征空间优化是基于各个类别的样本，以及初始特征集，找到类别之间区分的最大平均、最小距离的特征组合，作为分类的最优特征集，从而避免分类过程中盲目使用多种特征所导致的计算量急剧增大、分类精度降低、分类特征冗余等问题。图3分别为2015年到2019年的特征优选结果，主要展示了随着特征数量的变化，样本之间的区分距离的变化。

其中，图3(a)和图3(e)为2015年和2019年的特征优选结果图，图中直线一直上升没有出现变缓的趋势；图3(b)、图3(c)和图3(d)为2016年、2017年和2018年的特征优选结果图，图中特征维数达到5样本区分距离增加的趋势开始变缓，考虑两种情况，确定将特征维数设置为5。根据对图3中各年份优选特征的结构分析，可知图中随着特征维数的增加，直线都呈现了一直上升的趋势，兼顾特征数量与分类效率，得到表1的优选特征组合。

3.2影像对象的分类

基于上述分割后的分割影像，也即构建的新影像，利用监督分类方法对分割的影像对象进行分类。在面向影像对象的分类过程中，可以使用特征优选后的影像特征，影像特征包括像元的光谱信息，以及所获得的对象的相关特征和对象自身的特征，如对象的大小、形态信息和对象之间的空间关系信息等。

本实施例中针对2015年至2019年五年的密云区影像采取监督分类的方法，在构建好的影像对象中选择训练样本，为每年影像对象选取了600个分布均匀且具有代表性的样本作为分类依据。分别采用随机森林RF、支持向量机SVM和最近邻分类KNN分类器对分割对象进行处理，为每年影像对象确定精度最高的分类方法。

遥感影像在面向影像对象分类和变化检测之后，都需要验证和评估其分类和检测的准确性，以便与实际地类和变化情况进行比较，从而判定结果的有效性。混淆矩阵法是评价分析结果准确性的常用方法，可用于用户精度、生产者精度、总体精度和Kappa系数。

其中，总体分类精度反映分类结果的总体准确性，指所有分类类别的正确数之和占所有参考样本总数的百分比；生产者精度反映样本分类结果的准确性；用户精度指某一类别正确分类个数占该类别的分类总数的比例；Kappa系数是指实际地物类别与分类结果类别的匹配程度。用户精度、生产者精度、总体精度和Kappa系数的计算均为现有技术，在此不进行详细赘述。对于Kappa系数，学者Monserud和Leemans提出，小于0.4的值代表较差或非常差的一致性，0.4到0.55的值代表一般一致性，0.55到0.7的值代表良好一致性，0.7到0.85的值代表非常好的一致性，大于0.85的值代表图像间有着极好的一致性。

对2015年RF、KNN和SVM的分类结果进行精度评价，总体精度分别为89％、79％和75％，RF分类器的分类精度要高于其他两个分类器，其结果的混淆矩阵如表2：

表2 2015年影像对象的RF分类精度评价

对2016年RF、KNN和SVM的分类结果进行精度评价，总体精度分别为88％、84％和84％，RF分类器的分类精度要高于其他两个分类器，其结果的混淆矩阵如表3：

表3 2016年影像对象的分类精度评价

对2017年RF、KNN和SVM的分类结果进行精度评价，总体精度分别为89％、84％和81％，RF分类器的分类精度要高于其他两个分类器，其结果的混淆矩阵如表4：

表4 2017年影像对象的分类精度评价

对2018年RF、KNN和SVM的分类结果进行精度评价，总体精度分别为87％、78％和84％，RF分类器的分类精度要高于其他两个分类器，其结果的混淆矩阵如表5：

表5 2018年影像对象的分类精度评价

对2019年RF、KNN和SVM的分类结果进行精度评价，总体精度分别为94％、79％和54％，RF分类器的分类精度要高于其他两个分类器，其结果的混淆矩阵如表6：

表6 2019年影像对象的分类精度评价

综上，从2015年到2019年影像对象不同分类方法的对比分析与验证，采用RF分类器对影像对象进行分类的分类精度更高，因此，本发明中采用RF分类器对影像对象进行分类处理，其分类结果如图4所示，具体，图4(a)为2015年的影像对象分类结果，图4(b)为2016年的影像对象分类结果，图4(c)为2017年的影像对象分类结果，图4(d)为2018年的影像对象分类结果，图4(e)为2019年的影像对象分类结果。

步骤四、混合光谱分析

本发明中耕地变化类型中的撂荒耕地和退耕耕地都呈现出耕地到草地或林地的演变状态，增加了变化类型直接辨别的难度；另一方面，当非林草端元(如：水体、建筑用地、耕地等)与林草端元混淆在一起计算时大大增加了林草端元提取的工作量。

根据遥感信息图谱计算的理论，通过多光谱遥感定量化反演地物目标所蕴含的物理属性及生物化学属性的高精度数值描述，能够从地物波谱特性及其动态演化角度揭示地表特征信息场的分布规律。

本发明为了建立耕地变化相关的谱相特征，提出利用光谱混合分析(SMA)对林草覆被区进行端元分解，即利用线性光谱分解方法将影像对象图谱的混合像元分解成多个端元进行建模，能够得到混合光谱模型，混合光谱模型也可称为纯像元覆盖类型(称为端元)的光谱组合，计算像素内的各土地类型的覆盖率。覆盖率是指端元的丰度，其体现了某像元中各端元的占比，从而实现混合像元的分解。

对林草覆盖区进行端元分解，最终得到林地、草地和裸地3种端元，利用完全约束最小二乘法(Fully Constrained Least Squares，FCLS)计算每个像元内3种端元的丰度，从而利用撂荒耕地的林草混杂的中间状态来区分撂荒状态下的林草地和退耕状态下的林草地。

4.1端元选择

对2018年和2019年的影像数据进行主成分变换，取前四个成分作为光谱空间的内在维度，然后利用分类结果的矢量图掩模掉非林草区域。

基于几何顶点的端元选取，将主成分变换后的前四个成分波段两两组合，构造二维散点图，散点图呈现三角形，三角形的3个顶点为纯净像元，二维散点图的内部为混合像元，其中主成分分析在端元分解前进行，取主成分变换后的前四个主成分进行端元分解，前四个主成分包含原影像信息的90％并去除原有影像波段之间的冗余度。

本发明只对有林草覆盖的端元进行混合像元分解，二维散点图并不呈三角形，因此实际选取过程中只选择散点图凸出区域作为纯净像元“初步筛选”。然后将初步筛选的结果保存为矢量格式，并放置在N维可视化散点图中，选取像元集中区像元作为纯净像元，该过程可以看作是纯净像元“纯化”过程，得到的端元更具有代表性。通过基于几何顶点的端元“初步筛选”和“纯化”，结合影像对比查看，最终得到3种端元，分别为林地端元、草地端元和裸地端元。3种端元波谱反射率曲线如图5所示。

研究区三种地物反射率值有明显不同，其中林地和草地光谱曲线变化趋势相似，波峰、波谷出现的位置基本一致。在可见光范围地物反射率表现为裸地大于林草地，在近红外及植被红边波段反射率表现为林地大于草地。

4.2丰度计算

对于得到的林地、草地和裸地三种端元，利用FCLS计算每个像元内3种端元的丰度，最终得到研究区内林草地中的3种端元丰度混解图如图6和图7所示。图6(a)为2018年林地丰度图，图6(b)为2018年草地丰度图，图6(c)为2018年裸地丰度图，图7(a)为2019年林地丰度图，图7(b)为2019年草地丰度图，图7(c)为2019年裸地丰度图；其组分丰度介于0-1之间，丰度越高，说明该像元内该类型的相应地物覆盖度越高。

通过研究区域的训练样本进行统计分析，确定2018年影像以林地端元丰度G1小于0.7且草地端元丰度G2小于0.6为划分条件提取得到该年份的林草混杂区域，确定2019年影像以林地端元丰度G1小于0.8且草地端元丰度G2小于0.8为划分条件提取得到该年份的林草混杂区域，林草光谱解混结果如图8所示，具体图8(a)为2018年林草光谱解混结果图，图8(b)为2019年林草光谱解混结果图。

线性光谱混合模型(LSMM)是混合像元分解的常用方法。它被定义为：像元在某一光谱波段的反射率(光谱亮度值)是由构成像元的端元(基本组分单元)的反射率以其所占像元面积比例为权重系数的线性组合。可表达为

式中：R_iλ为第λ波段第i混合像元的光谱反射率；f_ki为对应于第i像元的第k个端元所占的分量值，又称丰度；C_kλ为第k个端元在第λ波段的光谱反射率；ε_iλ为残余误差值，是指实际观察值与模型拟合值之间的差值；n为端元的数量；m为波段的数量，1≤λ≤m。

通过残余误差值ε_iλ或用均方根误差RMS的方法来评价分解结果，RMS方法的计算公式如下：

LSMM模型从混合像元R_iλ中分离出各端元的平均光谱响应C_kλ，通过求解线性方程来反解各端元在像元中所占的面积比例f_ki，从而将所有像元分解成这些端元的丰度。模型计算的结果表现为各端元的丰度(图像)和以均方根误差表示的残余误差图像。RMS越小，效果越好。

步骤五，耕地变化检测；

分类后的变化检测基于对多期分类矢量结果的比较，能够最小化大气变化和传感器差异带来的影响。这个过程是基于所开发的知识规则而自动化实现的，知识规则结构的建立需要依据所研究的耕地变化类型的演化特征。然后利用叠置分析进行耕地变化检测的判断，具体是将两层地图要素进行叠置，按照指定的规则产生一个新要素层的操作。

以下以研究变化年份跨度5年为例，即起始年份为第i年，终止年份为第i+4年，给出知识规则的结构和耕地变化检测的流程如图9所示。

“新增耕地”的知识规则结构：如果第i年分类层中的“非耕地”是第i+4年分类层中的“耕地”，则对应的第i年至第i+4年间的耕地变化类型为新增耕地。

“退耕耕地”的知识规则结构：如果第i年分类层中的“耕地”是第i+4年分类层中的“水体”，则对应的耕地变化类型为退耕还湿地；如果第i年分类层中的“耕地”是第i+4年分类层中经过混合光谱分解和“提纯”后的林地和草地，且第i+1年至第i+3年间该区域未转变为建筑用地，则对应的耕地变化类型为退耕还林还草。

“撂荒耕地”的知识规则结构：如果第i年分类层中的“耕地”是第i+4年分类层中的林草混杂或裸地和第i+3年分类层中的林草混杂或裸地，且第i+1年至第i+2年间该区域未转变为建筑用地。

其中，林地和草地指的是对应年份光谱分解后提纯得到的林地和草地，其他类型中包括了部分退耕还湿地的区域，流程图主要是为了体现撂荒耕地和退耕还林还草区域的区别，所以退耕还湿类型没有被特别标出。

根据建立的耕地变化类型知识规则结构，对影像分类和端元分解后的结果进行耕地变化检测，计算变化检测后所得斑块的平均大小，将变化检测到的斑块大小设置为大于等于0.0005km2。分类后耕地变化检测结果如图10所示。

在研究区内随机布设130个样点作为真实验证数据，参考天地图北京中的高分辨率影像，统计耕地的真实变化类型，以此建立误差矩阵，如表7所示，对变化检测结果进行精度验证。

表7耕地变化检测误差矩阵

从精度评价指标来看，本次变化检测整体精度达到82.31％，Kappa系数为0.7753，说明检测效果较好。

对耕地变化进行分析：

研究统计了2015年至2019年间新增耕地、退耕耕地和撂荒耕地3种与耕地相关的变化类型，新增耕地在变化类型中占有最重的比例，面积达到了147.9602km2，具体主要是由于林地和未利用土地转变为了耕地，二者共占总新增耕地面积的85.79％；退耕耕地面积为119.7531km2，其中退耕还林占比55.33％；撂荒耕地面积为86.2996km2，其中96.63％处于林草混杂的状态。

表8研究区五年间耕地各变化类型面积统计表

由表8可知，2015年至2019年五年间，密云区的耕地在各地类间均发生着转变。

(1)根据当地实施的耕地补偿制度，新增耕地来源于土地的整理和土地开发，由研究结果可知五年间新增耕地主要来源于林地和未利用土地，其中来源于林地的新增耕地主要发生在山地丘陵区域，来源于未利用土地的新增耕地主要分布于建成区周边的平地，并且新增耕地的面积整体上小于同期耕地减少的面积，即密云区在五年间耕地面积呈现衰减的趋势。

(2)密云区作为北京市重要的生态屏障和水源地，近些年来推出实施退耕禁种等一系列的治理工程，要求密云区特别是水库上游地区的土地利用必须以生态维持和水源保护为核心，对土地利用的类型、方式、强度和布局有极为严格的环境要求和限制，研究结果显示了五年间退耕还湿地和还草主要集中于密云水库周边；密云区境内80％以上的用地为各级水源保护区，为了确保密云水库的水质清洁，其土地利用方式逐步改变为以水土保持和水源涵养林为主地土地利用结构，研究结果显示退耕还林的区域分散在四周的山地区域。

(3)研究显示撂荒耕地主要处在林草混杂的状态，处于裸地状态的较少，且撂荒地多分布于山区丘陵，平原分布较少。分布于平原区处于裸地状态的撂荒多为主动撂荒，具备复耕潜休耕地，待地力恢复后继续耕种；分布于山区丘陵的撂荒多为被动撂荒，有很大概率会持续撂荒。耕地撂荒的主要原因包括：水土流失，造成土地生产能力低、生产条件差不再适宜于耕种，从而撂荒；山区机械化种植的难度大，农作物产量降低，导致撂荒地增多；区域人口大量迁出、村庄衰落，导致农业完全废弃。

本申请结合了耕地对象图谱的差异特征和混合光谱计算的分异特性，提出一种高效、高精度的耕地类型变化检测方法。该方法首先基于遥感影像的面向对象进行分类，获取不同耕地对象的图谱差异特征；然后计算并确定混合光谱分异指标，并将其纳入耕地变化类型识别体系；最后综合耕地变化类型转移矩阵、不同耕地类型端元丰度值、完全约束最小二乘法高效和高精度地辨识三类耕地，实现耕地变化检测。

结果表明，基于耕地对象图谱的差异特征和混合光谱计算的分异特性的耕地变化检测方法是一种高效、高精度的方法。申请区耕地的整体检测精度达到82.31％，而且还具有以下优势：(1)在进行面向对象分类时，综合使用RF、SVM和KNN分类器后发现，RF作为一种非参数的机器学习方法，能够在提高分类精度方面提供很大的帮助。(2)混合光谱分解的端元指数能够对于混合地物提供直接的解释性，能够高效地对混合地物进行分类，在变化检测方面具有很大的潜力。(3)综合各地物要素和覆被类型的互动变化，建立起相应的知识规则结构，能够充分挖掘耕地变化的要素特征。

本申请的研究区在五年间耕地面积呈现衰减的趋势，新增耕地在变化类型中占有最重的比例。该方法能够更高效、快速地识别耕地变化，可望成为耕地变化检测申请的一种重要方法。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均在本发明待批权利要求保护范围之内。

Claims (5)

1.一种基于对象图谱和混合光谱的遥感耕地变化检测方法，其特征在于，包括如下步骤：

步骤1，获取至少五年同期遥感影像数据；

步骤2，根据所述遥感影像数据构建影像的对象图谱；构建的过程包括数据预处理、影像分割及特征提取和影像分类；

所述数据预处理为进行辐射定标、大气校正、影像融合和几何校正的处理；

所述影像分割采用多尺度分割方法，所述多尺度分割是基于异质性最小算法，通过特征选择获取分割影像，所述特征包括光谱特征、几何特征和纹理特征；对所述特征进行优选，得到5维特征空间，对特征空间优化，选择类别之间区分的最大平均、最小距离的特征组合作为分类的最优特征集；在进行多尺度分割处理时，分割参数为50，形状因子为0.3，紧致度因子为0.5；

所述影像分类是利用监督分类方法对特征优选后的影像特征进行分类；

步骤3，对所述对象图谱进行混合光谱分析：

将对象图谱的混合像元分解成多个端元进行建模，得到线性混合光谱模型；所述对象图谱的混合像元分解是仅对林草覆盖的端元进行混合像元分解；混合像元分解得到的端元为林地、草地和裸地；

通过对多个所述端元的丰度计算实现对混合光谱的分析；步骤4，通过叠置分析实现耕地变化进行检测，所述实现耕地变化进行检测是基于建立的知识结构规则利用叠置分析进行耕地变化检测的判断；

所述建立的知识结构规则的结构包括：

“新增耕地”的知识规则结构：如果第i年分类层中的“非耕地”是第i+4年分类层中的“耕地”，则对应的第i年至第i+4年间的耕地变化类型为新增耕地；

“退耕耕地”的知识规则结构：如果第i年分类层中的“耕地”是第i+4年分类层中的“水体”，则对应的耕地变化类型为退耕还湿地；如果第i年分类层中的“耕地”是第i+4年分类层中经过混合光谱分解和“提纯”后的林地和草地，且第i+1年至第i+3年间该耕地未转变为建筑用地，则对应的耕地变化类型为退耕还林还草；

“撂荒耕地”的知识规则结构：如果第i年分类层中的“耕地”是第i+4年分类层中的林草混杂或裸地和第i+3年分类层中的林草混杂或裸地，且第i+1年至第i+2年间该耕地未转变为建筑用地。

2.根据权利要求1所述的基于对象图谱和混合光谱的遥感耕地变化检测方法，其特征在于，所述步骤2中的影像分类为多层次分类，包括支持向量机、最近邻分类和随机森林。

3.根据权利要求2所述的基于对象图谱和混合光谱的遥感耕地变化检测方法，其特征在于，所述步骤2的构建的过程还包括对影像分类结果进行精度评价，所述精度评价方式为混淆矩阵法，能够获得用户精度、生产者精度、总体精度和Kappa系数。

4.根据权利要求3所述的基于对象图谱和混合光谱的遥感耕地变化检测方法，其特征在于，通过所述精度评价确定采用随机森林分类得到的分割影像。

5.根据权利要求1所述的基于对象图谱和混合光谱的遥感耕地变化检测方法，其特征在于，

线性光谱混合模型为：

式中：R_iλ为第λ波段第i混合像元的光谱反射率；0＜λ≤m，m为波段的数量；C_kλ为第k个端元在第λ波段的光谱反射率；

n为端元的数目，f_ki为对应于第i混合像元的第k个端元所占的分量值，又称丰度；ε_iλ为残余误差值。

Images