基于多时相环境卫星的冬前油菜种植面积估算

魏传文，黄敬峰，杨玲波

（1. 中国气象科学研究院灾害天气国家重点实验室，北京100081；2. 浙江大学遥感与信息技术应用研究所，杭州310058）

0 引言

油菜是重要油料作物之一，中国油菜种植面积和产量均占世界20%以上，及时有效地监测油菜空间分布情况具有重要意义［1］。根据油菜品种属性、区域气候特征等因素，中国油菜可划分为冬油菜和春油菜两大产区。其中，冬油菜种植面积及产量占全国总面积和总产量的90%以上，是中国油菜种植的主要品种［2］。冬油菜主要分布于我国长江流域湖北、四川、江苏、安徽、湖南等省份。由于冬季寒潮影响，该区域油菜易受冻害影响，据文献记载，长江流域发生过多次不同程度的油菜冻害［3-6］。监测油菜种植面积是油菜冻害等自然灾害预防和受灾评估的基础。

遥感技术具有大范围、近实时、高空间分辨率的特点，是油菜面积监测的可靠技术手段。目前油菜面积估算主要集中于开花期［7-8］，主要是因为开花期油菜与其它地物有显著的光谱差异，易于区分。如Pan等［8］使用Hyperion高光谱卫星数据估算了油菜开花期的种植面积。与光谱角匹配算法和植被指数（Normalized Difference Vegetation Index，NDVI）方法相比，基于多光谱特征拟合的方法能得到更精确的估算结果（RMSE＜0.06）。除了使用开花期影像，为了区分与冬油麦易混淆的冬小麦，有研究利用其它生育期的影像数据提取油菜种植面积［9-11］，如王凯等利用油菜特有的MODIS-NDVI时间序列变化特征，基于决策树算法估算了湖北省多年的油菜种植面积。然而，其中一些阈值仍然采用了开花期的影像数据。由于油菜冻害常发生于冬季，开花期的油菜面积是遭受冻害之后的面积。因此，为了准确获得油菜受冻害面积，冬前油菜面积估算是油菜冻害面积监测的关键。文章利用高空间分辨率、高时间分辨率的环境与灾害监测预报小卫星（HJ-1A/B）数据，采用决策树方法提取荆州市江陵县2009—2015年（不包括2011—2012年生长季）的冬前油菜种植面积，为准确获取油菜冻害面积奠定基础。

1 材料与方法

1.1 研究区概况

研究区位于湖北省荆州市江陵县，属江汉平原油菜产区。湖北省是长江流域冬油菜的主产区，油菜种植面积一直位于全国前列。通过查阅中国气象数据网［12］中国农业气象灾情旬值数据集，湖北省油菜有冻害记录的年份为1993、1994、2004和2008年，有冻害记录的站点10个，共15条记录。江陵县以两熟制耕作制度为主，包括水稻—油菜、水稻—棉花两种。研究区主要油菜品种为中早熟和中熟甘蓝型油菜，种植方式包括育苗移栽和直接栽培。其中育苗移栽一般在9月5—10日进行育苗播种，并在10月中旬进行移栽；直接栽培一般在9月20日播种，翌年5月上中旬收获。研究区地势平坦，属长江冲积平源和四湖滨湖平原并列地带，海拔高程为25.3～40 m，属亚热带季风性湿润气候，全年日照时数1 827～1 897 h，全年平均气温16～16.4℃，年均降雨量900～1 100 mm。研究区位置如图1所示。

1.2 研究方法

1.2.1 遥感数据及处理

该文选取HJ-1A/B卫星影像为数据源。HJ-1A/B数据具有高空间分辨率（30 m）、重访周期短（2 d）、数据免费等优点，更适合于对作物种植比较分散的区域进行遥感监测研究。HJ-1A/B CCD传感器具有蓝、绿、红、近红外4个波段，其中CCD传感器又包含了CCD1和CCD2。该研究拟对2009—2015年江陵县冬前油菜种植面积进行研究，但由于2011—2012年生长季影像数据质量不高，不能支持面积估算研究，因此剔除该生长季，仅对其中的5个生长季进行分析。影像数据从中国资源卫星应用中心下载，为经过相对辐射校正和系统几何校正的2级卫星数据产品，剔除含云量较高的影像数据，共获得可用影像数据32景，其中2009—2010年生长季共7景，2010—2011年生长季共11景，2012—2013年生长季共5景，2013—2014年生长季共4景，2014—2015年生长季共5景。

图1 油菜面积提取研究区Fig.1 Study area of oilseed rape planting extraction

对影像数据进一步校正处理。绝对辐射定标使用定标系数将DN值转换为辐射亮度值。由于影像数据覆盖范围较大，在大气校正之前，要对影像进行裁剪。采用MODTRAN模型的FLASSH模块进行大气校正。采用Landsat-8 OLI数据和二次多项式模型进行几何校正。均匀选取控制点并采用最近邻法进行重采样，保证影像校正位置误差小于0.5个像元。得到反射率后，计算NDVI指数，用于油菜面积估算。

1.2.2 样本选取

训练样本的选取是监督分类重要步骤，它的选取质量直接影响分类精度。该文采用研究区实地调查和参考高分辨率影像数据相结合的方式选取训练样本。根据江陵县农业生产状况和实地考察结果，将地物分为油菜、小麦、居民点和水体等4个类别，其它地物如大麦、林地和蔬菜等由于种植面积较小，予以忽略。实地考察时间为2014年3月31日和2014年5月22—24日。采用差分GPS对主要地物进行定位，将矢量数据加载到Google Earth高分辨率卫星影像上，观察4类地物在影像上呈现的特征，选取训练和验证样本用于油菜面积提取。为了衡量各地物训练样本之间的光谱可分离性，分别计算之间的J-M（Jeffries-Matusita）距离，值越大表明地物之间的可分离性越高［13］。训练样本的选取结果如表1所示。

1.2.3 估算方法

表1 2009—2015年（不包括2011—2012年）江陵县地物训练样本像元数Table 1 Number of training samples for different land covers from 2009 to 2015（excluding 2011—2012）

冬前油菜面积提取采用决策树算法（Decision Tree Classification，DTC）。决策树算法采用地物属性阈值分割原理对地物进行二叉树或多叉树划分，阈值通过样本训练（自动决策树）或人工经验获取［14］。采用决策树算法进行油菜面积提取的步骤是先将水体和居民点分离，再区分油菜和冬小麦。NDVI值具体到日期，表示为NDVIyyyy-mm-dd。为了对决策树算法提取的冬前油菜面积进行定性验证和对比，采用应用较广的最大似然分类法提取开花期油菜面积。最大似然分类法（Maximum Likelihood Classification，MLC）假定各类分布函数为正态分布，并通过样本数据获取目标及其他类别数据的后验概率，通过贝叶斯公式转换为先验概率，根据像元值将各个像元分类为最可能的地物类别。

1.2.4 评价方法

首先以油菜面积统计数据（来源于2011—2016年《湖北农村统计年鉴》）为参照数据对油菜估算面积进行精度验证，计算相对误差。公式：

式（1）中，AR是统计面积，AC是影像估算面积。

其次，对提取的面积进行空间位置验证。采用混淆矩阵评价遥感影像分类精度。评价指标主要有用户精度（User Accuracy，UA）、生产者精度（Producer Accuracy，PA）、错分误差（Commission Errors，CE）和漏分误差（Omission Errors，OE）。它们的计算公式：

式（2）～（5）中，Pi是每类被正确分类的像元数，Phi是混淆矩阵每行的像元总数，Pvi是混淆矩阵每列的像元总数。

2 结果与讨论

2.1 不同地物NDVI时间序列特征

江陵县主要冬季作物为冬小麦及油菜，其光谱特征较为相似，易于混淆，需着重探讨油菜及冬小麦的NDVI变化特征。考虑到2010—2011年生长季内有充足的影像数据，NDVI时间序列较为完整，该文以2010—2011年生长季为例说明油菜和冬小麦的NDVI特征变化。从图2发现，油菜的NDVI值呈现明显的“M”型特征，从2010年10月NDVI值开始增加，到12月达到第一个极大值，随后开始减小，并在2011年开花期（3月底）达到极小值。油菜角果期（4月）NDVI值又逐渐增大，而在角果期之后，NDVI值又逐渐降低。在2011年1—2月，由于油菜处于苗后期以及蕾苔初期，其NDVI值变化较小。

冬小麦的NDVI特征则较为单一，从2010年10月初开始至11月初其NDVI值先降低后逐渐增加，在2011年4月达到最高值，之后再次逐渐降低；水体的NDVI则一直处于低位，与其它地物区分显著；居民点则由于存在一定的植被混合像元，在2010年10月至12月其NDVI值逐渐降低，而开春（3月初）后其NDVI值又逐渐增加，基本小于0.3。

图2 江陵县2009—2015年5个生长季（不包括2011—2012年）4类地物NDVI时序曲线Fig.2 NDVI time series curves of 4 types of ground objects in 5 growth seasons in Jiangling County from 2009 to 2015（excluding 2011—2012）

从图2可知，研究区主要冬季作物冬小麦和油菜的生育期高度重叠，但是其NDVI时序曲线有一定差异。根据统计年鉴数据可知，油菜播种期一般为9月中下旬和10上旬，冬小麦的播种期一般为10月中下旬到11月上旬。油菜的播种期普遍要早于冬小麦，因此在生长早期油菜的NDVI值增长较冬小麦快，相同时期油菜NDVI值要高于冬小麦；其次，在冬季油菜基本停止生长，而冬小麦则继续生长，在NDVI曲线上表现为油菜的NDVI值保持稳定甚至下降，而冬小麦NDVI值继续增加；最后，在油菜开花期受到油菜花的影响，油菜的绿度下降，NDVI值下降，而冬小麦正处在拔节期，其NDVI值快速增加；而过了开花期后，油菜的NDVI值开始增加随后减小，而冬小麦NDVI则在拔节期后持续降低。因此，冬小麦NDVI峰值的出现要早于油菜峰值出现的时间。

2.2 冬前油菜面积估算

根据江陵县不同地物每个生长季NDVI的变化特征，构建相应的决策树。

（1）2009—2010年生长季：根据2009年10—11月油菜和冬小麦NDVI增大而其他地物NDVI变化较小的特征，设置NDVI2009-11-23-NDVI2009-10-24＞0.04，剔除水体和居民地；同时设置 NDVI2009-10-24-NDVI2009-09-15＞0，NDVI2009-09-15＜0.50，NDVI2009-11-23＞0.6，以此区分油菜和小麦。根据初步分类结果发现有一部分油菜的NDVI2009-10-24-NDVI2009-09-15＜0，因此，再次设置NDVI2009-09-15＜0.55，NDVI2009-11-23＞0.60来提取该部分油菜。

（2）2010—2011年生长季：根据2010年11—12月油菜NDVI增长缓慢甚至不变的特征，通过设置NDVI2010-12-21-NDVI2010-11-12＞-0.05，实现油菜提取；之后，根据2010年10—11月油菜NDVI缓慢增长以及10月底油菜NDVI较高的特征，设置NDVI2010-11-12-NDVI2010-10-04＞-0.03 and NDVI2010-12-21＞0.54，剔除上一步油菜提取结果中的其他地物；最后，考虑到部分油菜存在NDVI2010-11-12-NDVI2010-10-04＜0的情况，设置NDVI2010-10-04的范围为0.4～0.65 and NDVI2010-12-21＞0.5，提取剩余部分油菜。

（3）2012—2013年生长季：首先根据冬季油菜NDVI增长缓慢且变化较小的特征，结合实际影像分析，设置NDVI2013-01-01-NDVI2012-11-26＞-0.05，剔除居民点和水体，提取油菜和冬小麦混合像元。其次，根据油菜在早期生长较冬小麦快的特点，设置NDVI2012-11-26-NDVI2012-10-02＞0.25，NDVI2012-11-26＞0.58，从而剔除冬小麦，完成油菜识别。

（4）2013—2014年生长季：根据越冬前油菜生长较快而冬季油菜生长缓慢乃至停止 的 特 点， 设 置 NDVI2013-11-06-NDVI2013-10-12＞0，NDVI2013-11-19-NDVI2013-11-06＞0.15，NDVI2013-11-19＞0.39，初步识别油菜；其次考虑到部分油菜NDVI存在波动的特点，设置NDVI2013-11-06-NDVI2013-10-12＜0，NDVI2013-11-19＞0.43，NDVI2013-10-12＜0.45，剔除冬小麦像元，完成全部油菜提取。

（5）2014—2015年生长季：根据冬前油菜和冬小麦NDVI缓慢增大的特点，设置NDVI2014-12-19-NDVI2014-11-14＞0，剔除水体及居民地，获取油菜和冬小麦。同时，为了防止部分居民地存在绿地的影响，设置NDVI2014-12-19＞0.40，进一步剔除混淆的居民地。之后根据油菜在冬季生长基本停止的特征，设置NDVI2014-11-14-NDVI2014-10-07＞0，NDVI2014-12-19-NDVI2014-10-07＞0.25，剔除冬小麦，提取油菜。

综上研究根据决策树算法获取了2009—2015年5个生长季的油菜种植区域分布图（图3），结果可知江陵县冬小麦主要集中在西部和南部，而东北部则分布较少。，通过人工构建决策树的方式可以获取大部分的油菜和冬小麦，但是依然存在部分冬小麦和油菜及居民地的混淆情况，难以单纯使用NDVI进行识别，需要结合机器学习、深度学习等新技术进一步完善提取算法。

图3 2009—2015年（不包括2011—2012年）江陵县决策树算法油菜种植面积提取结果Fig.3 Classification results of winter oilseed rape based on decision tree method in Jiangling County from 2009 to 2015（excluding 2011—2012）

2.3 决策树算法提取结果精度评价

2.3.1 基于开花期油菜面积的定性评价

使用油菜开花期的遥感影像，基于最大似然分类方法获取江陵县油菜种植面积，并以此为基础评价决策树算法的分类识别精度［15-16］。基于最大似然分类方法提取的油菜种植面积最终结果如图4所示，从图4可知，对于2009—2010年和2010—2011年油菜生长季，油菜的面积提取结果与决策树提取较为接近，而其它时期，油菜面积提取结果差异较大。考虑到未获得江陵县2012—2013年和2013—2014年生长季油菜开花期的遥感影像，因此使用J-M距离方法计算各个时期影像地物的光谱可分离性，结果表明，2012—2013年生长季的地物识别最佳影像为2013-04-15的影像，而2013—2014年生长季地物识别最佳影像为2013-11-19的影像。此外，2014—2015年生长季油菜开花期的遥感影像存在少量云污染，可能对分类结果精度产生影响。

图4 2009—2015年（不包括2011—2012年）江陵县油菜种植面积最大似然提取结果Fig.4 Maximum likelihood classification results of oilseed rape planting area in Jiangling County from 2009 to 2015（excluding 2011—2012）

2.3.2 基于统计数据的油菜面积提取精度评价

参照2011—2016年《湖北农村统计年鉴》油菜种植面积统计数据，采用相对误差计算方法对遥感监测方法提取的油菜面积进行精度评价，结果如表2所示。从面积的相对误差而言，2009—2010年和2010—2011年油菜生长季中，最大似然分类方法油菜面积提取的误差要小于决策树算法提取误差。而由于2012—2013年、2013—2014年生长季油菜的开花期影像缺失，以及2014—2015年生长季开花期影像存在云污染，使得最大似然分类方法提取的油菜面积误差要高于决策树算法提取误差。总体而言，决策树方法提取的油菜种植面积相对误差都低于15%，表现较好。

表2 2009—2015年（不包括2011—2012年）江陵县油菜种植面积估算结果Table 2 Accuracies validation of oilseed rape planting area estimation from 2019 to 2015（excluding 2011—2012）

续表2

2.3.3 基于空间位置的精度评价

根据获取的样方数据，对两种方法的油菜提取面积进行评估，结果如图5所示，在2009—2010、2010—2011、2012—2013年油菜生长季中，最大似然分类方法提取的油菜面积在生产者精度和用户精度两方面都高于决策树方法提取结果。在2013—2014年生长季油菜面积提取中，决策树方法的生产者精度和用户精度更高。在2014—2015年生长季，虽然遥感影像含云，但是最大似然方法的生产者精度和用户精度更高。总之，最大似然分类方法和决策树方法在生长季油菜面积提取上都具有较好的精度，两者在不同年份表现各有优劣，影响其精度的具体原因有待进一步查证分析。

图5 2009—2015年（不包括2011—2012年）江陵县油菜种植面积估算结果的空间位置精度Fig.5 Accuracies validation of oilseed rape planting area estimation in Jiangling County from 2009 to 2015（excluding 2011—2012）

3 结论

冬前油菜面积提取是油菜冻害面积监测的关键。该文基于多时相HJ-1A/B NDVI时间序列特征，采用决策树方法提取了江陵县冬前油菜面积。冬前油菜面积的用户精度达到80.40%～95.56%，生产者精度达到82.56%～91.43%，相对误差低于15%。精度评价结果证明了该方法是可行的。影响油菜面积提取的主要因素是云污染以及其它地物（特别是小麦）与油菜的光谱相似性。为了降低云污染，增加可利用的影像数量，采用新的卫星数据是解决途径之一，例如静止卫星以及合成孔径侧视雷达数据具有较高的时间和空间分辨率，可以用于油菜面积估算研究。为了减小光谱相似性的影响，可以与其它光谱特征或非光谱特征（比如纹理特征）相结合进行油菜面积的估算研究。此外，人工智能等先进算法的应用为作物面积提取开辟了新的思路。后续研究将采用新数据，基于更多特征以及新型分类算法，提高冬前油菜种植面积提取精度，为油菜冻害面积监测提供更精准的支持。