农业生产是国家社会经济的基础，粮食产量对于制定国家和区域社会经济发展规划、确保国家粮食安全和社会稳定、指导和调控宏观种植结构等均有重要意义(任建强 等，2011；吴炳方，2004)。农作物种植面积，尤其是粮食作物种植面积已成为影响我国粮食产量的关键问题。准确、及时、全面的农作物种植种类信息对于各级政府管理生产和制定相关政策至关重要(刘扬 等，2015)。遥感技术具有宏观、迅速的大范围监测能力，可在短时间内获得客观、准确的农情信息，在农作物识别及其面积估算方面已经得到了重要的应用。

目前大多农作物遥感分类研究中，农作物遥感识别特征相对单一，较多地局限于光谱、微波散射等信息，难以充分利用农作物之间的多元化、全方位信息，农作物遥感识别和面积估算精度难以得到提高。多光谱特征是农作物遥感识别应用最广泛和最基本的识别特征，但由于遥感数据空间分辨率、时间分辨率等方面的限制，以及广泛存在的“同物异谱、异物同谱”现象，致使单纯依靠多光谱特征时产生诸多分类混淆，对于优势的大宗作物(如冬小麦)的分类精度基本维持在80%左右，即使综合多时相特征难以提高到90% (贾坤和李强子，2013)。不同的农作物具有不同的高度、孔隙率、覆盖度等几何特征，从而对微波信号表现出不同的散射特性(Chakraborty和Panigrahy，2000)，且在农作物不同的物候期表现出不同的散射特性。但微波散射特征目前主要用于水稻的识别，由于雷达数据噪声的影响，在旱田农作物分类方面，雷达图像对不同农作物虽具有指示作用(Jia 等，2012；Soares 等，1997；丁娅萍，2013)，但农作物分类尚需进一步挖掘。

综合多源遥感信息，使其得到整合互补(陈启浩，2007)，改善单一类型特征或者较少类型特征组合造成的分类精度较低和分类结果破碎化的现象。在农作物遥感识别过程中，依据不同的分类系统、空间尺度、时间尺度、分类精度等需求，光谱特征、时相特征、极化特征、纹理特征、DEM等多种特征变量都可以被选用于农作物遥感分类(Lu和Weng，2007；Wang 等，2010)，改善分类效果。若以上特征均用于遥感识别，一定程度上会造成维数灾难和信息的冗余，降低分类精度和效率(程希萌 等，2016)。特征优选降低维数、使模型泛化能力更强，减少过拟合，同时增强特征和特征值之间的理解。选择并利用合适的特征变量是提高农作物遥感分类精度的有效途径之一(贾坤和李强子，2013)。因此在众多的特征中进行特征选择获取最优的特征和特征组合作为农作物识别的典型识别特征是迫切需要的。特征选择算法较多应用于高光谱数据中(Kaya，2013；Persello和Bruzzone，2009；Zortea和Haertel，2004)。对于农作物分类应用最为广泛的多光谱数据，目前应用较多的特征选择算法是特征递归消除法RFE (Recursive Feature Elimination)(Atkinson 等，2014；Tuia 等，2009)和随机森林算法RF (Random Forest)(杨珺雯 等，2015；姚登举 等，2014)。杨珺雯等人(2015)利用基于随机森林的RF-RFE去除价值低的特征得到最佳波段组合，实现高光谱数据降维，但特征递归消除算法的稳定性很大程度上取决于在迭代时底层所用的模型(Jarvis，2016)。利用随机森林算法进行特征选择可以在有效降低数据维度的同时保证农耕区土地利用类型分类精度(马玥 等，2016)。随机森林存在两个主要问题，一是特征关联导致重要特征有可能得分很低，二是偏向问题致使特征变量类别多的特征较为有利(Jarvis，2016)。单变量特征选择与复杂的多元特征选择方法相比精度同等甚至更高(Drotár 等，2015)。该算法能够很好理解数据及其特点和结构，利于特征排序优选来优化模型。目前单变量特征选择方法在遥感影像分类方面，尤其是农作物遥感识别方面研究较少。

多源特征综合和识别特征优选是提高农作物分类精度的两个重要方面。本文基于多源特征综合和识别特征优选进行以下两个方面的研究：首先，基于随机森林分类算法，综合多时相光谱特征、波段差值特征，植被指数特征和纹理特征对主要农作物类型进行分类，分析不同类型特征及其组合对作物分类精度的影响程度。其次，利用单变量特征选择方法对以上多种特征构建的特征空间进行特征优选，根据变量的重要性得分进行选择，舍弃重要性低的特征，降低变量维度，在应用随机森林分类的基础上，提高农作物识别精度。本文实验结果表明了单变量特征选择方法的实用性和优越性。

本文以江苏省泗洪县为实验区，如图1所示。泗洪县位于江苏省西北部，淮河中游，紧邻中国5大淡水湖之一的洪泽湖，地理坐标为33°08′—33°44′N、117°56′—118°46′E。泗洪县属于东亚季风区，四季分明，雨热同季，光能充足，降水充裕，适合种植水稻和玉米等粮食作物。

图 1 实验区示意图及野外调查样方布设图 Figure 1 Schematic diagram of experimentation area and of quadrat layout in the field investigation

本实验以泗洪县的水稻和玉米作物为研究对象。泗洪水稻和玉米分别是一季稻和春玉米，前者5月中上旬播种，10月中上旬收获；后者是五月中下旬播种，9月中上旬收获。2015年10月9、10日开展了地面调查。均匀选择泗洪县20个村庄作为调查样区，每个样区内布设5个200×200的样方，调查并记录每个样方内作物分布情况，主要包括作物类型、作物种植面积、地块经纬度及拍摄的作物照片。在上述获取的100个样方中，共有812个地块，种植类型主要有水稻、玉米、林地、草地、休耕地及菜地等。野外调查数据为目视解译提供了可靠的先验知识和验证信息。野外调查样方布设如图1所示。

根据实验区遥感数据的质量和主要农作物的物候信息，本文选取了3景高分一号(GF-1)数据，时相分别为2015-07-13、2015-08-03、2015-10-12和1景环境一号(HJ-1A)数据，时相为2016-08-28。数据信息见图2和表1。遥感数据预处理过程包括：辐射定标、云处理、大气校正、几何校正及重采样(将HJ-1数据重采样为16 m)等。首先，利用ENVI5.1中的辐射校正模块对四景影像进行辐射定标，需注意不同传感器对应的增益和偏移参数不同；其次，同样利用辐射校正模块下的FlAASH大气校正模块对辐射定标后的实验影像做大气校正；然后，采用ERADS2010软件的图像几何校正模块对3景GF-1和1景HJ-1A影像进行正射校正。最后，运用多幅图像插值方法对HJ-1A影像进行云像元处理。本实验利用配准后的8月和10月等时相的多景HJ影像与8月28日HJ-1A影像建立空间映射关系，对其有云区域进行插值，达到去云目的。

表 1(Table 1)

表 1 遥感实验数据列表 Table 1 Table of experental data 影像编号 时相 传感器 数据质量 1 2015-07-13 WFV2 良好 2 2015-08-03 WFV4 良好 3 2015-08-28 CCD1 有云 4 2015-10-12 WFV2 良好

表 1 遥感实验数据列表 Table 1 Table of experental data

图 2 遥感影像示意图 Figure 2 Schematic diagram of remote sensing image

本文的辅助数据除了上述的野外调查数据外，还包括2015年现场勾画的泗洪耕地数据、泗洪县及周边县区的矢量边界数据、水稻和玉米的物候历等。

本文以随机森林分类算法为分类器，基于多时相光谱、植被指数和纹理、波段差值等特征及其组合，对泗洪县主要农作物进行识别和提取；利用单变量特征选择法对所有特征构建的特征空间进行特征选择，然后利用随机森林分类器对优选特征组合分类，达到高精度区分主要农作物的目的。图3为文章的总体技术流程图。

随机森林算法是由美国科学家(Breiman，2001)提出的由多棵CART决策树构成的新型机器学习方法，是一种高效的分类器集成策略。因具备较高的分类准确率、较强的抗噪、抗异常值的能力，同时能处理大量数据，且不易出现过度拟合的优点(杨珺雯 等，2015)，随机森林算法已成为目前遥感影像分类算法中炙手可热的分类器。

图 3 总体技术流程 Figure 3 Flow chart of technical process

随机森林算法的基本思想是通过自助法(boot-strap)重采样技术，不断生成训练样本(约为原始数据的2/3)和测试样本(约为原始数据的1/3)，由训练样本生成多个分类器组成随机森林，测试数据的分类结果按分类树投票多少形成的分数而定(王书玉 等，2014)。分类过程首先从给定的训练集通过多次随机的可重复的采样得到多个bootstrap数据集；然后，对每个bootstrap数据集构造一棵决策树，构造是通过迭代将数据点分到左右两个子集中实现的，这个分割过程是一个搜索分割函数的参数空间以寻求最大信息增量意义下最佳参数的过程。最后，在每个叶节点处通过统计训练集中达到此叶节点的分类标签的直方图经验地估计此叶节点上的类分布。这样的迭代训练过程一直执行到用户设定的最大树深度或者直到不能通过继续分割获取更大的信息增益为止(雷震，2012)。

随机森林算法通过Python语言实现，设置两个参数：生长树的数目N1和特征参数的个数N2。本文经过多次实验，N1取值100时，泛化误差OOB最为稳定。N2则是根据各个试验方案中参与分类的特征参数个数而定。

特征变量重要性和特征选择在遥感影像分类占据着重要的地位(Yang 等，2013)。好的特征选择方法能够提升模型的性能，更能帮助我们理解数据的特点、底层结构，这对进一步改善分类模型及算法都有着重要作用。特征选择是利用一系列的规则，得到特征重要程度的相对关系，自动的选择出对分类过程最重要的特征子集的过程。特征选择模型分为过滤式(filter)、封装式(wrapper)和嵌入式(embedded) 3类(Blum和Langley，1997；Liu和Yu，2005)。过滤式特征选择方法评价特征的主要依据是数据内在的特征，独立于分类算法。封装式特征选择方法是将选择算法作为分类算法的一部分，通过特征的分类性能作为特征选择的评价指标。嵌入式特征选择则是上述两种方法的折中。

本文选择的是基于第一种思想的单变量特征选择算法。该选择算法是一种统计测试方法，能够对每一个特征进行测试，衡量该特征和响应变量之间的关系，得到每个特征的重要性得分，去除得分不高的特征。单变量特征选择有多种特征选择方法，其中较为常用的是SelectKBest方法，该方法计算特征间相关性得分，选择重要性排名排在前K的变量。

本实验利用python的数据分析工具包scikit-learn进行特征优选。首先，设置单变量特征选择方法的参数—需提取的特征个数M；其次，应用单变量特征选择方法计算训练样本数据中每个特征的重要性并排序；然后，利用SelectKBest方法选择前K的特征；最后，结合随机森林分类器和测试样本数据进行分类并计算每个M对应的分类精度，选择最高分类精度对应的子集特征，即是优化的特征子集。

像元尺度的分类精度评价采用的主要参数是基于精度检验样本的混淆矩阵(误差矩阵)计算得到的，通过对混淆矩阵建立的各种统计参数进行的。应用混淆矩阵分析的主要参数有总体分类精度、Kappa系数、错分误差、漏分误差、生产者精度、用户精度等目前普遍采用总体分类精度或者Kappa系数评价整个分类图的精度，以条件Kappa系数评价单一类别的精度。

本文利用部分野外调查数据作为验证样本数据，基于混淆矩阵，对各个实验方案进行评价。

本文选取的农作物识别特征包括4类，分别是原始光谱、波段差值、植被指数以及纹理等特征，共76个特征参数(表2)。

表 2(Table 2)

表 2 特征信息列表 Table 2 Characteristic parameter information table 特征变量 特征时相 特征参数 光谱特征 7、8、9、10 R、G、B、NIR 波段差值特征 7、8、9、10 NIR_10-NIR_9、R_10-R_9、G_10-G_9、B_10-B_9 NIR_10-NIR_8、R_10-R_8、G_10-G_8、B_10-B_8 NIR_10-NIR_7、R_10-R_7、G_10-G_7、B_10-B_7 NIR_9-NIR_8、R_9-R_8、G_9-G_8、B_9-B_8 NIR_9-NIR_7、R_9-R_7、G_9-G_7、B_9-B_7 NIR_8-NIR_7、R_8-R_7、G_8-G_7、B_8-B_7 植被指数特征 7、8、9、10 NDVI、PVI、DVI、SAVI 纹理特征 7、8、9、10 variance、entropy、second moment、correlation、dissimilarity 注：数字10，8，7分别代表时相为2015-10-12，2015-08-03，2015-07-13的影像。为避免混淆，数字9则表示2015-08-28的影像。故NIR_10-NIR_7表示2015-10-12和2015-07-13影像的近红外波段差值。

表 2 特征信息列表 Table 2 Characteristic parameter information table

光谱特征：时相为2015-07-13、2015-08-03和2015-10-12的GF-1及2015-08-28的HJ-1A数据的4个波段反射率。

波段差值特征：任意组合两个时相影像，分别对4个波段反射率做对应的差值。本文共有6组影像组合，24个波段差值特征。

植被指数特征：分别提取上述4个月份影像的归一化植被指数(NDVI)、垂直植被指数(PVI)、差值植被指数(DVI)及土壤调整植被指数(SAVI)。

纹理特征：计算4景影像第一主成分的5个纹理特征，方差(variance)、信息熵(entropy)、二阶距(second Moment)、相关性(correlation)和相异性(dissimilarity)。

本实验波段差值、植被指数特征均由ENVI5.1的Band Math工具实现。由ENVI5.1的主成分分析的正向主成分旋转获取影像的第一主成分，并利用灰度共生矩阵(GLCM)获取上述5种纹理。

为了得到最佳的分类方案，本文共设计6组实验，如表3所示。

表 3(Table 3)

表 3 实验方案组合信息表 Table 3 The combination information table of experimental programs 实验序号 组合方法 A 4个时相的4波段光谱特征 B 4个时相的4波段光谱特征+波段差值特征 C 4个时相的4波段光谱特征+植被指数特征 D 4个时相的4波段光谱特征+纹理特征 E 4个时相的4波段光谱特征+波段差值特征+植被指数特征+纹理特征 F 优化特征子集

表 3 实验方案组合信息表 Table 3 The combination information table of experimental programs

本实验利用单变量特征选择方法计算实验E中76个特征变量的重要性得分，再与随机森林分类器结合，在分类过程中选择最高的分类精度对应的特征个数及特征子集。各个特征得分见图4，不同特征个数对应的总体分类精度和Kappa系数如图5。

图 4 特征重要性得分 Figure 4 Scores of feature’s importance

图 5 基于单变量选择法的优选特征数目与分类精度关系(红点表示精度最高的点) Figure 5 Relationship between number of preferred features and classification accuracy based on this algorithm (Red dot represents the point of highest accuracy)

图4中，band1、band2、band3和band4表示GF-1和HJ-1A影像的蓝光、绿光、红光及近红外波段；7_ndvi表示2015-07-13影像的NDVI，10-9_band1则表示2015-10-12影像和2015-08-28影像的蓝光波段做差值，其他依次类推。由图4可知，不同特征的重要性差异较大，有的特征得分低至2.43，有的特征得分高达1216.7。结合图5(a)(b)可知，随着参与分类的特征个数的增加，分类精度处于不断波动的状态，表明分类精度不随特征个数的增加而提高。特征变量个数取值19时，分类的总体精度和Kappa系数均达到最大值，分别是0.9707和0.96，表明优选特征子集中包含19个特征。从基于重要性得分的优选组合得出，19个特征中，光谱特征有5个，植被指数特征有6个，波段差值特征占7个，纹理特征有1个。光谱特征和植被指数特征在分类中的作用较为显著，纹理特征则不明显。通过单变量特征优选后得到的优化特征共有19个，每个特征重要性得分如图6(a)。不同类型特征个数所占比重如图6(b)。

图 6 优选特征重要性得分及不同类型特征所占个数比重图 Figure 6 Scores of preferred feature’s importance and numeric proportion of different kinds of features

6组实验的分类效果图如图7所示。对比原始影像与图7(a)、(b)及(d)可以看到，影像中的地物类型被错分或漏分的情况较多，如影像的右下部分本应是水的区域被分成了玉米或者建筑用地，说明实验A、B、D的特征构建不理想。基于光谱特征和植被指数特征综合的实验C以及基于优选特征组合的实验F的分类效果较好，极大了减少了错分或者漏分的几率，且后者的分类结果优于其他的实验方案。

应用随机森林分类器，对以上A—F 6个实验方案分类。计算不同方案的分类结果的混淆矩阵，选择混淆矩阵中的总体分类精度、Kappa系数、错分误差和漏分误差等4个指标，对比分析各个分类结果间的差异，评价方案的优劣。误差对比图如图8所示。

由图8(a)(b)可知，随着不同类型特征变量的加入，总体精度和Kappa系数的走势相同，但并不呈单一趋势。实验方案A和B的总体精度分别为0.8822和0.8923，Kappa系数则分别是0.8358和0.845，精度基本持平。纹理特征的加入使得分类的总体精度和Kappa系数均有小幅度的提高，其中，总体精度为0.8876，Kappa系数为0.8427。对于加入植被指数特征，分类精度大幅度上升，总体精度较实验A提高了4.5%，Kappa系数则提高7.6%。对于所有特征综合的实验E，总体精度和Kappa系数均得到了提高，但是低于方案C。特征优选后，分类精度和Kappa系数达到了本次实验精度的最高峰，分别是0.9707和0.9603。

由图8(c)(d)可得，对于错分误差，各种类型特征的加入几乎改善了所有地物类别的精度，表明不同特征综合有效地减少了地物错分的概率。在漏分误差方面，波段差值和纹理特征对裸地和其他两种地物类别的精度改善不明显，地物漏分现象未得到有效减轻。除了上述两个类别，其他4类地物识别精度提高。由本实验中错分和漏分误差两种分类指标可知，植被指数特征和优选特征组合对分类贡献显著。

由上可知，光谱特征中加入不同类型特征，分类效果均得到了不同程度的改善，特征按分类贡献度排列依次是植被指数特征、波段差值特征、纹理特征。此外，优选特征组合分类贡献最大。

图 7 各个实验分类效果 Figure 7 Classification result of each experiment

本文共设计6组实验方案，分别应用多时相光谱特征、光谱+波段差值特征、光谱+植被指数特征、光谱+纹理特征、光谱+波段差值+植被指数+纹理特征以及优选特征组合，主要得到以下两个结论。

图 8 误差对比 Figure 8 Error comparison

(1) 不同类型特征的分类贡献。不同类型特征的加入对分类精度均有一定的贡献。根据各个实验的精度可知：在光谱特征的基础上，单独引入波段差值、植被指数或纹理等特征，分类精度均得到了不同程度的提升，其中植被指数特征的贡献最大，分类效果最好，说明植被指数特征在分类中作用显著。当综合所有类型特征时，4个精度指标均有了较大的提高，玉米和水稻区分效果显著，说明多源特征综合有利于分类效果的改善，但是该方案的分类精度低于光谱特征和植被指数特征综合方案，表明特征过多，数据间的相关性太强，一定程度上会降低分类精度。

(2) 单变量特征选择可获取优选特征子集。多种类型特征组合是分类精度提升的重要因素。但该方案不能确定对分类做出较大贡献的具体的特征变量。本文利用单变量特征选择方法计算特征对分类的贡献度，优选出能获得最佳分类精度的特征变量组合。本实验优选特征组合中光谱特征分别是8月影像的蓝、绿、红3个波段，10月影像的红、近红外波段；波段差值特征分别是10月和8月底的红、近红外波段各自对应的差值，10月和8月影像的近红外波段差值，10月和7月影像的近红外波段差值，8月和7月的蓝、绿、红波段各自对应的差值；植被指数特征分别为10月份影像的垂直植被指数、差值植被指数、土壤调整植被指数、归一化植被指数以及8月影像的归一化植被指数和土壤调整型植被指数；纹理植被指数是10月的二阶距。优选特征组合在6组实验中，分类精度最高，在有效地区分耕地和非耕地的同时，玉米和水稻的识别效果最佳。优选特征组合高于其他两两类型特征组合或者单一特征分类结果，原因在于前者是各个类型特征的优选子集的组合，利用了多源信息的优势，使得有用信息最大化；所有特征综合的分类精度低于优选特征组合的主要原因是后者特征多、数据维度高、数据间相关性强，而后者去除了重要性低的特征，保留了贡献度大的特征，降低了维度，减少了特征之间的信息相关。同时，特征优选降低了工作量，减少了运作时间。

此外，从优选的特征子集中不同类型特征参与分类的贡献度可知以下几点。首先，在时相上，10月是适宜的分类时相，其次是8月初。主要是因为该时间点上玉米已经收割，对应影像的主要农作物为水稻，水稻、玉米的区分明显，识别效果较好；时相0803影像上玉米和水稻等作物的光谱信息最为丰富，信息量最大且差异表现明显，易于两者的区分；其次是波段差值较为重要，在优选特征组合中的个数最多。不同时相的波段差值反映了作物的4个波段的光谱时序差异，与作物的生长周期契合，减弱甚至消除“异物同谱”的现象。虽然波段差值特征与光谱特征融合方案的分类精度逊于其他一些方案，但是在优选的各类特征融合分类中作用显著。三是纹理特征的贡献最为薄弱。该分类影像是16 m和30 m的数据，其纹理信息不如高分辨率影像丰富及明显，对分类影响较小。

针对目前农作物遥感识别中存在的问题，综合多时相光谱信息、多时相波段差值、植被指数和纹理信息，利用随机森林算法，实现了主要农作物的识别，克服了农作物遥感识别特征较为单一的不足，说明了多源特征在分类中的重要性；基于单变量特征选择算法，获得特征重要性，与随机森林分类器的结合有效地选取了特征子集，不但有效地降低了数据维度，而且改善了分类效果。可见基于随机森林分类，单变量特征选择方法既保证分类精度又提高了计算效率，是一种可行且高效的特征优选方法。同时为本文实验区秋季主要农作物识别和后续的作物面积估算提供了有效依据。优选特征中包含适于玉米和水稻识别的适宜时相，优选特征以及优选特征组合，可有效避免盲目选择识别特征而带来的工作量大和精度不理想的问题。

本研究的优势在于(1)引入不同时相间影像波段差值。不同作物在不同物候期其光谱信息不一，尤其是在光谱信息丰富的典型物候期，作物间光谱差异更为显著。波段差值则可充分地利用不同时相间光谱的差异性来识别提取作物。(2)运用单变量特征选择。目前有研究应用特征递归消除和随机森林算法进行特征优选，但对单变量特征选择方法的应用甚少，尤其是农作物识别研究。本文实验验证了该方法的可行性。

本研究在以下两个方面还需要在今后的研究中进一步补充：(1)本实验只应用了上述4种类型特征的部分变量，后续的工作将重点放在其他变量组合和优选的研究上。(2)本文的研究对象是苏北地区的主要农作物，本文算法有待于在其他区域或者更大范围以及不同种植结构的实验区内进行验证。