科学、准确的粮食产量预报信息可为我国农业生产防灾减灾、合理调整作物结构、保障粮食安全等提供重要的基础数据支持。因此，产量预报工作一直是我国科研领域内的一项重中之重的研究课题^[1]。

作物产量预报方法开始以农学预报方法和统计模型预报方法为主^[2-3]，其预报时效多数在农作物将要收获之前的大约1~2个月。农学预报方法在应用上比较简单，但预报结果准确性差。统计预报方法虽然所用资料易获取，方法也相对灵活，是较有效的产量预报方法，预报准确率也较高，但其预报时效短，且有在短时间内进行筛选因子难等缺点^[4-5]。近几年，随着农业气象业务发展的需要，相关学者对作物产量动态、跟踪预报技术进行了较为深入的探讨研究，提出了基于气候适宜指数^[6-11]、关键气象因子影响指数^[12-16]、产量历史丰歉气象影响指数^[17-19]等农业气象产量动态预报方法，发展到现阶段则主要以作物生长模拟模型作为最新的研究方向之一^[20-23]。这些预报方法均建立在一定的科学理论基础之上，具有较高的预报准确率和较长的预报时效，并在省级农业气象产量预报业务服务中得到广泛应用，取得了显著成效。

但上述每种单一动态预报方法的预报结果均存在一定不稳定性，且这些研究均是在省级或国家级基础上，而随着业务发展需求的进一步提高，地区级精细化产量预报已成为我国农业气象业务中一项亟待拓展的重要服务内容。特别是山东省地区是我国重要的经济和农业大省，是冬小麦的重要生产基地，其区域冬小麦产量的及时、准确预报信息，对保障我国的粮食产量安全具有重大作用。因而，建立一种预报稳定且准确率高的地区级精细化产量动态集成预报技术方法迫在眉睫。

因此，本研究在新型统计检验聚类方法^[24-27]对山东省冬小麦种植区进行客观分区的基础上，建立各个地区级的基于作物产量历史丰歉气象影响指数、关键气象因子影响指数、气候适宜度指数以及作物生长模拟模型的冬小麦产量逐旬动态预报模型。通过分析各单一预报模型的历年预报结果平均准确率，采用加权方法，开展基于动态集成预报技术的冬小麦产量预报方法研究，以期提高各个区域冬小麦产量预报的稳定性和准确率。

以山东省17个地级市冬小麦为主要研究对象，产量资料取1980—2011年的冬小麦总产、种植面积以及平均单产。根据山东省冬小麦的种植分布，考虑到资料的完整性和连续性，选取山东省14个分布基本均匀、代表性较好的气象站：惠民、济南、莘县、长岛、龙口、威海、海阳、成山头、青岛、兖州、莒县、日照、沂源、潍坊，1980—2011年的逐日最高气温、最低气温、降水量、日照时数、风速和水汽压等要素，以及1992—2011年逐旬20 cm土壤墒情资料。发育期资料选取代表站冬小麦2001—2005年多年发育期平均值。文中所用WOFOST (world food study) 作物生长模型参数主要是参照邬定荣等^[28-29]对华北平原的研究结果，在此基础上进行参数调整。

新型统计检验聚类分析 (New Cluster Analysis of Statistical Test，CAST) 方法需要通过一定标准的统计显著性检验，因而具有明确的数理统计基础，分区结果客观唯一，可为精细化产量预报提供基础依据，具体分区方法参见文献[27, 30]。图 1为应用CAST方法分区结果，由图 1可见，C1区位于山东省的西北部，包括东营、滨州、德州、聊城和济南；C2区位于山东省中部，包括泰安、莱芜、淄博和潍坊；C3区位于山东省南部，包括菏泽、济宁、枣庄、临沂和日照；C4区包括东部沿海的青岛、烟台和威海。因此，在CAST分区基础上，以14个代表站逐日最高气温、最低气温、降水量和日照时数、风速和水汽压等气象要素为基本资料，计算区域日气象要素值、区域旬气象要素值以及区域平均积温、累积降水量、标准化降水量、累积日照时数等要素，具体计算方法参见文献[27, 30]。

图 1

图 1. 山东省冬小麦分区 Fig 1. Regionalization of winter wheat in Shandong Province

将区域冬小麦实际单产进行产量分离，分解为趋势产量和气象产量。趋势产量采用二次曲线方法模拟，即

(1)

从实际产量中分离气象产量并计算

(2)

其中，Y_t为年趋势产量，单位：kg·hm^-2；X为时间序列，即1981，1982，…，2011；a，b，c均为模拟系数；Y_w为年相对气象产量，单位：%；Y_d为年实际产量，单位：kg·hm^-2。

具体原理方法参见文献[31-33]，利用1980—2011年山东省14个气象站逐日气象资料，分别计算区域积温、标准化降水量、累积日照时数3类因子的相似距离、相关系数、综合诊断指标；在3类因子中分别选出3个最大综合诊断指标所对应的历史相似年；通过分析确定9个历史相似年冬小麦产量丰歉气象影响指数与预报年产量丰歉气象影响指数之间的关系，最终采用大概率方法^[30]计算预报年冬小麦产量丰歉气象影响指数。

利用历年冬小麦相对气象产量与旬气象因子的关系广泛筛选气象因子^[30]。判断各区域旬平均气温、旬降水量、旬日照时数与相对气象产量之间的相关性大小是否通过显著性检验。如果通过显著性检验，即可初步确定为关键气象因子；根据冬小麦生长发育的生物学特性判断该因子是否有意义，如果有意义，即确定为关键气象因子。将确定为关键气象因子的旬平均温度、旬降水量、旬日照时数按温度、降水、日照用下列方法进行因子归类：

(3)

式 (3) 中, M_tps为选定时段内同类关键因子的加权组合，R_i为选定时段内某个关键因子与气象产量的相关系数，E_i为选定时段内的关键气象因子。建立相对气象产量与综合关键气象因子影响指数的线性回归模型，进行冬小麦产量动态预报。

适宜度模型的具体构建方法参考文献[27, 30]。根据山东省冬小麦各生育期所需的最低温度、最高温度和最适温度^{[27, 34]}，在计算各站逐日温度适宜度的基础上，采用算数平均法分别计算逐站逐旬温度适宜度和区域逐旬温度适宜度。

采用算术平均方法对14个气象站旬降水量进行处理求得区域旬降水量 (区域旬日照时数处理方法类似)，旬降水距平百分率大于30%时，区域旬降水适宜度值为1981—2010年区域旬平均降水量与旬降水量的比值；旬降水距平百分率小于-30%时，区域旬降水适宜度值为旬降水量与1981—2010年区域旬平均降水量的比值；其他情况下区域旬降水适宜度值为1。单站20 cm土壤相对湿度小于冬小麦各个发育阶段适宜相对湿度值时，土壤墒情适宜度为20 cm土壤相对湿度与适宜相对湿度的比值；其他情况下土壤墒情适宜度为1；区域土壤墒情适宜度为区域内所有站点土壤墒情适宜度的算数平均值。在降水适宜度和土壤墒情适宜度的基础上，确定冬小麦水分适宜度：越冬期没有土壤墒情观测记录，用越冬前最后一次土壤墒情适宜度；冬小麦抽穗前用土壤墒情适宜度；冬小麦抽穗后，当区域某旬降水距平百分率大于30%时，用降水适宜度，其他情况用土壤墒情适宜度。

当区域旬日照时数不小于区域旬平均日照时数 (1981—2010年平均值) 时，区域旬日照适宜度的值为1；当区域旬日照时数小于区域旬平均日照时数时，区域旬日照适宜度的值为区域旬日照时数与区域旬平均日照时数的比值。

综合考虑气温、日照、水分适宜度的冬小麦气候适宜度为

(4)

式 (4) 中，F(c)_k，F(t_k)，F(m_k)，F(s_k) 分别为第k旬的气候、温度、水分和日照适宜度。逐旬适宜度的加权集成构成了不同时段的适宜度指数，权重系数通过计算区域各旬气候适宜度与该区域冬小麦相对气象产量的相关系数及其所占比例确定。将各区域冬小麦不同生育时段气候适宜度指数与相对气象产量作相关分析 (表 1)，两者的相关性均达到0.01显著性水平。因此，利用冬小麦相对气象产量和不同生育时段的气候适宜度指数进行回归，建立不同生育时段的冬小麦产量动态预报模型。

表 1

表 1 冬小麦不同生育时段气候适宜度指数与气象产量的相关系数 Table 1 Correlation coefficient of meteorological yield and climatic suitable index in different developmental stages of winter wheat 区域 播种至1月下旬 播种至2月下旬 播种至3月下旬 播种至4月下旬 播种至5月下旬 C1 0.713 0.704 0.704 0.721 0.714 C2 0.703 0.747 0.776 0.797 0.770 C3 0.716 0.736 0.737 0.765 0.781 C4 0.566 0.586 0.743 0.766 0.797

表 1 冬小麦不同生育时段气候适宜度指数与气象产量的相关系数 Table 1 Correlation coefficient of meteorological yield and climatic suitable index in different developmental stages of winter wheat

WOFOST模型通过不断变换气象资料初始场进行实时滚动预报，模型的气象初始场分别是作物播种后各旬的最后一天。气象初始场之前的数据采用实际气象资料，此后采用平均气候条件下的气象资料驱动作物模型。依次类推，通过驱动模型得到历年和预报年的潜在和雨养产量，这里假设模型模拟得到的潜在产量为灌溉地产量，雨养产量为非灌溉地产量。

采用冬小麦相对气象产量与最近5年平均产量通过式 (5) 的方法对实际冬小麦产量加以订正，去除产量长序列中由于科技进步的影响，得到历年订正后的实际产量

(5)

式 (5) 中，Y为2007—2011年5年平均产量, 单位：kg·hm^-2；Y_w为通过式 (1) 和式 (2) 二次曲线方法处理得到的历年相对气象产量。订正后的历年实际产量Y、灌溉地产量Y₁、雨养地产量Y₂三者之间的关系：

(6)

式 (6) 中，A为灌溉地面积比例系数，且0≤A≤1。

利用最小二乘法，拟合历年订正后的实际产量Y与灌溉地产量Y₁、雨养地产量Y₂之间的关系式 (6)，可得到灌溉地面积比例系数A。通过模型输出的预报年灌溉地产量和雨养地产量以及灌溉地面积比例系数A通过式 (6) 得到预报年订正后的实际产量Y，带入式 (5) 得到预报年的相对气象产量Y_w，回代至式 (1) 和式 (2) 可得到预报年的预报单产。

记作物产量历史丰歉气象影响指数、关键气象因子影响指数、气候适宜度影响指数、WOFOST作物生长模型4种方法在第j年的预报单产分别为Y_{1, j}，Y_{2, j}，Y_{3, j}和Y_{4, j}。集成方法第j年的预报值为

(7)

式 (7) 中，Y_j为基于集成预报方法第j年的预报值，Y_{k, j}为第k种预报方法在第j年的预报值。λ_k为第k种方法的权重，为该种方法的历年预报平均准确率占所有方法历年 (预报年第j年往前推5年，下同) 预报平均准确率之和的比值。如果某区域某一旬某种预报方法的历年预报平均准确率低于90%，则加权时不考虑该方法。如果4种预报方法的历年预报平均准确率均低于90%，则用4种方法的预报结果进行集成。另外，要考虑气候适宜度指数预报方法中气候适宜指数应与产量呈正相关，若不符合，则集成时不考虑该方法 (表 2)。

表 2

表 2 4种方法不同预报时间的权重系数 Table 2 The weight coefficient of four methods at different forecasting time 区域 预报时间 产量历史丰歉 关键气象因子 气候适宜度指数 WOFOST模型 C1 1月下旬 0.4827 0.0000 0.0000 0.5173 2月下旬 0.4803 0.0000 0.0000 0.5197 3月下旬 0.4833 0.0000 0.0000 0.5167 4月下旬 0.4832 0.0000 0.0000 0.5168 5月下旬 0.5047 0.0000 0.0000 0.4953 C2 1月下旬 0.2492 0.2422 0.2495 0.2591 2月下旬 0.3301 0.3188 0.0000 0.3510 3月下旬 0.3261 0.3279 0.0000 0.3460 4月下旬 0.3312 0.3319 0.0000 0.3369 5月下旬 0.3370 0.3338 0.0000 0.3292 C3 1月下旬 0.0000 0.3346 0.3329 0.3325 2月下旬 0.4849 0.0000 0.0000 0.5151 3月下旬 0.3353 0.0000 0.3181 0.3466 4月下旬 0.3393 0.0000 0.3286 0.3321 5月下旬 0.3351 0.0000 0.3277 0.3372 C4 1月下旬 0.3361 0.3307 0.0000 0.3332 2月下旬 0.3275 0.3363 0.0000 0.3363 3月下旬 0.0000 0.4977 0.0000 0.5023 4月下旬 0.3218 0.3437 0.0000 0.3345 5月下旬 0.3233 0.3392 0.0000 0.3375  注：以2009年为例。

表 2 4种方法不同预报时间的权重系数 Table 2 The weight coefficient of four methods at different forecasting time

本研究利用动态滚动预报的方法，即利用预报年之前的历史资料预报当年的产量，预报年改变，则历史资料年代随之滚动改变，而不是基于一段固定年份建立不变的模型。基于产量历史丰歉气象影响指数和集成技术的产量动态预报方法没有固定模型，其历史资料回代检验结果即为外推预报结果。从其他几种方法不同建模年份产量动态预报模型历史回代检验平均准确率 (表 3) 可以看出，基于关键气象因子影响指数方法普遍达93.0%以上；基于气候适宜度指数方法，平均准确率相对较高，普遍在95.0%以上；而基于WOFOST模型方法，平均准确率相对较低，但均在90.0%以上。总体上看，不同方法不同建模年份产量动态预报模型历史回代统计检验结果较好。

表 3

表 3 不同建模年份产量动态预报模型回代检验平均准确率 (单位：%) Table 3 Average accuracy of return test for the dynamic yield forecast model based on different years (unit:%) 预报方法 建模年份 C1区 C2区 C3区 C4区 关键气象因子 1980—2003 95.4 92.5 96.6 94.7 1980—2004 95.2 92.5 96.6 94.1 1980—2005 95.2 92.4 96.0 94.2 1980—2006 95.7 93.3 95.7 94.5 1980—2007 96.0 93.4 96.5 95.1 1980—2008 95.6 93.0 95.5 95.0 1980—2009 94.9 93.5 94.4 94.7 1980—2010 94.9 93.6 94.8 94.4 气候适宜度指数 1992—2003 98.8 97.0 97.0 97.1 1992—2004 98.4 97.0 96.6 96.9 1992—2005 97.7 96.7 96.5 95.8 1992—2006 97.3 96.4 96.1 95.5 1992—2007 97.9 96.2 96.1 95.7 1992—2008 97.1 95.8 95.6 95.8 1992—2009 96.5 95.5 95.7 95.9 WOFOST模型 1980—2003 91.2 91.6 90.8 91.9 1980—2004 91.4 92.0 92.3 90.7 1980—2005 91.6 92.2 92.4 90.6 1980—2006 91.9 92.2 91.6 92.9 1980—2007 94.2 93.5 92.9 90.9 1980—2008 92.0 92.8 92.6 91.3 1980—2009 91.6 92.9 92.8 91.8 1980—2010 93.7 93.7 92.2 91.9  注：预报为5月下旬。

表 3 不同建模年份产量动态预报模型回代检验平均准确率 (单位：%) Table 3 Average accuracy of return test for the dynamic yield forecast model based on different years (unit:%)

本研究建立的冬小麦产量动态预报模型，可从冬小麦播种开始逐旬做产量预报。但限于文章篇幅，文中仅列出2004—2011年1—5月各个月下旬的产量预报结果。

由图 2可知，C1区2004—2011年冬小麦产量动态预报准确率普遍较高，除了2007年的1月、2月和3月下旬未超过90.0%，其他均超过90.0%，预报准确率最高值出现在2004年5月下旬，达到98.8%。C2区的预报准确率不如C1区，除2004年、2006年、2009年和2010年个别旬超过90.0%，多数未达到90.0%，2005年1月下旬和2011年5月下旬预报准确率甚至低于80.0%。C3区的预报准确率也普遍较高，多数为90.0%以上。C4区冬小麦产量预报准确率相对较低，波动幅度较大，准确率最大值达97.1%，最小值低至73.2%，且没有达到90.0%的情况较多。说明基于产量历史丰歉气象影响指数的产量动态预报方法在C1区和C3区比在C2区和C4区更适用。

图 2

图 2. 基于产量历史丰歉气象影响指数方法预报结果 Fig 2. Forecasting result of method based on meteorological influence index of yield historical harvest conditions

由图 3可知，各区域基于关键气象因子影响指数方法的冬小麦产量动态预报准确率波动幅度较大。C1区在同一年各时段预报准确率差异相对较小，2004年、2008年和2009年各个时段预报准确率均较低，普通低于90.0%，最低值为79.1%，出现在2004年2月下旬。C2区在同一年各时段预报准确率差异也较小，在2008年之前以及2011年各个时段预报准确率普遍较高，均在90.0%以上，最高值达100.0%；但2008年和2010年各时段预报准确率较低，没有达到90.0%。C3区各时段预报准确率差异稍大，最高值在2007年3—5月下旬，达到99.9%，最低值在2008年3—5月下旬，为85.5%。C4区在2006年、2007年、2009年以及2011年各时段预报准确率均较高，最大值达98.6%，在2005年2月下旬预报准确率最低，为83.9%。说明这种方法在各区域的特定时段比较适用。

图 3

图 3. 基于关键气象因子影响指数方法预报结果 Fig 3. Forecasting result of method based on influence index of key meteorological factors

各区域各时段基于气候适宜度指数的产量预报准确率普遍在80.0%以上 (图 4)。C1区在同一年不同时段的预报准确率差异较小，2010年和2011年各个时段预报准确率均在96.0%以上，在2010年2月下旬几乎达到100.0%，其他时段预报准确率相对稍低。C2区在2006年、2008年以及2009年的个别时段预报准确率偏低，最低值为84.9%。C3区预报准确率呈波浪式变化，预报准确率最大值在2010年1月下旬达99.9%，最小值在2008年1月下旬，为83.9%。C4区在2005年、2006年和2009年各时段预报准确率普遍偏低，均低于90.0%。

图 4

图 4. 基于气候适宜度影响指数方法预报结果 Fig 4. Forecasting result of method based on climatic suitable index

2004—2011年1—5月下旬冬小麦产量动态预报结果中 (图 5)，总体上，C1和C2区的预报效果比C3和C4区相对更好。C1区只有2004年的1月下旬和2月下旬，2007年的3月下旬和4月下旬以及2010年5月下旬预报准确率低于90.0%。C2区预报准确率相对比较稳定，波动幅度较小，最高值在2005年1月下旬为99.8%，最低值在2008年1月下旬为88.7%。C3区在2007年3月下旬预报准确率最低为67.8%，其他各时段预报准确率普遍较高。C4区各时段的预报准确率在2006年之前的变化趋势一致，到2007年1月下旬、2月下旬和5月下旬的预报准确率较高，而3月下旬和4月下旬的预报准确率偏低，最高值、最低值相差10.0%左右，2009年以后，除了5月下旬其他各个时段的变化趋势恢复一致，且普遍较高。

图 5

图 5. 基于WOFOST作物模型预报结果 Fig 5. Forecasting result of WOFOST crop model

由图 6可知，集成预报方法的冬小麦产量预报稳定性与其他4种方法相比，均明显提高，预报准确率波动范围较小，且准确率普遍较高。C1区冬小麦产量预报准确率最低值出现在2009年1月下旬和2月下旬，分别为94.4%和93.4%，其他大多数预报准确率普遍在95.0%以上，准确率在98.0%以上的情况较多。C2区冬小麦平均单产预报准确率多数在95.0%以上，最高值达到100.0%，出现在2004年1月下旬；最低值在2008年3月下旬，为91.3%。C3区预报准确率普遍较高，几乎均在95.0%以上，其中，在2005年3月下旬和2007年5月下旬达到100%。C4区除了2004年1月下旬预报准确率未达到90.0%，其他各时段普遍在90.0%以上，最高值达到100.0%，出现在2005年4月下旬。说明这种集成预报方法在各区域大多数时段均适用，可达到较好的预报效果。

图 6

图 6. 基于集成预报技术方法预报结果 Fig 6. Forecasting result of method based on integrated forecasting technology

2004—2011年各区域不同时段5种产量动态预报方法预报平均准确率如表 4所示。整体上，集成预报方法在各区域预报平均准确率均达到95.0%以上，普遍比其他4种方法高。C1区WOFOST作物模型和产量历史丰歉气象影响指数方法的预报结果较好。C2区产量历史丰歉气象影响指数方法的预报结果最差，而WOFOST作物模型的预报结果次于集成预报方法。C3区各种方法预报准确率相差不明显。C4区产量历史丰歉气象影响指数方法预报效果在1月下旬仅为85.1%，预报效果相对较差。

表 4

表 4 2004—2011年山东冬小麦产量预报平均准确率 (单位：%) Table 4 Average accuracy of yield forecast of winter wheat in Shandong Province from 2004 to 2011(unit:%) 区域 预报时间 产量历史丰歉 关键气象因子 气候适宜度指数 WOFOST模型 集成预报方法 C1 1月下旬 93.7 92.0 90.6 94.7 96.9 2月下旬 93.9 90.4 90.1 96.0 97.1 3月下旬 93.3 89.9 90.0 95.4 97.6 4月下旬 94.2 89.3 90.0 94.9 97.3 5月下旬 94.5 89.1 90.0 93.0 98.5 C2 1月下旬 88.2 92.9 93.4 95.4 98.2 2月下旬 86.8 93.4 92.1 96.6 97.7 3月下旬 87.3 93.8 92.5 96.6 97.5 4月下旬 86.7 93.2 92.7 95.0 96.5 5月下旬 85.8 93.4 93.1 94.4 96.6 C3 1月下旬 90.2 93.6 94.4 93.4 96.5 2月下旬 92.2 94.2 93.3 95.8 98.0 3月下旬 93.5 93.2 92.4 92.8 98.1 4月下旬 94.6 93.9 92.7 94.6 98.4 5月下旬 94.6 94.1 93.0 92.9 97.9 C4 1月下旬 85.1 93.0 90.3 92.2 95.1 2月下旬 86.4 91.3 90.3 93.6 96.2 3月下旬 86.9 91.7 90.1 93.6 95.9 4月下旬 87.0 92.1 90.7 93.3 96.3 5月下旬 85.7 92.0 91.3 94.2 96.3

表 4 2004—2011年山东冬小麦产量预报平均准确率 (单位：%) Table 4 Average accuracy of yield forecast of winter wheat in Shandong Province from 2004 to 2011(unit:%)

本文在CAST分区基础上分析了作物产量历史丰歉气象影响指数、关键气象因子影响指数、气候适宜度影响指数、WOFOST作物生长模型以及集成预报方法不同冬小麦产量动态预报方法在山东省各区域的预报效果。研究表明：

1) 单一产量预报方法在进行滚动预报时，其预报结果在某些区域某一时段的预报效果相对较好，预报准确率可达90.0%以上，而在某些情况下预报效果较差。单一预报方法对因子的选择和参数的考虑具有一定局限性，预报结果不稳定，波动性较大。

2) 集成预报方法利用历年平均预报准确率加权集成法，在某一特定区域特定时段，如果某一种产量预报方法的预报效果不好，则集成时不考虑该方法。其单产预报准确率相对于其他4种方法均有所提高，准确率普遍在95.0%以上，且预报结果稳定性较好。说明该方法克服了单一预报方法预报结果稳定性差的不足，可广泛应用。

本研究将产量预报区域降到省级以下，并实现了在冬小麦播种后任意时段进行动态滚动预报，且采用了稳定性较好的集成预报方法。但本研究仅针对山东省冬小麦，对于其他省份其他作物，该方法是否实用，需进一步探讨。