CN104794424A

CN104794424A - 一种新的中低分辨率遥感数据复合的耕地识别方法

Info

Publication number: CN104794424A
Application number: CN201410028272.2A
Authority: CN
Inventors: 秦晓俊; 李玉婷; 雷燕飞; 谭遵泉
Original assignee: BEIJING TIANHE DIGITAL TECHNOLOGY Co Ltd
Current assignee: BEIJING TIANHE DIGITAL TECHNOLOGY Co Ltd
Priority date: 2014-01-20
Filing date: 2014-01-20
Publication date: 2015-07-22

Abstract

本发明提供一种新的中低分辨率遥感数据复合的耕地识别方法——物候光谱复合方法(Phenological and Spectral Combination Model，PSCM)，包括以下步骤：步骤一、数据预处理；步骤二、低分辨率地物分布概率图提取；步骤三、物候光谱复合影像(Phenological and Spectral Combined Image，PSC Image)生成；步骤四、选择分类样本；步骤五、利用PSC图像进行耕地识别；步骤六、耕地识别精度评价。本发明通过分析典型地物的物候特征和光谱响应特征，由低分辨率植被指数时间序列数据所反映的植被物候信息提供像元分布概率信息，中分辨率多光谱数据提供像元光谱信息，复合物候、光谱两种信息共同识别耕地，以达到提高耕地识别精度，且保持识别结果为中分辨率的目的。

Description

一种新的中低分辨率遥感数据复合的耕地识别方法

【技术领域】

本发明涉及遥感图像分类领域，特别是一种中低分辨率遥感数据复合的耕地识别方法。

【背景技术】

耕地是农业发展的基础，是保障粮食安全生产的重要资源之一。及时获得准确、详实的耕地面积及其空间分布信息，对合理利用和保护耕地资源，保障粮食安全生产，实现社会、经济、环境和生态协调发展至关重要。同一耕地地块，由于轮作复种方式的不同，在不同时期可以表现为裸地或植被信息；不同耕地地块，由于种植和管理方式不同，在同一时期也可以表现为裸地或植被信息。因此，耕地是年内年际变化最为复杂的一种土地利用／覆盖类型，其光谱会随时间变化而变化，对于遥感识别而言也最为复杂。

目前，耕地遥感识别主要以中分辨率(10-30m)和低分辨率(250——1000m)遥感数据源为主。中分辨率遥感识别耕地通常需要选用作物生长季内单期或多期关键期影像。然而，单期影像不能反映农作物生长动态变化信息，无法良好地表达耕地上地物覆盖强烈的光谱变化特征，耕地识别难度大；多期影像虽然可以良好地表达地表信息变化，充分利用耕地特殊光谱特征识别耕地，但由于受到卫星重访周期、作物物候以及天气等诸多因素影响，大范围适合时相的多期影像数据获取十分困难。低分辨率遥感卫星覆盖范围广，重访周期短，可以及时获取大范围影像数据，为耕地识别提供了数据保证。但由于存在大量混合像元，区域内部详细的耕地空间分布信息难以反映，识别结果精度低且可靠性不高，很难满足实际需求。

因此，复合中、低空间分辨率遥感数据进行耕地识别，借助中分辨率数据丰富的空间细节和低分辨率数据较高的时间维信息，在一定程度上解决区域耕地识别中存在的难题同时提高耕地的识别精度，满足实际应用的需求。但目前的中低分辨率遥感数据复合的耕地识别方法，多以低分辨率数据为主要研究数据，中分辨率数据仅起辅助作用(用于选择训练样本和精度检验样本或求解研究方法的模型参数或阈值，进而指导低分辨率数据分类)，耕地识别结果的空间分辨率仍然维持低空间分辨率，识别精度近似于仅采用低分辨率遥感数据的识别结果，无法与中分辨率耕地识别结果相比较，没有发挥出两种数据的各自优势。

【发明内容】

针对前述方法的各种不足，本发明提出一种新的中低空间分辨率遥感数据复合的耕地识别方法——物候光谱复合模型(Phenological and SpectralCombination Model，PSCM)，具体包括以下步骤：

步骤一、数据预处理；

步骤二、低分辨率地物分布概率图提取；

步骤三、物候光谱复合影像(Phenological and Spectral CombinedImage，PSC Image)生成；

步骤四、利用PSC图像进行耕地识别；

步骤五、耕地识别精度评价。

本发明通过分析典型地物的物候特征和光谱响应特征，复合低分辨率植被指数时间序列数据所反映的地物物候信息及单期中分辨率数据的高空间分辨率、多光谱信息共同识别耕地，达到提高耕地识别精度，同时保持了耕地识别结果的高空间分辨率特征(与所采用的中分辨率数据空间分辨率相同)，为实际应用提供更为准确的耕地空间分布信息。

【附图说明】

图1为地理位置及TM影像；

图2为总体流程图；

图3为TM影像8类地物光谱反射率曲线；

图4为种主要土地利用／覆盖类型的MODIS NDVI时间序列曲线；

图5为低分辨率地物分布概率图提取流程图；

图6为PSC影像生成流程图；

图7为分类结果对比图：(a)TM影像分类结果，(b)复合影像分类结果，(c)、(f)为GE高分影像，(d)、(g)分别是对应(c)、(f)的TM影像分类结果细节图，(e)、(h)分别是对应(c)、(f)的复合影像分类结果细节图。

【具体实施方式】

下面结合附图说明及具体实施方式对本发明进一步说明。

研究区和数据

研究区地处黑龙江省牡丹江市与鸡西市交界处内，地理坐标为北纬45°14′-45°4′，东经130°14′-130°22′，区域大小为50X50km，如图1所示。该地区地处山区，地势中部高、南北低，地貌复杂、地形多样。全区主要粮食作物为水稻、玉米和大豆等一年一熟型作物，耕种模式相对简单，符合研究需要。水稻4月中上旬开始育秧，5月中下旬开始移栽，9月下旬成熟收获；玉米和大豆5月上旬开始播种，分别在9月下旬和10月上旬成熟收获。自然植被从4月中上旬开始进入叶芽开放期，9月下旬进入落叶期。

在本实验中，选取2007年8月24日的一景作物生长关键期的LandsatTM5影像，空间分辨率为30m，条代号为115-029；NASA USGS提供的MODIS16天合成250m空间分辨率NDVI产品构建植被指数时间序列数据集，时间从2007年1月1日至2007年12月31日，共23景影像，条代号H27v04；2007年至2011年的高分辨率Google Earth(以下简称GE)数据作为精度验证辅助参考；2007年8月野外调查获取的137个GPS野外详查样本，辅助训练样本和精度评价样本。

方法流程

本发明研究流程主要包括数据预处理、PSCM构建、耕地识别以及精度评价，具体流程如下图2所示：

1数据预处理

由于耕地有别于其他土地利用／覆盖类型，其在不同时期可以表现出的不同的光谱特征，年内变化十分复杂。因此，本发明以研究区内作物生长物候规律作为主要依据来选取适合时相的中分辨率影像。本发明选取了2007年8月24日的一景作物种植关键期Landsat-5TM影像(除去第6波段热红外波段)进行耕地识别，条带号为115-029，UTM投影，WGS84坐标系。在该时期，农作物生长旺盛，植被特征明显，较易识别。采用ENVI4.7软件对Landsat TM5影像进行辐射定标，大气校正，以减少或消除大气效应对影像光谱反射率的干扰。同时，本发明选取16天最大合成的MOD13(陆地2级标准数据产品)中的NDVI产品构建植被指数时间序列数据，空间分辨率250m时间从2007年1月1日至2007年12月31日，共23景数据，条带号H27v04，Sinusoidal投影，WGS84坐标系。利用Erdas9.2软件将MODIS NDVI时间序列数据投影转换成与TM影像相同的投影，即UTM投影，WGS84坐标系。

2PSCM构建

2.1研究区分类体系建立

研究区作物种植模式为一年一熟类型，基本土地利用／覆盖类型有以下5种类别：耕地、自然植被、建筑用地、水体和未利用地。基于TM影像及高分辨率GE数据辅助，耕地可分为旱田、水田和休耕地三个子类；自然植被可分为向阳植被和背阴植被两个子类(二者光谱反射率曲线差异较大，为保证分类精度，故区分为两个子类)。具体分类类别信息如表1所示。

表1分类体系

Table l Classification scheme

2.2研究区地物光谱特征分析

图3显示了TM影像上8种土地利用／覆盖类型的平均光谱反射率曲线，可以发现，休耕地、建筑用地、水体和未利用地4种地物之间光谱反射率曲线差异十分明显。向阳植被与背阴植被虽然同属于自然植被类别，但是由于地形的影响和光照条件的差异，导致了这两类地物光谱反射率曲线在第4至第7波段差异明显。旱地和水田的光谱反射率曲线形状虽然相似，但在第4至第7波段上，旱地光谱反射率曲线幅度高于水田曲线，两者较易区分。因此，除向阳植被的光谱反射率曲线以外，旱地与其他7种地物的光谱反射率有着较为明显的差异。即在本研究区内，旱地与向阳植被光谱特征十分相似，区分度低，极易混淆。

2.3研究区地物物候特征分析

图4显示了研究区内4种基本土地利用／覆盖类型(耕地、自然植被、建筑用地及水体)的平均NDVI时间序列曲线(未利用地由于占地面积小且分布分散，在MODIS NDVI时间序列数据多以混合像元形式存在，无法提取较纯净端元，因此未选取此类别)。由于研究区作物种植模式为一年一熟类型，因此植被物候曲线理论上讲只有一个波峰。自然植被的生长季从每年4月开始至10月结束，其NDVI曲线生长季跨度最长且波峰相对平坦。农作物的生长季从每年4月底开始至9月底结束。耕地NDVI值在4至5月份缓慢上升，从6月左右开始急剧上升，峰值出现在8月份左右；农作物收割以后，NDVI值又迅速下降。耕地特殊的关键物候期形成了其独特NDVI曲线。从图中可以发现，每种典型地类的NDVI曲线形状和幅度都不相同且差异明显。因此，利用NDVI时间序列曲线可以有效地区分出几种典型地物。

2.3PSCM构建

基于以上分析，本发明提出一种新的中低分辨率遥感数据复合耕地识别方法——物候光谱复合模型(Phenological and Spectral CombinationModel，PSCM)。主要包括低分辨率地物分布概率图提取和物候光谱复合影像生成两个步骤。

2.3.1低分辨率地物分布概率图提取

利用MODIS NDVI时间序列数据所反映的地物物候信息获取研究区内主要土地利用／覆盖类型的分布概率图，流程图如图5所示。

首先以低分辨率植被指数时间序列数据代替地物光谱反射率作为混合像元的分解因子，然后利用中分辨率影像辅助选取合适的端元(端元类别由研究区具体情况决定)，最后采用不带约束条件的线性光谱混合模型进行混合像元，获得基于地物物候信息的低分辨率地物分布概率图。

2.3.2物候光谱复合影像生成

PSC影像的生成PSCM构建的关键步骤，目的是将低分辨率地物分布概率图中包含的主要土地利用／覆盖类型的分布概率信息(通过物候信息得到)与中分辨率数据包含地物光谱反射率信息相复合，使其既具有光谱信息又包含地物物候(空间分布)信息，流程图如图6所示。

首先，需要在TM影像上找出研究区内容易产生“同物异谱”和“异物同谱”现象的土地利用／覆盖类型，仅选取与这些类别相对应的低分辨率地物分布概率图进行物候光谱复合。这是为了避免在非光谱混淆区域低分辨率数据可能带入的误差。由于研究区内只有旱地与向阳植被的光谱反射率混合较为严重，其余地类光谱反射率均能被良好的区分开。然后，将选定土地利用／覆盖类型的低分辨率地物分布概率图重采样为与TM影像相同的空间分辨率。最后，进行物候和光谱信息的复合，即采用多波段叠加的方法将重采样后的地物分布概率图与中分辨率影像复合，生成PSC影像，公式如下所示：

\{\begin{matrix} I = Σ_{i = 1}^{a} F_{i} + Σ_{j = 1}^{m} T_{j} \\ 1 \leq a \leq n \end{matrix}

式中I为PSC影像，共包含a+m个波段；F_i为低分辨率NDVI时间序列数据混合像元分解后所得到地物分布概率图(集)F(f₁，f₂，f₃，...，f_i，ε，i<n)中第i个端元所对应的分布概率图；T_j为中分辨率影像T(t₁，t₂，t₃，...，t_i，j<m)的第j个波段；n为端元类别数；a为参与影像复合的地物类别；m为中分辨率影像波段数。

3耕地识别

耕地识别过程分为训练样本选择和分类器选择两部分。训练样本的质量是影响分类器分类精度的重要原因之一。通过分析PSC影像光谱反射率曲线，以GE高分辨率数据为真值，结合部分野外调查样本，对每类地物(共8类)分别选择200个有代表性的像元作为训练样本。选择支持向量机模型(Support Vector Machine，SVM)，核函数选取径向基函数(Radial BasisFunction，RBF)，对PSC影像进行分类，进而得到耕地识别结果。

4精度评价

为了与PSCM识别结果进行对比，本发明采用同一套训练样本和相同的分类方法对TM影像进行耕地识别，对仅依赖于光谱信息的耕地识别结果(以下简称光谱分类)与复合物候信息和光谱信息的PSCM耕地识别结果进行定性(目视分析)和定量(误差矩阵)对比。

误差矩阵是用来进行精度评价的一种标准格式，可以计算出Kappa系数、总体分类精度(OA)、用户精度(UA)、生产者精度(PA)等参数。由于Kappa与OA在表达分类精度时有等效的效果，本发明只对比分析SVDD分类结果的后三种分类精度。

4.14.1.目视评价分析

图7为光谱分类结果和PSCM结果及局部细节展示。从图7(a)可以发现，光谱分类整体效果一般。虽然休耕地、建筑用地、水体、未利用地和大部分水田被良好的区分出来，但是旱地与向阳植被相互混淆的现象比较严重。图7(c)所示为山区密林地区，从原始TM影像以及GE数据中可以确定，此区域绝大部分为自然植被，几乎没有耕地。但是由于旱地与向阳植被光谱信息十分相似，导致光谱分类将部分自然植被错分为旱地，见图7(d)。图7(f)所示为研究区北部作物集中连片种植地区，同样，由于向阳植被与旱地光谱信息易混淆难以区分，导致光谱分类将此区域大部分旱地错分为向阳植被，见图7(g)。相反，对于加入了物候信息的PSCM来说，其分类结果中绝大部分旱地与向阳植被都被比较准确的区分开，如图7(b)(e)(h)所示。

从整体分类结果或局部细节结果都可以发现，PSCM无论在光谱信息单一的区域(如建筑用地、水体)或是光谱信息易混的区域(如旱地与向阳植被)都可以产生更加准确的分类结果。而光谱分类在光谱信息单一的区域分类效果良好，在光谱信息易混的区域分类结果比较破碎，错分、漏分现象严重。

5.2误差矩阵分析

本发明利用分层采样的方法对每类地物随机选取0.035％的样本点(由于建筑用地和水体样本点过少，因此人工添加了若干样本点)，基于GE高分辨率数据将这些样本进行目视解译(剔除类别不确定的点)。将解译后的样本点与未参与训练样本选择的野外样本点一起作为精度评价的参考真值。利用检验参考真值，对研究区内5种基本土地利用／覆盖类型的光谱分类结果和PSCM分类结果进行混淆矩阵分析。采用总体精度、用户精度、制图精度、错分误差、漏分误差5个指标以及Kappa系数来评价分类结果，

如表2、表3所示。

表2TM影像分类结果混淆矩阵

Table2Confusion matrix for classification of Landsat5 TM data

*Cul-耕地，Veg-自然植被，Built-建筑用地，Water-水体，Und未利用地

表3复合影像分类结果混淆矩阵

Table2Confusion matrix for classification of integrated data

通过混淆矩阵分析可知，分类总体精度由光谱分类的78.78％(表2)提高到PSCM的90.53％(表3)，Kappa系数由光谱分类的0.71提高到复合分类的0.87。对比表2、3可以发现，PSCM在耕地信息识别上具有较强的优势：耕地的漏分误差从光谱分类的30.4％大幅下降到PSCM的10.00％，错分误差从光谱分类的18.31％下降到PSCM的4.26％；耕地的制图精度及用户精度从光谱分类的69.60％和81.69％大幅度提升为复合分类的90.00％和95.74％。同时，自然植被的分类精度也相应提高。

今后需要在以下三个方面进行深入的研究：其一，在可以选用空间和时间分辨率更高的时间序列数据，并根据研究区的实际情况选用其他类型的植被指数，如增强型植被指数(EVI)进行方法研究；其二，尝试使用不同的软分类模型和硬分类模型；其三，深入探讨训练样本质量和数量与PSCM模型之间存在的内在联系

本发明通过分析典型地物的物候特征和光谱响应特征，由低分辨率植被指数时间序列数据所反映的植被物候信息提供像元分布概率信息，中分辨率多光谱数据提供像元光谱信息，复合物候、光谱两种信息共同识别耕地，以达到提高耕地识别精度，并且保持耕地识别结果的空间分辨率为中分辨率的目的。

以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。

Claims

1.一种新的中低分辨率遥感数据复合的耕地识别方法——物候光谱复合方法(Phenological and Spectral Combination Model，PSCM)，其特征在于包括以下步骤：

步骤一、数据预处理；

步骤二、低分辨率地物分布概率图提取；

步骤四、选择分类样本利用PSC图像进行耕地识别；

步骤五、耕地识别精度评价。

2.一种新的中低分辨率遥感数据复合的耕地识别方法——物候光谱复合方法(Phenological and Spectral Combination Model，PSCM)，其特征在于包括以下步骤：

步骤一、数据预处理。对中分辨率多光谱影像数据进行辐射校正以及大气校正，以获取地物真实的光谱信息；对低分辨率植被指数时间序列数据进行滤波操作，并将其投影转换成与中分辨率多光谱影像数据相同的投影。

步骤二、低分辨率地物分布概率图提取。首先以低分辨率植被指数时间序列数据代替地物光谱反射率作为混合像元的分解因子，然后利用中分辨率影像辅助选取合适的端元(端元类别由研究区具体情况决定)，最后采用不带约束条件的线性光谱混合模型进行混合像元，获得基于地物物候信息的低分辨率地物分布概率图。

步骤三、物候光谱复合影像(Phenological and Spectral Combined Image，PSC Image)生成。首先通过分析研究区中分辨率遥感影像上容易产生“同物异谱”和“异物同谱”现象的土地利用／覆盖类型，仅选取与这些类别相对应的低分辨率地物分布概率图进行物候光谱复合。然后，将选定土地利用／覆盖类型的低分辨率地物分布概率图重采样为与中分辨率遥感影像相同的空间分辨率。最后，进行物候和光谱信息的复合，即采用多波段叠加的方法将重采样后的地物分布概率图与中分辨率影像复合，生成PSC影像，如公式3所示：

\{\begin{matrix} I = Σ_{i = 1}^{a} F_{i} + Σ_{j = 1}^{m} T_{j} \\ 1 \leq a \leq n \end{matrix}

式中I为PSC影像，共包含a+m个波段；F_i为低分辨率NDVI时间序列数据混合像元分解后所得到地物分布概率图(集)F(f₁，f₂，f₃，...，f_i，ε，i＜n)中第i个端元所对应的分布概率图；T_j为中分辨率影像T(t₁，t₂，t₃，...，t_i，j<m)的第j个波段；n为端元类别数；a为参与影像复合的地物类别；m为中分辨率影像波段数。

步骤四、对图像进行分类。首先，选择合适的监督分类器。然后通过目视判定的方式选择出适合数量(30p，p代表用于分类遥感影像的波段数目)的纯净像元构成训练样本集合。最后，对利用选定的分类器和训练样本对PSC影像进行监督分类，获得耕地的空间分布范围。

步骤五、利用误差矩阵对耕地识别精度进行评价。利用检验样本集对PSCM耕地识别结果计算误差矩阵，并比较分析采用相同分类器和训练样本对中分辨率遥感影像分类的耕地识别结果的分类总精度、用户精度和生产者精度。