4种水稻气象产量分离方法的对比研究

葛道阔，曹宏鑫，张文宇,2，张伟欣，金之庆，李秉柏，唐普传,2

(1.江苏省农业科学院 农业信息研究所，江苏 南京 210014；2.江苏大学 农业装备工程学院，江苏 镇江 212013)

一个多世纪以来，全球气候变化已经是一个不容置疑的事实和不可阻挡的自然趋势[1-2]。自改革开放起，我国作物生产系统的生产力，虽然随着持续的科技进步得到不断提升，但全球气候变化始终威胁其稳定性，进而影响到国家的粮食安全[3]。以气候本身及其变率变化为特征，造成的农业气象灾害，对包括江苏在内的长江中下游省份的水稻生产产生了严重的影响[4]。采用科学、可靠的方法，评估气象灾害的影响范围与程度，因地制宜地制定和采取粮食生产的应对策略，对保障我国优势产区的粮食安全具有十分重要的意义[4-7]。为了深度解析气象灾害对作物产量的影响，通常采用以下3种方法：一是作物模拟试验方法，二是生长模型方法，三是生物与数理统计方法[6-8]。后者是把作物产量分解为3个部分，即作物趋势产量、作物气象产量和随机误差。作物趋势产量，是综合反映一个历史时期不同技术水平以及农业政策、投入、品种等的长周期分量，作物气象产量则是由农业气象灾害主导的短周期分量。一般通过数学模型或统计方法将作物气象产量波动从实际产量中分离出来，然后分析其与灾害性气候因子的关系，这对于趋利避害、规划作物生产意义重大。现阶段作物趋势产量的模拟与提取方法，主要有单一的3a、5a滑动平均法[9-10]和Logistic拟合法[11-13]，以及采用二次曲线、三次指数平滑法[8、14]。近年来，有学者采用HP 滤波等方法开展了多年长序列趋势产量拟合[15-16]，但其序列长度长达 50～60年之多。考虑到1979年之前与之后粮食产量统计技术和方法的截然不同，本研究以 1979年为起点，以保证数据的可用性和均一性，并首次采用一致性相关系数与一般相关系数以及年代分段趋势符合度等相结合的分析方法，对上述4种方法进行了综合评比，以优选出最为科学、合理、可信的作物气象产量分离方法。

1 数据来源

苏中、苏北地区是江苏的大粮仓，近40年来，其水稻总产占全省比重连年攀升，已由上世纪80年代初的约5～6成上升到 2018 年的近9成。通过不同途径收集不同年份的江苏省苏中、苏北地区代表性市(县)的水稻统计资料，其中1979～2009年的资料来自江苏省农村统计年鉴[17]，2010～2018年的资料来自所属地级市统计局官网以及相关县市年鉴或年度统计公报；以上资料均采用了统一统计口径。考虑地理分布及资料的可用性，选取了睢宁、灌南、涟水、盱眙、兴化和如东等6个市(县)统计资料开展研究，并在灌南、兴化实施了多年的水稻试验，以对统计资料加以印证。

2 研究方法

为了从技术上逐步量化气象因子与区域水稻单产之间的关系，首先要剥离历史时期长周期下不同农业生产力水平下的单产分量，求得水稻气象产量。长周期趋势产量通常肇因于强农惠农政策的落实、农业产业的投入、涉农科技进步等因素[16]。而江苏的水稻生产又有高投入、高产出的特点，例如化肥投入与水稻单产，均在国内名列前茅[18]。任一实际产量的单产时间序列，均可分解为以下3个部分，即趋势产量(社会因素导致)、气象产量(自然、灾害因素导致)和随机产量(随机噪声)，表达式如下：

(1)

2.1 HP滤波法

HP滤波法由Hodrick R J等[19]首创，可被理解为高通滤波器，它认为经济变量(如产量)既不是永远不变的，也不是随机变动的，其趋势是缓慢变动的。将产量的时间序列看作是高频成分(气象产量)与低频成分(趋势产量)的叠加。其优点是气象产量与原实际产量年份同步，故序列首尾年无缺失。水稻实际产量序列、长期趋势产量序列以及短期气象产量序列分别以ot、yt、yw表示，可将yt定义为求取公式(2)的最小化问题的解[17]，即：

(2)

yt=(100F+I)-1ot

(3)

公式(3)中：F和I分别为系数矩阵和单位阵。

2.2 二次指数平滑法

指数平滑法不同于一般的加权滑动平均法，而是遵循“厚近薄远”的原则，根据距离预测期远近对不同的历史数据给予不同的关注度；在一次指数平滑的序列上再做一次指数平滑，同时考虑偏差值以及滞后性的影响，然后建立新的趋势预测模型[20]。其优点是适应性强、结果稳定。二次指数平滑法属于再平滑，其基础为一次指数平滑法，后者的数学模型为：

S1=ao1

S2=ao2+(1-a)S1

……

St=aot+(1-a)St-1

(4)

公式(4)中，在t时，St为平滑水稻产量值，ot为水稻实际产量值，St-1为上期水稻产量平滑值，a为水稻产量平滑系数(0≤a≤1)。二次指数平滑实际是一次指数平滑的简单递归形式，其定义为：

Dt=aSt+(1-a)Dt-1

(5)

公式(5)中，St是一次指数平滑水稻产量序列，Dt是二次指数水稻平滑产量序列。

2.3 一致性相关系数与一般相关系数相结合的分析方法

首先，本文采用一致性相关系数rc[21]对区域内各市(县)的趋势产量进行测度，对不同拟合方法进行初选。

(6)

理论上讲，不同方法拟合的40年各市(县)趋势产量应与40年本区域平均趋势产量具有一致性，rc的取值范围在-1和1之间，3个具体数值范围：<0.5、0.5≤0.85和>0.85分别代表一致性差、一致性较好和一致性极好3个评价等级。本文所选市(县)均属于苏北平原代表性产粮县，同属于暖温带半湿润季风气候，社会发展和技术进步历程差异较小，趋势产量变化特征与走向应相对一致。其次，开展拟合方法决选，鉴于严重的气象灾害往往对应于较大范围的天气系统，其影响往往遍及周边较多市(县)。因此，本研究对于初选入的市(县)，对经过剥离趋势产量产生的气象产量、两两分别计算一般相关系数，以最终判断拟合方法是否最优。

2.4 水稻趋势产量与农业发展实际趋势符合度分析方法

趋势产量不仅受到品种更新、技术更新等影响，一系列制度创新和强农惠农政策的推出，也能立竿见影地解放生产力，加快趋势产量的增长。王桂芝等[15]对1961～2012 年中国粮食单产逐年资料的趋势分析表明：全阶段普遍增长，但分阶段增速差异明显。1978～1990 年间，受“家庭联产承包”等制度创新的推行以及农业技术水平等提高的影响，中国的粮食单产增长迅速，增速最大；而在1991～2012年，粮食单产随着技术不断进步保持持续增长。趋势产量是科技进步和政策调整在农业上的集中反映，趋势产量模拟曲线必须与不同社会发展阶段的农业生产实际相吻合[8]。基于全国粮食生产的历史大背景，本研究为判明4种趋势产量提取方法在本研究区域孰优孰劣，将不同方法提取的趋势产量按1979～1990年、1991～ 2003年和2004～2018年共3个年代时段，分段计算其增长速率并进行相互比较。

2.5 数据处理

分别采用 Eviews 10、SPSS 19.0 以及Excel 2016软件，对数据进行分析、处理以及作图。

3 结果与分析

3.1 水稻趋势产量、气象产量的双重相关分析

由表1所列的4种方法拟合的水稻趋势产量序列的一致性相关系数可见，所有数据数值范围为0.859～0.990，均>0.85，表明4种方法拟合的6市(县)水稻趋势产量序列与江苏区域平均趋势产量序列间的一致性极好。因此，4种方法均基本适用于宏观的粗线条的气象产量与趋势产量分离。

表1 4种方法拟合的各市(县)水稻趋势产量序列与江苏省平均趋势产量序列间的一致性相关系数

表2列出了4种方法分别拟合的不同市(县)水稻气象产量序列之间的相关系数，按平均而言，相关系数为0.49～0.51，4种方法并无明显差异，且每一种拟合方法的所有站点间的气象产量序列之间均存在不同程度的相关，而且从4种方法总体上看，临近站点之间相关系数的平均值一般大于其他站点。表2中水稻气象产量序列之间的相关系数变化基本符合区域农业气象灾害影响规律和水稻产量形成规律。但各个试点显著性检验表明：方法1、2、3中，在农业政策和技术甚为接近、天气气候条件相似的部分临近站点之间(睢宁与灌南、灌南与涟水、涟水与盱眙以及兴化与如东两两之间)，均表现为微弱正相关(P<0.1)或低度正相关(P<0.05)，这显然与生产实际不完全相符。而方法4在所有临近站点之间均达到显著相关(P<0.01)或极显著相关(P<0.001)。因此，方法4即二次指数平滑法值得推荐。

表2 4种方法拟合的不同市(县)水稻气象产量序列之间的相关系数

3.2 水稻趋势产量与社会发展实际趋势相符合状况分析

表3显示，就研究区域水稻而言，1990年前与全国趋同，而2003年是政策、价格因素导致全省进入低谷之年，而随后的2004年是水稻良种补贴、直补政策实施元年。尽管就6个站点平均而言，本研究 4种方法拟合水稻趋势产量，总体上体现了研究区域第一阶段产量迅猛增长、第二阶段产量平稳增长、第三阶段量质并重+平衡增长的特点，与王桂芝等[15]利用全国累计资料研究的粮食增长趋势大体一致，可见4种方法均能基本体现较大范围的平均趋势产量变化。但不同方法之间对研究区域不同时段水稻产量的增速的拟合结果差异明显。但是，3个阶段趋势产量分析结果还表明，前3种方法计算的部分市(县)的一部分时段出现了负值，在10年到15年的跨度内，趋势产量出现负值，这显然是不合理的，而指数平滑法模拟的上述3个阶段区域平均增速分别为279.21、88.76和74.94 kg/hm2，并且符合江苏水稻单产的区域变化规律，即1991年之前，苏中的平均单产增速明显低于苏北，而1991年之后，苏中的平均单产增速明显高于苏北[22]，与区域社会发展实际趋势相符。因此，指数平滑法更能精准地捕获区域性水稻趋势产量长期变化特征，且拟合效果最佳。

表3 4种方法计算的不同市(县)3个年段划分下的水稻趋势产量序列的增速 kg/hm2

3.3 水稻趋势产量拟合结果分析

通过三点滑动平均、二次函数、HP滤波和二次指数平滑4种方法，分别对江苏6个市(县)40年的水稻产量实施了趋势提取。结果显示，以兴化市苏北的灌南县和苏中的兴化市为例(图1)，20世纪80年代初起，由于国家一系列农业政策的改革和制度创新，极大地解放了农业从业者的生产力，同时若干主推的科技创新如高产杂交水稻的大面积推广应用，使农村劳动生产率大幅提升，因此，本阶段江苏水稻产量增速较快；而从90年代至最近5年，随着国家政策对水稻高产、优质并重的要求及其相应技术的跟进，水稻产量虽稳定增长，但增长速度呈现放缓、再放缓的态势。图1表明，4种方法提取的水稻趋势产量均能在总体上较好地体现1979～2018年区域粮食产量的增产特征。由于4种方法各自原理和计算方法各异，尽管各自拟合的6市(县)趋势产量与其实际产量皆达显著相关水平(P<0.01)，总体上趋势一致，但拟合结果也有一定差异。各地趋势产量表现出3个类型，第一个类型是迟钝型：即二次函数、HP滤波2个方法，虽然模拟的峰谷略有差异，总体上均呈现平缓向上的态势，但对剧烈产量波动年份反应迟钝。第二个类型是敏感型：即三点滑动平均，同时捕捉了较小波动年份的细微特征，显得过分敏感。而第三个类型表现为适宜型：即二次指数平滑方法，显现出气象因素以外的诸因素叠加、但波动向上的良好趋势。

图1 4种方法拟合的灌南县和兴化市趋势产量与其实际产量

4 讨论

本研究利用江苏省5个地级市6个站点40年的水稻历史产量数据，利用多项评估指标分析评价了三点滑动平均、二次函数、HP滤波和二次指数平滑4种方法分别分离出的水稻的气象(灾害)产量序列的合理性和普适性，结果显示:相对而言，4种方法中均以二次指数平滑方法最为合理、可信。有学者[8，15]分别以全国粮食单产、5个产棉大省的棉花产量为研究对象，均采用了本文4种分离方法中的某3种方法，分别取得的结果的优选顺序为：二次函数>滑动平均、HP滤波>滑动平均，从一致性相关看，这与本文的结果一致。赵东妮等[16]对辽宁省62年以及不同时段的水稻趋势产量拟合和增长速率计算结果表明：HP滤波法与指数平滑法均表现较好，所提取的水稻趋势产量均与不同时段农业生产实际相符，这与本文的结果稍有不同，这可能与资料的区域范围、农业生态条件以及产量序列长度有关。二次指数平滑方法可能更适用于市(县)级地理尺度的气象灾害的精细化评估。总体上看，通过6个站点的趋势产量分别与江苏省平均趋势产量对比可知，两两序列均表现为一致性极好，一致性相关系数达0.859～0.990。由此可见，本研究所选择的上述4种方法所拟合出的趋势产量序列，总体上与区域平均趋势产量序列一致，4种方法均基本适用于宏观的粗线条的趋势产量分离。而指数平滑法更能精准地捕获区域性水稻趋势产量的长期变化特征以及阶段性变化特征。

随着序列分段方法的成熟，未来可对作物产量资料序列进行科学分段，选取最适宜的产量分离方法分段拟合求取趋势产量，进而得到气象产量[23-24]，以及多种分离方法并举且混合使用，从而面向精细化应用，更客观地反映气象灾害对市县地理尺度作物产量的影响。

另外，由于中国农业区地域差异和气候变化显著，影响粮食产量的因素纷繁复杂，为了更准确地从趋势产量中分离出气象产量，未来必须加强以下几个方面的工作：依赖于更为详细的病、虫害影响资料，以便于在统计资料中进一步区分气象灾害的影响和病、虫灾害的影响；将基于过程的作物模拟模型引入，以相互印证；从政策上，加强统计资料的严谨性；从技术上采用更为科学的方法，以确保统计资料的真实性。