安徽省(29°41'N—34°38'N, 114°54'E—119°37'E)地处我国南北交界处, 地形种类齐全, 皖西、 皖南以丘陵山地为主, 皖北为平原地区, 总面积约14万km²。 该区域年均温为15 ℃左右, 年均降水量介于700~1 600 mm之间, 全年光热条件良好, 适宜多种农作物的生长, 为粮食产业提供了优越的发展条件。

遥感数据来自2012年—2018年3月—5月MOD15A2H8d合成的叶面积指数产品, 每年12景, 空间分辨率为500 m。 该产品已经过辐射定标、 大气校正以及几何校正等处理, 在MRT(MODIS Reprojection Tool)中完成区域拼接与投影转换。

逐日气象数据来自中国科学院青藏高原研究所的《中国区域地面气象要素驱动数据集》^{[8, 9]}, 时间分辨率为3 h, 空间分辨率为0.1°×0.1°, 包含7个气象要素: 近地面气温、 降水、 气压、 比湿度、 10 m风速、 下行短波和长波辐射。 经过格式转换与要素计算得到模型所需的气象要素。 饱和水汽压差(vapor pressure deficit, VPD)由近地表2 m气温数据与露点温度计算得到, 数据来自欧洲中期天气预报中心(European Centre for Medium-Range Weather Forecasts, ECMWF)发布的ERA5数据集。

研究选用安徽省冬小麦生长关键期多幅遥感影像, 以县为单位, 使用空间分辨率16 m左右的GF1和Landsat8影像, 每县至少2幅影像, 关键时相主要为农作物未播种前期上年10月—12月上旬影像(北部主要为10月中旬影像, 中南部主要为11月—12月上旬影像), 生长中期本年3月中下旬—4月上旬影像, 基于物候数据及地面调查数据选取训练样本, 主要利用最大似然法(少部分县使用神经网络方法)进行监督分类; 采用2 m空间分辨率影像提取耕地图斑, 与分类结果求交集, 获得安徽省冬小麦空间分布。 训练和精度验证样本为全省均匀分布的1 198个200 m×200 m样方压盖的自然地块。 主要使用2 m分辨率GF1影像, 补充使用0.8 m GF2和BJ2影像重采样为2 m影像, 人工勾画地块及目视解译判读地块类型。 自然地块的作物属性由人工野外调查获得。 全部样方中的90%用于训练, 10%用于精度验证, 采用误差矩阵评价分类精度。 经验证各县的分类精度可达85%。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 研究区概况Fig.1 Overview of the study area

实割实测产量数据来自国家统计局, 为地块尺度数据, 处理成村级尺度的冬小麦单位面积产量。 2012年—2018年间, 全省共包含1 265个可用的村级站点/年组合, 其中皖北650个, 皖中506个, 皖南109个。

SCYM是Lobell团队提出的一种新的多尺度农作物估产框架, 旨在利用作物生理知识去解释卫星观测结果, 可快速应用于不同田块、 不同种类作物的估产, 估产精度高并且可扩展性强。 基本流程如下: (1)针对一系列合理的管理措施和环境条件进行作物生长模拟。 (2)训练模拟叶面积指数(leaf area index, LAI)、 季节性天气变量与模拟产量的随机森林模型。 (3)基于遥感观测LAI和随机森林模型进行产量估算。

WOFOST模型由荷兰Wageningen农业大学和世界粮食研究中心(CWFS)共同开发研制, 以作物品种、 土壤属性、 气候资料及田间管理措施为驱动, 逐日模拟作物出苗到成熟的生理过程, 目前在诸多国家和地区得到了广泛应用^{[10, 11, 12]}。 LAI是模型中最重要的变量之一, 其影响着作物的干物质积累, 进而决定了最终产量, 并将贮藏器官干重(WSO)的最大值作为冬小麦模拟产量指标。

随机森林由多棵决策树组成, 本质是将若干个弱分类器集成为强分类器的集成算法。 利用bootstrap有放回地抽取与样本集相同大小的训练集, 构建分类回归树, 针对每棵树的分类结果, 采取投票(取众数)或算术平均的方法得到分类或回归结果, 误差则基于未被抽取的样本组成的数据集来计算。 本研究中将季节性天气变量(播种—越冬期积温、 全生育期降水量、 5月太阳辐射、 12月—2月平均气温)及不同观测日期下的模拟LAI作为随机森林的特征变量集合, 模拟产量作为输出结果训练模型。

1.9.1 WOFOST模型参数优化

WOFOST模型中作物品种、 土壤属性、 田间管理等参数同天气变化一起影响着最终产量的形成, 模型模拟结果对不同参数变化的响应程度不同。 搜集了WOFOST模型模拟冬小麦生育动态相关的大量文献, 在这些文献中完成了参数敏感性分析以及本地化率定工作, 并提供了全部参数的取值范围(表1)。

表1

Table 1

表1(Table 1)

表1 WOFOST模型中冬小麦主要参数标定范围 Table 1 Range of calibrated values of maincrop parameters of WOFOST model 参数名称 定义 单位 取值范围 TSUM1 出苗至开花有效积温 ℃·d 600.0~1 500.0 TSUM2 开花至成熟有效积温 ℃·d 500.0~1 200.0 IDSOW 播种日期 d 273~313 TDWI 初始生物量 kg·hm-2 100~210 SPAN 叶片衰老系数 d 21~41 LAIEM 出苗时叶面积指数 hm-2·hm-2 0.136 5 CVL 叶生长同化物转化效率 kg·kg-1 0.485~0.685 CVO 贮存器官生长同化物转化效率 kg·kg-1 0.729~0.81 CVR 根同化物转化效率 kg·kg-1 0.694 CVS 茎同化物转化效率 kg·kg-1 0.662~0.762 EFFTB 对应生育期单片光能利用效率 kg·hm-2·h-1·J-1·m-2·s-1 0.35~0.5 AMAX01(DVS=0.00) AMAX02(DVS=1.00) AMAX03(DVS=1.30) AMAX04(DVS=2.00) 叶片最大CO2同化速率 kg·hm-2 35.83 38.83 30~35.83 4.48 SLA01(DVS=0.00) SLA02(DVS=0.50) SLA03(DVS=2.00) 比叶面积 hm2·kg-1 0.001 25~0.002 12 0.001 35~0.002 75 0.0~0.001 5 注: 在WOFOST模型中物候用无量纲状态变量DVS来描述, 对于大多数一年生作物, DVS在出苗期设为0, 开花期设为1, 成熟期设为2 Note: In WOFOST phenology is described by the dimensionless state variable DVS. For most annual crops, DVS is set to 0 at seedling emergence, 1 at flowering and 2 at maturity

表1 WOFOST模型中冬小麦主要参数标定范围 Table 1 Range of calibrated values of maincrop parameters of WOFOST model

鉴于田间实验数据难以获取这一重要原因, 基于大量前人研究成果^[10], 挑选出高敏感性的生长参数: 从出苗到开花的积温TSUM1、 叶片衰老系数SPAN、 从出苗到开花的比叶面积SLATB2以及播种日期IDSOW, 综合考虑安徽省气象条件和作物生长物候期, 分别在其适宜范围内取等间隔的10个值, 对四个参数的取值排列组合, 以各站点气象数据为驱动, 分别在皖北、 皖中、 皖南区域进行模拟。

1.9.2 产量预测模型建立

基于前人研究成果, 结合安徽省冬小麦的各生育阶段, 设计出以下指标作为待选变量(表2)。

表2

Table 2

表2(Table 2)

表2 待选择的季节性天气变量 Table 2 Seasonal meteorological weather indicators to be selected 天气变量 待选指标 平均气温/℃ 12月平均气温 1月平均气温 2月平均气温 12月—2月平均气温 3月—5月平均气温 平均最低气温/℃ 12月平均最低气温 1月平均最低气温 2月平均最低气温 12月—2月平均最低气温 降水量/mm 越冬前降水量 返青期降水量 拔节期降水量 抽穗扬花期降水量 灌浆期降水量 全生育期总降水量 太阳辐射/(kJ·m-2) 灌浆期太阳辐射 3月太阳辐射 4月太阳辐射 5月太阳辐射 3月—5月平均太阳辐射 饱和水汽压差/kPa 3月平均饱和水汽压差 4月平均饱和水汽压差 5月平均饱和水汽压差 3月—5月平均饱和水汽压差 有效活动积温/℃ 播种—越冬期积温 返青—成熟期积温 全生育期积温 越冬期负积温

表2 待选择的季节性天气变量 Table 2 Seasonal meteorological weather indicators to be selected

利用安徽省每年10月至次年6月的再分析资料, 结合冬小麦生育期计算出以上指标。 利用待选指标对相应站点/年的模拟产量逐步回归, 最终选取播种—越冬期积温、 全生育期降水量、 5月太阳辐射及12月—2月平均气温作为构建产量预测模型的天气变量。

将返青至成熟期(3月—5月)按8 d为单位分段, 提取每个时段的最大模拟LAI, 选择以8 d为单位是便于后续将模拟LAI与LAI遥感观测值一一对应。 将上文选取的季节性天气变量与两个不同时段的LAI值排列组合, 作为随机森林的输入特征变量, 并将不同日期组合训练成的随机森林模型储存起来, 供后续使用。

1.9.3 产量估算

根据MOD15A2H产品的质量控制文件及冬小麦空间分布数据筛选出总体质量好且属于冬小麦种植区的像元。 产量验证数据为村级尺度, 而MODIS产品分辨率为500 m×500 m, 需要对LAI数据进行均值滤波处理。 付立哲等研究表明, 由于农作物种植范围相对较小, MODISLAI数据中不可避免地会出现大量混合像元, LAI被低估的现象普遍存在^[13]。 利用2012年—2018年3月—5月均值滤波后的LAI数据与相应时段的模拟LAI进行线性拟合, 以修正遥感观测值。

对不同时相LAI组合下模型的解释能力按照从高到低排序, 提取出每个站点/年所对应像元的最佳可用观测日期组合。 将每个站点/年的季节性天气变量及两个时段的LAI修正值代入到训练好的随机森林模型中预测产量。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 WOFOST模拟产量与站点实测产量的比较Fig.2 Comparison of WOFOST simulated yield with site-measured yield

将全部参数组合下各站点/年的模拟产量与地面实测单产线性拟合, 剔除异常值点, 确定拟合优度最高的参数组合。 最终WOFOST模型估算产量与站点实测产量的R²达到0.421, 皖北地区为0.452, 皖中地区为0.376, 皖南地区为0.668。 随机抽取50次最优参数组合下的模拟结果, 绘制冬小麦生育期LAI及WSO时序曲线[图3(a)和(b)]。 可以看出LAI最大值普遍在3~7之间, 冬小麦单产在5 000~8 000 kg·ha^-1波动。 在不同种植场景下, 安徽省冬小麦的生长情况呈现出广泛的变异性。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 50次随机模拟的LAI(a)和WSO(b)时间序列曲线Fig.3 LAI (a) and WSO (b) daily outputs from 50 random simulations of the WOFOST

图4表示不同时段LAI组合作为特征变量时模型的解释能力。 可以看出, 基于LAI模拟值和季节性天气变量的产量预测模型对模拟结果变异性的解释能力取决于LAI的具体观测日期, 但总体上决定系数较高, 均在0.68以上。 相较于前期(3月初), 在冬小麦发育后期(4月中旬以后)进行LAI观测时模型的解释能力更强。 在冬小麦播种来年的5月9日至16日、 5月17日至5月24日之间进行LAI观测, 可以获得该区域最准确的产量估计值。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 基于LAI模拟值和季节性天气变量的随机森林模型的决定系数(R2)Fig.4 The coefficient of determination (R2) for the Random Forests models based on simulated LAI and monthly weather

2.2.1 站点尺度估产精度分析

将全部站点/年组合视为一个整体时, SCYM的产量估算值与实测值总体相关性为0.758(R²为0.575)(图5), RMSE为790.92 kg·ha^-1。 其中皖北地区冬小麦产量预测效果最佳, 相关系数可达0.604(R²为0.365), RMSE为807.77 kg·ha^-1, 皖南地区预测效果不佳, 相关系数仅为0.210(R²为0.044), RMSE为786.56 kg·ha^-1。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 SCYM估算产量与站点实测产量的比较Fig.5 Comparison of SCYM estimated yield with site-measured yield

从SCYM估产的空间分布差异来看(图6), 估产精度的空间异质性显著, 估产精度较高的站点主要分布在淮北平原, 相对误差普遍低于1%, 而高误差区域集中于皖南丘陵地带, 最大相对误差为60.6%。 可能的原因是皖北种植规模较大, 产量验证点分布较为密集, 而皖南多为山地丘陵, 站点分布稀疏, 为SCYM框架解释低产信息带来了一定的困难。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 SCYM估算产量与站点实测产量的相对误差空间分布Fig.6 Spatial distribution of relative error between the SCYM estimated yield and the site-measured yield

2.2.2 安徽省冬小麦估产精度分析

以季节性天气变量及每个像元最佳观测日期组合下的LAI为驱动, 对2012年—2018年安徽省冬小麦的产量进行预测, 多年平均产量分布情况如图7。 图中白色为非耕地区域, 小麦单产在4 531.91~6 997.58 kg·ha^-1。 总体而言, 小麦单产由北向南逐渐减少, 高值区出现在皖北的淮北平原, 低值区则在皖中、 皖南均有分布。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 2012年—2018年SCYM估算产量平均值Fig.7 Average SCYM estimated yield in 2012—2018

将SCYM估产结果与站点实割实测数据逐年比较(图8), 2014年后估算产量持续下降, 分析七年的天气变量发现, 2014年后5月太阳辐射总量下降11.92%, 全生育期降水量增加22.94%, 多雨寡照的气候条件对冬小麦灌浆十分不利, 从而导致冬小麦减产。 自2015年开始, 播种—越冬期积温、 12月—2月平均气温不断升高, 研究发现, 气候变暖导致小麦发育前期旺长, 若气温回升而不及时调控, 后期易早衰、 倒伏, 从而产量大幅下降。 但总体来看, 2012年—2018年实测单产平均值为6 058.00 kg·ha^-1, SCYM估算单产平均值为5 984.95 kg·ha^-1, 两条时序曲线的相关性达到0.822, RMSE为189.96 kg·ha^-1, 每年估产的相对误差均不超过6%, 说明该方法具有一定的可行性, 可为安徽冬小麦的产量估算提供技术支持。

图8

Fig.8

	Figure Option ViewDownloadNew Window
	图8 安徽省冬小麦单产变化趋势Fig.8 Trend of winter wheat yield in Anhui Province