基于长时间序列遥感数据反演NPP的耕地质量评价

马佳妮 张 超 吕雅慧 高璐璐 郧文聚 朱德海

(1.中国农业大学信息与电气工程学院， 北京 100083； 2.中国农业大学土地科学与技术学院， 北京 100083； 3.自然资源部农用地质量与监控重点实验室， 北京 100035)

0 引言

耕地质量管护对耕地数量、质量、生态三位一体保护以及耕地占补平衡具有重要意义，而掌握耕地质量空间分布是耕地质量管护的重要前提。目前，耕地质量监测与评价技术方法仍以区划—布点—调查—评价为主[1-2]，但在样点布设科学性、数据与结果时效性、空间连续性等方面尚存在局限。遥感数据具有覆盖面积大、获取数据快、时空连续的特点，可有效弥补传统方式的不足，为耕地质量监测与评价提供有效手段。

基于遥感的耕地质量监测与评价研究主要包括耕地质量指标遥感监测和耕地质量遥感反演。在耕地质量指标监测中，黄健熙等[3]利用多年蒸散发产品结合年度有效灌溉数据评价区域灌溉能力；吕雅慧等[4]采用0.5 m的GEO Eye-1 影像，构建归一化差值植被指数(Normalized difference vegetation index, NDVI)与二维熵的决策树，实现了耕地质量重要内容农田林网的自动识别；杨建锋等[5]通过建立多光谱遥感影像波段与有机质含量、地形坡度、表层土壤质地的反演模型，获得了耕地质量自然等。以上通过遥感数据信息提取与波段信息反演解决了部分指标获取困难的问题。此外，还有通过遥感数据直接监测和评价耕地质量相关研究，林晨等[6]建立了基于MODIS数据的耕地质量自然等反演模型。以上耕地质量监测与评价多以一个时间断点，忽略了时间断点数据不稳定性与评价结果的时效性。

已有学者将多期遥感数据应用于土地监测与评价[7-11]。研究表明，净初级生产力(NPP)较NDVI更能全面反映作物长势，已广泛应用于农田生产力评价[12-13]。以上研究采用时序数据弥补了时间断点数据的不稳定性，但利用植被指数表征全域各类土地质量状况，忽略了不同作物监测指标的差异，造成结果的不可比性。耕地质量各构成要素的特点和相互间的影响，决定了耕地质量的外在表现[14]。耕地质量包含自然质量、利用质量与经济质量等[15]。而作物长势是耕地质量利用的外在表现，产量与产能是反映耕地质量的重要指标。依据吉林省气象站点物候统计，吉林省玉米与水稻生长期为5—9月。农田NPP在作物生长期变化明显[16]，耕地质量影响因素如田间施肥、灌溉等人为因素和土壤理化性质在作物生长期变化明显。多年作物生长期的NPP体现了多年影响耕地质量的各因素对作物的贡献程度，间接评价了耕地质量。因此，本文选择吉林省为研究区，采用2000—2010年5—9月500 m的MODIS MOD09A1遥感数据，结合气象站点数据，构建省域耕地质量利用等评价方法，考虑不同作物植被指数差异，分别计算反映不同作物长势的指标，经距平分析法消除量纲，利用反映多年作物长势的NPP评价耕地质量。

1 研究区和数据

1.1 研究区概况

选择吉林省为研究区，地理位置为东经121°38′～131°19′，北纬40°50′～46°19′，如图1所示。吉林省属于温带大陆性季风气候，四季分明，雨热同季。多年平均气温2～6℃，降水量400～900 mm。以中部大黑山为界分为东部山区与西部平原。全省耕地共553.5万hm2，占吉林省面积28.98%。70%耕地分布在中、西部平原地区。全省耕地主要为旱地和水田。耕作制度为一年一熟制，种植结构稳定，主要农作物有玉米和水稻等，粮食单产地区间差异大，粮食平均单产在全国处于中等偏下水平。耕地质量利用等别分布在8～13等别。耕地后备资源相对充足，有31万hm2耕地后备资源。

图1 研究区位示意图Fig.1 Geographical distribution map of study area

1.2 数据源

1.2.1遥感数据

采用2000—2010年玉米与水稻生长期5—9月的8 d合成500 m分辨率MODIS地表反射率产品MOD09A1(www.gscloud.cn)，共215景。经投影转换、镶嵌、重采样预处理后，形成2000—2010年长时间序列遥感数据集。

1.2.2气象数据

气象数据来源于中国气象科学数据共享网，包括吉林省内与周边共51个国家气象站点，包括2000—2010年间逐日最大温度、最小温度、总辐射量、日照时数等数据。

1.2.3光合有效辐射数据等

光合有效辐射数据(Photosynthetically active radiation, PAR)是利用晴空指数、太阳高度角与日照时数建立的估算模型而得[17]。其他数据包括2013年吉林省农用地等别成果中的国家利用等别数据以及吉林省行政边界数据。

2 研究方法

NPP是植物光合作用所固定的有机物总量扣除植物自养呼吸后的剩余部分。单年NPP值易受气候等因素影响，用于评价长期稳定的耕地质量误差较大。计算2000—2010年的NPP均值评价耕地质量可减少气候等因素引起的不确定性。基于影像计算每8 d的NPP值，经累加求得年NPP值，通过计算多年NPP均值，评价耕地质量。由于NPP受作物类型的影响[18-20]，需对研究区主要作物水稻(C3)与玉米(C4)进行分类，利用不同作物生长期对指数响应程度不同，构建基于时序指数的分类特征，并采用决策树主要区分水稻与玉米，剩余少量大豆与水稻均为C3作物。不同作物类型(主要分C3和C4)分别利用植被光合作用模型(Vegetation photosynthesis model，VPM)计算NPP，并经累加计算2000—2010年NPP年累积量。后经距平分析法去除作物类型间差异，结果用于评价耕地质量利用等别。相应的技术路线如图2所示。

图2 技术路线图Fig.2 Flow chart of method

2.1 基于时序遥感数据的水稻与玉米识别

由于水稻、玉米的光合作用过程、光能利用率、积温等差异，导致两种作物类型计算NPP不同，所以首先利用多时相进行农作物分类。通过构建增强植被指数(Enhanced vegetation index, EVI)、地表水指数(Land surface water index, LSWI)与NDVI时序数据集，选取水稻移栽期、移栽期后40 d及收获后为关键时期，采用决策树分类区分水稻和玉米。

2.1.1遥感指数计算

利用MODIS数据计算NDVI、EVI和LSWI，计算公式为

(1)

(2)

(3)

式中RB——蓝波段波长(459～479 nm)

RR——红波段波长(620～670 nm)

RNIR——近红外波段波长(841～876 nm)

RSWIR——短波红外波段波长(1 628～1 652 nm)

尽管8 d合成MODIS数据最大限度削弱低云、云阴影、气溶胶等对图像质量的影响，然而6、7月存在长时间云覆盖情况，故利用S-G滤波消除云的影响[21]。

2.1.2基于时序指数的水稻玉米分类

依据水稻与玉米物候期时序变化差异(图3)，采用关键生育期LSWI和EVI进行作物分类。选择水稻移栽期、移栽期后40 d及收获后作为分类关键期。移栽期水稻田因灌水存有2～15 cm积水[22]， LSWI时序曲线出现突增，EVI出现突降。采用EVI-LSWI≤0.05作为特征参量。移栽40 d后，水稻进入返青期，分蘖数逐渐增加，EVI迅速增大，该时段利用EVI>0.5EVImax作为分类特征。在收获后期，利用NDVI<0.6，采用决策树实现作物分类。

图3 水稻与玉米时序图Fig.3 Rice and corn timing chart

2.2 NPP计算

VPM是基于遥感数据进行NPP估算的光能利用率模型[23]。首先计算

GPP=εgFPARchlPAR

(4)

其中

FPARchl=aEVI

(5)

εg=ε0WscalarTscalarPscalar

(6)

(7)

式中εg——光能利用率，g/MJ，受温度、叶龄与水分的影响

FPARchl——植被光合的部分吸收光合有效辐射比例

a——系数，取1

ε0——植被最大光能利用率，g/MJ，C4植物最大光能利用率高于C3植物[24]，玉米与水稻的ε0分别为3.65 g/MJ和1.15 g/MJ[25-26]

Tscalar——温度系数

Wscalar——水分系数

Pscalar——叶龄系数

tlswimax——叶面完全伸展时期

Tmin、Tmax、Topt——作物的最低温、最高温、最适温[27]，玉米、水稻的最低温分别为0、-3℃，最高温为45、42℃，最适温为23、16℃

2.3 耕地质量利用等别评价

计算水稻与玉米年NPP累积量后，通过距平分析法消除作物类型间的量纲，即

(8)

式中NPP——11年NPP均值

D——像元距平百分率

最终根据D值参照农用地等别的分级标准得到耕地质量的评价结果。

3 实验结果与分析

3.1 多时相遥感作物分类

由于水稻和玉米是吉林省主要作物，多年种植结构稳定，故采用2006年分类结果作为11年后续研究的基础数据。分类样本来源于实地采样结合现有耕地图斑，分类结果如图4所示，总体分类精度为80.56%，Kappa系数为0.7。吉林省玉米面积为70.75万hm2，水稻面积32.82万hm2。玉米种植区集中在北部白城市、中部松原市。水稻种植区位于中东部通榆市、德惠市与中西部四平市，处于松花江与东辽河两侧。

图4 2006年吉林省水稻玉米分布图Fig.4 Rice and maize distribution map of Jilin Province in 2006

3.2 VPM模型计算NPP

利用VPM模型计算8 d NPP。后经累加得到年NPP值。并生成2000—2010年吉林省NPP空间分布图，如图5所示，全域NPP累积量在[0，507.249]之间。

3.3 耕地质量评价

利用距平分析法将归一化后的D值作为表征耕地质量。并与2013年农用地分等中吉林省耕地质量利用等对比验证。为与吉林省国家利用等别等级一致，按照吉林省国家利用等别各等别面积占比，将耕地质量反演结果分为9个等别，分别对应吉林省国家利用等的6～14等别，6为最高等别，按序等别降低。从图6b可见，全域耕地质量呈阶梯状分布。北部白城市与松原市耕地质量处于11～14等别。中部长春市、四平市耕地质量处于9～10等别，少量分布8、12等别。南部地区在9～14等别都有分布，主要集中在10～13等别。北部与南部质量低，中部质量高。由图6可得，反演结果与国家利用等别在整体空间分布较为一致，北部反演结果集中在9～12等别，说明南部整体反演等级偏高，主要是由于南部耕地零碎分散并存在于耕地与林地混合的像元导致反演结果偏高。对比发现，两种结果在吉林省北部均处于9～12等别，但2013年吉林省农用地分等国家利用等别的结果位于北部的乾安县整体处13～14等别，与周边县等级差异较大。图6a中北部的评价结果空间连续，符合耕地质量实际分布。图6b中县域存在异常值的原因是现有耕地质量评价是以县域为评价单元，各县在指标获取与计算的差异造成县域结果整体偏高或偏低。本文方法可检验县域评价结果异常。

图5 2000—2010年吉林省累积NPP空间分布Fig.5 Cumulative NPP spatial distribution during 2000—2010 in Jilin Province

图6 本文评价结果与2013年吉林省农用地分等国家利用等别对比Fig.6 Comparison of evaluation results with utilization of agricultural land in Jilin Province in 2013

采用反演结果等级与农用地分等成果等别之差，评价构建方法的可靠性。空间分布如图7所示，级差大于2的区域主要集中在南部的敦化市、汪清县、柳河县、桦甸市，中部的公主岭市及北部镇赉县。说明本文的评价结果高于农用地分等结果，主要因为耕地零散破碎，林地NPP高于耕地NPP，耕地与林地混合像元造成反演结果较高。级差小于-2的区域在中部的蛟河市与北部的洮南市。

如图8所示，级差在-1～1级有73.72%耕地，级差-2～2级的耕地面积有17.01%，级差在3级及以上的耕地占9.27%，说明本文提出方法的评价结果与农用地分等中的分级结果具有较好的吻合度。

图7 本文评价结果与2013年吉林省农用地分等国家利用等别级差分布Fig.7 Inversion results and spatial distribution of difference in utilization of agricultural land in Jilin Province in 2013

图8 级差概率分布Fig.8 Level difference probability distribution

4 结论

(1)构建了利用长时间影像数据集耕地质量评价方法，通过监测多年反映作物长势的NPP表征耕地质量。在耕地集中连片的平原地区评价结果与现有数据结果一致，同时表明，500 m的遥感数据可用于平原耕地集中连片区域的耕地质量监测与评价。

(2)利用多时序植被指数作物分类，选择水稻

移栽期、移栽期后40 d及收获后作为关键生育期，移栽期未被稻苗遮挡的田间积水，在影像上可作为区别水稻田与旱地的特征。移栽期后40 d，水稻处于急速生长的分蘖期。收获后田间使用LSWI、EVI、NDVI指数特征参与分类，最终的分类精度达80.56%。

(3)采用2000—2010年MOD09A1影像数据，结合气象站点数据，经VPM模型计算水稻和玉米NPP值，并得到多年NPP累积量。使用距平分析法消除作物类型差异，以此表征耕地质量，并对其结果的合理性进行验证，结果表明，该方法适合评价吉林省中部与北部等平原地区，在耕地质量评价方面具有可行性。

(4)基于长时间序列遥感数据反演NPP的耕地质量评价，利用多年作物长势间接反映耕地质量，避免单年评价的偶然性，提高了准确性。同时利用覆盖全域的MODIS影像，有效解决了省域耕地质量不可比、评价由点及面与评价周期长等问题。在获取耕地质量空间分布的同时，可以获得耕地质量时空演变过程。