【期市月月谈】植被指数与美豆单产的相关性探讨

　　安如泰山

　　信守承诺

　　摘要

　　7 月末-8 月末的 NDVI、EVI 均与美豆单产有较强的正相关关系，其中 EVI 的相关 系数高于 NDVI。在以 NDVI、EVI 分别和美豆单产建立的一元回归模型中，7/27-8/12 时间段内的 NDVI 表现最优，且对单产变化具有一定的解释力度。因此在进行单产预测 时，可以考虑将 7 月末-8 月末这一时间段内的 NDVI 或 EVI 作为因子之一。而简单的一元线性回归拟合度仍有提升空间，且出现预测值普遍偏低的问题，因此在下一次的期市月月谈中，我们将拓展单产预测模型，并进行更多的讨论。

　　正文

　　北美大豆的收获正在进行中，虽然定产预估需等待明年 1 月 USDA 报告公布，但是新作产量的波动对市场的影响已在盘面得到反映。2022 年北美干旱问题贯穿始终，8 月下旬 Pro Farmer 的田野调查显示多数州的大豆、玉米单产值都显著低于 3 年平均水平，数据公布当周 CBOT 大豆、玉米、小麦合约均表现出较大程度的涨幅。单产的变化对市场的影响无疑，因此利用遥感数据来判断单产预估将是本文探讨的重点。

　　在上一期的遥感专题《另眼看风景——以植被指数为例初探遥感数据在农产品研究中的应用》中，我们谈到了关于 NDVI、EVI 的基本概念、季节性特征，以及利用植被指数进行土地覆盖类型分类的问题。在此基础之上，本文将研究聚焦到具体的品种（美国大豆），选取了 9 个主产州（Iowa, Illinois, Indiana, Ohio, Missouri, Minnesota, Nebraska, South Dakota, North Dakota）作为对象，进一步探讨植被指数和单产的相关性。

　　首先介绍文中使用的数据来源。第一，美国大豆种植区域的数据来自 USDA NASS CDL（Cropland Data Layer），该数据覆盖时间自 1997 年至 2021 年，结合卫星数据和实地数据，利用监督学习的方式给出某自然年度内不同农作物的种植区域（参见https://www.nass.usda.gov/Resear ... pland/sarsfaqs2.php）。因为 2022 年数据尚未提供，故分析时以 2021 年数据替代。第二，NDVI、EVI 数据来自 LP DAAC MOD13Q1 （ 6.1 版 本 ） 250 米 16 天 合 成 数 据 （ 相 关 介 绍 参 见 https://lpdaac.usgs.gov/products/mod13q1v061/）。第三，美国各州单产数据来自 USDA NASS 数据库。

　　大豆的生长期 NDVI 曲线特征

　　以 2021 年美国大豆主产州的 NDVI 时间序列为例，通过 NDVI 的变化来反映大豆生长期的过程。从图一中能总结出两个较明显的特点：一、从时间上看，自 5 月开始 NDVI 的增长速度变快，其最大值一般发生在 8、9 月份，这代表了大豆从开始种植到结荚鼓粒成熟的时期。从 9 月之后，NDVI 开始下降，且 9-10 月的下降速度较快，反映大豆已进入收获期。到 11 月收获近结束时，NDVI 下降速率已然趋缓。二、从空间上看，在 2021 年，位置偏北的各州，例如 North Dakota，South Dakota，Minnesota 等的 NDVI 整体上低于位置偏南的各州，例如 Indiana, Ohio, Missouri 等，且其到达峰值、开始下降的时间点更早一些。

　　从图 3 各州的 NDVI 季节性图上则能反映出更多细节问题。

　　首先，NDVI 值是否对单产的预测具备意义。诸多文献资料已经针对 NDVI 在预测单产上的作用方面有过研究并得出积极的结论。这里我们定性地观察一下，从 North Dakota 的 NDVI 季节性图上找到一条明显低于其他年份的曲线，即 2021 年的 NDVI。而 2021 年 North Dakota 的单产为 25.5 蒲式耳/英亩，是自 2005 年以来的最低值，这一点符合 NDVI 的表现。但这样较极端的情况比较少见，且有可能是异常点，因此我们将 9 个主产州 2000- 2021 年的单产和 NDVI 数据进行了相关性分析，具体讨论见下一节。

　　其次，生长期内不同阶段 NDVI 的重要性。图 3 中多个主产州在 2019 年的 NDVI 曲线都明显向右偏移，说明该年度大豆种植的时间处于过去 6 年的偏晚位置。结合 2019 年大豆生长期的 cropprogress 数据（图 4 下），可以看到，2019 年从种植到收获的进度都是明显低于 2018 年和过去 5 年平均值，这一点与 NDVI 的表现非常吻合。回顾 2019 年美国大豆的种植情况，降雨在其中的影响非常显著。从图 5 可以看出，从 5 月中下旬开始，2019 年的降雨量明显处于近 5 年高位，在种植窗口期多个地区的降雨量几乎翻倍。降雨对作物单产的影响可以从优良率反映出来（图 4 上），2019 年整个生长期的优良率都位于 5 年低位。再结合 NDVI 的曲线，可以发现，尽管曲线反映出生长期向后推迟的特征，但 NDVI 最大值的数值变化却不太明显，所以只使用生长期内最大 NVDI 和单产做相关性分析还不足以解释单产的变化。因此我们将生长期内各个月份的 NDVI 分别和单产做相关性和一元线性回归，试图找出关键时期的 NDVI，具体讨论见下一节。

　　植被指数与美豆单产的相关性

　　接下来讨论植被指数与大豆单产的关系。数据的时间为 2000 年至 2021 年，范围包括 9 个大豆主产州，选取大豆生长期内（4 月-10 月）的植被指数和当年单产为研究对象。因为数据量较为有限，此处将不同州、不同年份的数据混置在一起研究，并未考虑不同地区或者不同年份的影响。

　　从 NDVI 与单产的散点图（图 6）中，可以发现 7 月 27 日至 8 月 12 日、8 月 13 日 至 8 月 28 日这两段时间的 NDVI 与单产具有较为明显的正向相关关系，分别为 0.606、0.648。若以 EVI 作为观察对象，也得到类似的结果，只是相关系数较高的时间段多了 7 月 12 日-27 日。具体的统计结果见表 1。

　　在和单产的相关性比较中，NDVI 和 EVI 的表现一致，同时 EVI 与单产的相关性略强于 NDVI。接下来，我们利用 2000 年至 2021 年 7 月 27 日-8 月 12 日以及 8 月 13 日28 日的 NDVI、EVI 分别与当年的单产建立一元回归模型，并将 2022 年的数据作为测试集。从结果来看，以 8/13-8/28 时间段内的植被指数为解释变量的训练集拟合效果最佳，但从测试集结果来看，7/27-8/12 时间段内的效果略好，且以 NDVI 作为解释变量的效果好于 EVI。从图 7 测试集的真实值和拟合值来看，有一个点的拟合效果不佳（由于数据集有限，无法判断该点是否是异常值）。该点代表的植被指数显著低于其他数据点，但是单产 却没有明显的降低，由此仅依靠植被指数来拟合单产无法很好地解释这类数据。从测试集的效果来看，2022 年单产的拟合值有一定的解释度，但明显比真实值偏低。因为 2022 年的单产非定产数据，且 CDL 沿用 2021 年的数据，因此和训练集的处理方式有一定差异，这可能为结果带来了偏差。另外，回归的效果和训练集和测试集的分组也有关系，此处将 2022 年的数据作为测试组，2000-2021 年的数据作为训练组，得到的可能并不是最好的回归效果，可以考虑使用 k-fold 交叉验证进一步优化。

　　结论

　　7 月末-8 月末的 NDVI、EVI 均与美豆单产有较强的正相关关系，其中 EVI 的相关 系数高于 NDVI。在以 NDVI、EVI 分别和美豆单产建立的一元回归模型中，7/27-8/12 时间段内的 NDVI 表现最优，且对单产变化具有一定的解释力度。因此在进行单产预测 时，可以考虑将 7 月末-8 月末这一时间段内的 NDVI 或 EVI 作为因子之一。而简单的一元线性回归拟合度仍有提升空间，且出现预测值普遍偏低的问题，因此在下一次的期市月月谈中，我们将拓展单产预测模型，并进行更多的讨论。