基于多生育期MODIS-NDVI的区域冬小麦遥感估产研究

李 冰,黎世民,周 磊,黄灿辉,王来刚*

(1.河南省农业科学院 农业经济与信息研究所，河南 郑州 450002； 2.郑州澍青医学高等专科学校，河南 郑州 450000； 3.河南省农业遥感监测中心，河南 郑州 450002)

基于多生育期MODIS-NDVI的区域冬小麦遥感估产研究

李 冰1,2,黎世民1,周 磊3,黄灿辉3,王来刚1*

(1.河南省农业科学院 农业经济与信息研究所，河南 郑州 450002； 2.郑州澍青医学高等专科学校，河南 郑州 450000； 3.河南省农业遥感监测中心，河南 郑州 450002)

以河南省开封市为研究区，采用2005—2013年种植区域冬小麦生育期内16 d合成的250 m空间分辨率的归一化植被指数(NDVI)之和，与冬小麦产量数据进行相关分析，筛选得到最佳遥感估产时相，建立单变量回归模型；同时，基于冬小麦生长的各个关键时期NDVI构建多元回归模型，再通过主成分分析方法对多元回归模型参数进行改进，得到新的估产模型；最后使用开封市2014年的产量数据对估产模型进行验证，旨在构建具有较高精度的估产模型，从而更好地指导小麦生产。结果显示，3种估产模型的估产误差均控制在10.55%内，根据3个模型得到研究区内冬小麦最佳产量预测时段为3月下旬，即拔节期；3个模型中，主成分回归估产模型的产量拟合精度最高，达93.12%，具有一定的实用价值。

遥感估产； MODIS； NDVI； 回归模型; 冬小麦

粮食生产作为人类社会可持续发展最基本的物质条件，对于每一个国家和地区都举足轻重[1]。虽然我国每年粮食产量均处于上升趋势，但是由于人口不断增多，耕地面积逐渐减少，粮食问题依然非常严峻。小麦是我国主要的粮食作物，我国是世界上种植小麦面积最大、产量最高的国家。因此，对小麦长势进行监测，并准确地预估小麦产量，具有重要的现实意义，其可为政府相关部门科学地制定粮食政策、调整国家粮食储备以及制定粮食进出口计划提供重要信息[2]。

20世纪70年代以来，由于航天水平不断提高，遥感探测技术迅速发展成一种新兴的科学技术。卫星遥感技术以其探测周期短、覆盖范围大、资料丰富、时效性强、费用消耗低等优势[3]，被普遍应用于小麦的产量估测研究中，取得了大量的研究成果。王东伟[4]、闫岩等[5]通过遥感影像进行叶面积指数(LAI)的反演，建立作物生长模型，从而实现对冬小麦长势的监测。任建强等[6]构建冬小麦产量形成关键期3—5月的小麦干质量与植物净初级生产力(NPP)的关系模型，并对小麦产量进行预测，取得了较好的结果。黎锐等[7]基于多时相Landsat TM 影像得到的归一化植被指数(NDVI)，运用支持向量回归(SVR)方法构建模型，预测冬小麦产量，与常规多元回归方法相比，有效地提高了估产精度。朱再春等[8]利用光谱角聚类法(光谱角制图和K均值结合)划分试验区，通过回归分析建立了冬小麦产量和生育期内16 d合成的增强化植被指数(EVI)之间的关系模型。当前利用NDVI建立回归模型进行作物估产的研究中，多是基于某个特定区域。各个区域的气候条件不同，估产时期的选择也不同，若选择冬小麦生长的典型时相进行分析，信息会有缺失[9]，若选用整个生长关键期的NDVI信息，则忽略了因素之间存在共线性的影响。针对以上问题，利用2005—2013年多时相的中高时空分辨率的TM和MODIS多光谱影像数据，分析开封市小麦生育期内各个时相的NDVI数据与统计年鉴产量数据之间的相关性，选择最佳时期建模，再选取生长期内所有NDVI值与产量进行多元回归，丰富信息数据，最后，利用主成分回归分析方法对模型进行优化，旨在构建具有较高精度的估产模型，从而更好地指导小麦生产。

1 材料和方法

1.1 研究区概况

开封市地处黄河中下游平原东部，太行山脉的东南方向，位于河南省中东部，介于113°52′15"～115°15′42″E、34°11′45"～35°01′20″N，面积约6 444 km2，气候为暖温带大陆性季风气候，雨热同期，降雨集中，四季分明，年均气温为 14.52 ℃，平均降水量为 627.5 mm，多年平均日照时数为 2 267.6 h，年均无霜期为 221 d。开封市自然条件适宜，其所辖的兰考县、尉氏县、杞县等5县都是我国重要的商品粮和小麦生产基地。该区的粮食作物以小麦、玉米为主，实行一年两熟轮作制度。

1.2 数据来源及处理

1.2.1 数据来源 本研究所用的2005—2013年冬小麦关键生育时期，即2月下旬至5月下旬的0.2～0.8的MODIS-NDVI数据从美国航空航天数据中心共享网(https://ladsweb.nascom.nasa.gov/)获得，空间分辨率为250 m。2005—2013年3—5月份的TM影像从美国地质勘测局共享的Landsat 5(http://glovis.usgs.gov/)下载，空间分辨率为30 m。另外，本研究种植区的提取数据来源于河南省行政区划矢量图，开封市地区的冬小麦种植面积和产量数据来源于2005—2014年《河南省统计年鉴》。

1.2.2 数据处理

1.2.2.1 遥感数据预处理 本研究所选用的NDVI数据主要为TERRA/MODIS 的16 d合成植被指数数据产品(MODIS13Q1)，经过大气校正、几何校正等预处理。利用开封市的行政区划矢量图掩膜提取得到研究区内NDVI影像图。利用GIS的空间分析技术剔除NDVI值中的云值及噪声点。冬小麦种植区的提取主要是基于高空间分辨率的TM影像。本研究所用的遥感影像数据均采用Albers投影(正轴等积割圆锥投影)，全部重采样为250 m分辨率。

1.2.2.2 冬小麦种植区提取 通过分析比较研究区域内所有农作物的类别和物候期，研究区冬小麦在春季具有区别于其他农作物的时相特征，故选用覆盖开封市的2005—2013年3—5月的2景TM影像。以开封市土地利用类型图作为参考，运用非监督分类方法对研究区进行分类，再基于4、3、2波段和5、4、3波段进行人工解译分类修正，通过合并相似类别，获得研究区的冬小麦种植区图。

1.2.2.3 冬小麦种植区NDVI的提取 利用ArcGIS 10.1的操作平台，对获得的冬小麦种植区图进行处理，使小麦种植区栅格属性均为1，其余为0，通过与NDVI遥感影像进行叠加分析，得到研究区内冬小麦种植区的NDVI数据，并计算每年冬小麦生长期内该区域介于0.2～0.8的NDVI之和作为研究指标。通过分析研究区内冬小麦的物候期、耕作制度及水热等气候条件，在2—5月小麦进行间作或套种的作物较少，故可直接研究NDVI与冬小麦产量的关系。

1.3 冬小麦遥感估产模型构建

本研究选取2005—2013年开封市冬小麦产量数据与小麦生长关键时期的MODIS-NDVI数据构建模型。首先，将开封市2005—2013年2—5月的16 d合成的MODIS-NDVI影像图与每年的冬小麦种植区图进行叠加，得到冬小麦种植区的NDVI影像，再计算介于0.2～0.8全部NDVI数据之和。通过对每年不同时期的NDVI值与冬小麦产量之间的相关性进行分析比较，选择与产量相关性最高的NDVI指标作为自变量，得到以冬小麦产量为因变量的单变量的回归估产模型；然后基于多元回归分析方法，利用整个关键生育时期的NDVI值，建立估产模型，并运用主成分分析方法对冬小麦生长期内2月下旬到5月下旬NDVI值进行分析，得到优化的综合指标，进而研究综合指标与冬小麦产量的关系；最后，通过开封市2014年的冬小麦产量数据和提取得到的NDVI数据验证所构建模型的精度及稳定性。整体技术路线如图1所示。

图1 冬小麦遥感估产模型构建流程

2 结果与分析

2.1 冬小麦NDVI与产量的相关性分析

对2005—2013年小麦产量和其重要时期的NDVI进行相关性分析(表1)发现，冬小麦产量与其研究区内各个主要时相的NDVI均有一定的相关性。其中，3月下旬NDVI与产量之间的相关性最好，2月下旬次之，即从2月下旬开始，冬小麦产量与NDVI的相关性逐渐升高，3月下旬达到最好，之后相关性逐步下降。

2.2 冬小麦NDVI与产量的回归模型的构建

根据表1相关性分析结果以及开封市的气候条件和冬小麦物候历，选择使用3月下旬NDVI数据(介于0.2～0.8的NDVI之和，x)建立其与冬小麦产量(y)的回归模型,见式(1)。

y=a+bx

(1)

其中，a、b是回归参数。使用2005—2013年冬小麦产量数据和3月下旬拔节期的NDVI进行线性统计回归，得到的线性回归分析结果见式(2)。

y=30.273x+810 478.075 (2)表1 小麦生长关键期NDVI与产量的相关系数

该方程的R2为0.512，拟合效果较好；F为7.342，大于F0.05(5.59),说明该模型显著。冬小麦的生长是一个连续过程，应综合考虑冬小麦各关键生育时期的生长状况对产量的影响。因此，将每年2月下旬、3月上下旬、4月上下旬、5月上下旬NDVI(介于0.2～0.8的NDVI之和)为自变量x1、x2、x3、x4、x5、x6、x7，冬小麦产量为因变量y，建立多元回归方程见式(3)。

y=a+b1x1+b2x2+b3x3+b4x4+b5x5+b6x6+b7x7

(3)

其中，a、b1、b2、b3、b4、b5、b6、b7是多元回归模型的参数，运用最小二乘法对方程中各参数的值进行估计。使用2005—2013年冬小麦产量数据和2月下旬、3月上下旬、4月上下旬、5月上下旬的冬小麦种植区域内0.2～0.8的NDVI之和进行多元回归，得到的多元回归统计模型见式(4)。

y=-35.211x1-1.707x2+107.018x3-42.528x4+5.330x5+9.022x6-16.630x7+553 842.091

(4)

模型的R2是0.976，说明多元回归方程的拟合度为97.6%，具有一定的稳定性。

2.3 冬小麦NDVI与产量主成分回归优化模型的构建

主成分分析是主成分回归分析的关键内容，对原始变量进行线性变换从而得到较少的新的指标，所获得的新指标变量之间相互独立并且具有很高的代表性。对冬小麦关键生育时期NDVI进行主成分分析，并利用多元回归方法对得到的综合变量进行分析得到式(5)。

y=a+b1F1+b2F2+b3F3+…+bnFn

(5)

F1、F2、…、Fn是主成分分析得到的主成分，y是因变量，a、b1、b2、…、bn为方程的参数。通过计算回归方程中各主成分的载荷，以及原始变量与求取的综合变量之间的线性关系，获得7个初始变量在回归方程中各自的回归参数。

选用方差膨胀因子对构建的多元回归方程进行共线性分析，结果表明，除4月上旬NDVI的和x4，5月上旬NDVI的和x6及5月下旬NDVI的和x7这三者的方差膨胀因子小于10以外，其他几个指标的方差膨胀因子均大于10，所以拟合的回归方程存在有较强的多重共线性问题。对7个自变量进行主成分分析，提取得到2个主成分及其相应线性关系式见式(6)、(7)。

F1=0.110x1+0.261x2+0.632x3+0.816x4+0.756x5+0.885x6+0.644x7

(6)

F2=0.951x1+0.880x2+0.748x3-0.131x4-0.336x5-0.152x6-0.485x7

(7)

利用得到的主成分F1、F2作为自变量，冬小麦产量为因变量进行回归分析，得到方程见式(8)。

y=970 904.38+3.769F1+7.484F2

(8)

该方程的R2为0.884，表明该方程有较好的拟合精度。将上式中的2个主成分F1和F2分别转化为7个初始变量，得到用自变量表示的回归方程见式(9)。

y=7.376x1+4.491x2+8.903x3+5.979x4+5.931x5+8.707x6+2.686x7+970 904.38

(9)

2.4 冬小麦估产模型的验证

为了检验估产模型是否具有可行性以及各自的准确性，利用2014年开封市的冬小麦产量数据和其NDVI数据对建立的3个模型分别进行验证，得到结果如表2所示。从表2可以看出，3个遥感估产模型的误差均是负值，相对预测误差在6.88%～10.55%。模型1是基于小麦产量和与产量相关性最大的3月份下旬的NDVI数据构建的一元回归模型，其相对误差为-10.55%，对产量的预测相对较准确。分析开封市作物生长的物候期得知，3月下旬冬小麦正处在拔节期，这一时期对冬小麦穗粒数和穗数的多少起着关键作用，对最终产量的形成相当重要。

与模型1相比，综合利用冬小麦生长各个时期的NDVI值建立的多元回归模型具有更高的R2，预测精度更高，其产量的拟合精度为92.92%，可能是由于该方法系统考虑了不同时期的NDVI状态，涵盖更加丰富、全面的信息，有较好的灵敏度，综合、全面地反映冬小麦在各个生长时期的情况，因此，模型2的产量估测度较模型1有了很大提高。分析模型2中各个变量的回归系数可以发现，拔节期(3月下旬)的NDVI数据的回归系数最大，即其对产量的影响最大，且和产量之间呈正相关关系，与相关分析的结果吻合，开封市拔节期是冬小麦产量预估的最佳生育时期，其次为返青期和孕穗期。

构建的冬小麦生长关键期内NDVI值的主成分回归产量预测模型主要是对所建立的多元回归模型的优化，将原来的7个变量进行线性组合，转换为2个相互独立且代表性好的综合主成分因子。主成分分析得到的2个综合主成分因子的载荷矩阵如表3所示，比较分析所提取的2个主成分，其中，第一主成分能更好地反映冬小麦抽穗期至成熟期的信息，即包含4月上下旬和5月上下旬的NDVI数据，而第二主成分主要是2月下旬和3月上下旬的NDVI数据，处于冬小麦的返青期至拔节期。综上，通过主成分分析对变量进行科学整合，得到了代表冬小麦生育期的2个关键时期的NDVI，分别是代表冬小麦从抽穗期至成熟期内生长状态的第一主成分，及表示返青期至拔节期的第二主成分，2个主成分均有很高的代表性。主成分回归模型的估产水平最高，精度达到93.12%。对比分析该模型各个系数可知，3月下旬的NDVI的系数最大，从而再次验证3月下旬(拔节期)对产量的形成影响最大。

表2 冬小麦估产模型验证结果

表3 主成分因子载荷矩阵

3 结论与讨论

作物每个物候期的生长状况均与其最终产量的形成密切有关[10]。选择与冬小麦产量形成相关的重要物候期NDVI作为本研究的分析指标，用于最终产量的预测，建立最佳观测期NDVI与最终产量的回归模型及各生育时期NDVI和产量的多元回归模型，并在多元模型基础上，有效提取涵盖原变量的综合变量，经由回归统计分析，得到最佳遥感估产模型。

1)通过2005—2013年冬小麦关键生长时期(2—5月)的NDVI与产量的相关性分析发现，3月下旬(拔节期)为开封市冬小麦估产的最佳时期，且多元回归模型和主成分回归模型中拔节期(3月下旬)NDVI值的回归系数最大，证实了拔节期是研究区域内冬小麦最佳估产期。

2)比较基于最佳估产期和生长关键期NDVI分别建立的一元回归模型和多元回归模型，后者的拟合优度R2是0.976，估产精度为92.92%，均高于前者。据此可知，综合考虑冬小麦不同生育时期NDVI状态，对于提高估产模型的稳定性和预测精度均有帮助。

3)采用主成分回归法优化多元回归模型可以有效提高模型的估产精度，估测精度达到93.12%，在3个模型中估产效果最好。说明选用冬小麦所有生长关键时期的NDVI信息并科学地提取综合指标进行回归建模，能达到更高的预测精度。

本研究所构建的3种估产模型，均可以利用NDVI进行冬小麦估产，估产误差均控制在10.55%内，证明了利用NDVI进行冬小麦估产的可行性，但在研究过程中还存在有很多不足，有待进一步的论证和完善。

1)研究中冬小麦种植区域的提取受TM影像数据源的影响，每年所依据的影像在时相上并不完全统一，并且，研究区部分地区种植有蒜苗等其他蔬菜，易于混淆，造成遥感解译的误差，使得每年提取的冬小麦种植面积与实际面积的统计值存有差异。因此，可以选择更高时空分辨率的遥感数据或采用其他辅助数据，来提取获得冬小麦种植区域。

2)本研究选用NDVI数据进行冬小麦产量估测，主要是因为其应用较为广泛。而增强化植被指数EVI将抗大气植被指数和土壤调节植被指数进行耦合从而实现对植被信息优化[11]，并降低了植被冠层背景不稳定和大气衰减的影响[12]，可在今后的冬小麦估产研究中进行应用。

3)冬小麦产量受多个因素的影响，如水热条件、土壤质量等，可以以植被指数为基础引入气候等相关因素进行综合分析，这有待进一步的研究。

[1] 金涛.中国粮食生产时空变化及其耕地利用效应[J].自然资源学报,2014,29(6):911-919.

[2] 杜天昊.基于MODIS-NDVI的河南省冬小麦估产模型研究[D].重庆：重庆大学,2013.

[3] 赵庚星,余松烈.冬小麦遥感估产研究进展[J].山东农业大学学报(自然科学版),2001,32(1):107-111.

[4] 王东伟.遥感数据与作物生长模型同化方法及其应用研究[D].北京：北京师范大学,2008.

[5] 闫岩,柳钦火,刘强,等.基于遥感数据与作物生长模型同化的冬小麦长势监测与估产方法研究[J].遥感学报,2006,10(5):804-811.

[6] 任建强,陈仲新,唐华俊,等.基于植物净初级生产力模型的区域冬小麦估产研究[J].农业工程学报,2006,22(5):111-117.

[7] 黎锐,李存军,徐新刚,等.基于支持向量回归(SVR)和多时相遥感数据的冬小麦估产[J].农业工程学报,2009,25(7):114-117.

[8] 朱再春,陈联裙,张锦水,等.MODISEVI时间序列数据和光谱角聚类的冬小麦遥感估产分区方法研究[J].光谱学与光谱分析,2012,32(7):1899-1903.

[9] 任建强,陈仲新,唐华俊.基于MODIS-NDVI的区域冬小麦遥感估产——以山东省济宁市为例[J].应用生态学报,2006,17(12):2371-2375.

[10] 任建强,陈仲新,唐华俊,等.长时间序列NOAA-NDVI数据在冬小麦区域估产中的应用[J].遥感技术与应用,2007,22(3):326-332.

[11] 王长耀,林文鹏.基于MODISEVI的冬小麦产量遥感预测研究[J].农业工程学报,2005,21(10):90-94.

[12] Huete A,Didan K,Miura T,etal.Overview of the radiometric and biophysical performance of the MODIS vegetation indices[J].Remote Sensing of Environment,2002,83(1/2):195-213.

Regional Scale Yield Estimation of Winter Wheat by Remote Sensing Using MODIS-NDVI Data in Different Growth Stages

LI Bing1,2,LI Shimin1,ZHOU Lei3,HUANG Canhui3,WANG Laigang1*

(1.Institute of Agricultural Economics and Information,Henan Academy of Agricultural Sciences,Zhengzhou 450002,China; 2.Zhengzhou Shuqing Medical College,Zhengzhou 450000,China;3.Remote Sensing Monitoring Center of Agriculture of Henan Province,Zhengzhou 450002,China)

This study took Kaifeng city of Henan province as study area.The correlation between the sum of 16 d synthetic 250 m normalized difference vegetation index(NDVI) in winter wheat growth stages from 2005 to 2013 and winter wheat yield was analyzed to obtain the best phase of yield estimation by remote sensing to construct single variable regression model.Meanwhile,the multiple regression model based on NDVI of winter wheat in different growth stages was constructed.Then,the principal component analysis method was used to improve the parameters of multiple regression model.Finally,the yield data of Kaifeng city in 2014 was used to verify the yield estimation model,aiming to build the yield estimation model with higher precision,so as to better guide the wheat production.The results showed that the estimation errors of three models were controlled within 10.55%,according to the three models,the optimal period of yield estimation of winter wheat was in late March,ie,jointing stage,the principal component regression estimation model had the highest fitting precision of 93.12%,and had some practical value.

yield estimation by remote sensing； MODIS； NDVI； regression model； winter wheat

2017-01-15

河南省科技攻关重点项目(172102110090)；河南省农业科学院自主创新项目(2017ZC60)

李 冰(1980-)，女，河南开封人，讲师，硕士，主要从事遥感图像自动分类研究。E-mail:42383949@qq.com

*通讯作者：王来刚(1979-)，男，河南辉县人，副研究员，博士，主要从事农业遥感应用研究。E-mail:wlaigang@sina.com

S512.1；S126

A

1004-3268(2017)06-0150-06