基于灰色马尔科夫模型贵州省水果产量预测方法

Prediction Method of Fruit Yield in Guizhou Province Based on Grey Markov Model

摘要: 水果在人们的膳食中占据了很重要的地位,与人们的身体健康有着千丝万缕的关系，它的营养价值十分丰富，随着国内发展越来越好，家庭对水果的需求也在逐渐增加。贵州省水果产量会受到气候、政策、种植方式、储存方式等等方面的因素的影响，使得产量会具有极强的波动性和随机性。在本文中对比了灰色系统模型和灰色马尔科夫模型对产量预测的结果，发现灰色马尔科夫模型对水果产量预测会更加精准。本文将选取贵州省2004年到2019年的水果产量数据(数据来源：贵州省2020年统计年鉴)进行分析研究。

Abstract: Fruits occupy a very important position in people’s diet and are inextricably related to people’s health. It is very rich in nutritional value, and with the development getting better and better in China, the demand for fruits by families is gradually increasing. Fruit production in Guizhou province can be affected by factors such as climate, policies, planting methods, storage methods, etc., making the output extremely volatile and random. In this paper, we compare the results of the grey system model and grey Markov model for yield prediction, and find that the grey Markov model is more accurate for fruit yield prediction. In this paper, fruit yield data which is from 2004 to 2019 in Guizhou province (data source: Statistical Yearbook of Guizhou Province in 2020) will be selected for analysis and research.

1. 引言

中国的水果产量占据了亚洲水果产量的16%，由2020年统计年鉴可知，全国水果产量27,400.8万吨占全国主要农产品产量的15%，可知水果行业在中国农业经济里占据了重要的地位；根据数据显示，中国水果产量持续增长，这也从侧面体现了中国人民在不断追求美好生活。

通过统计年鉴可知，在贵州果园种植面积从2015年的307.65 (千公顷)增加到2019年的684.5 (千公顷)，在五年时间内水果种植面积翻了1.2倍；贵州省水果种植面积在近五年内提升速度特别快，从2015年时在全国排十八上升到现在排全国第七；水果产量上，2019年比2015年翻了1.03倍；从水果种植面积和产量迅速增长可知，人们对水果的需求日益增长。但是水果产业又是一个靠天吃饭的产业，需要科学地预测水果产量，再科学地规划和调整水果生产种植方式。运用科学的方法预测水果产量，可以使果农们提前做好准备，认清楚市场形式，起到有效的结构性调整作用。

国内有不少的学者对我国或者是各个省份的水果产量进行研究，季洪霄、许峰利用灰色系统模型对中国的主要生产的几种水果进行预测，表明灰色预测的误差较小 [1]；姚飞、王波等利用灰色马尔科夫模型预测了中国未来的水果产量，指出灰色马尔科夫模型比单纯的灰色模型的精准度要高 [2]；马创、袁野等利用灰色马尔科夫模型对中国的粮食产量进行预测的准确度比传统的灰色模型以及马尔科夫模型的预测精度要高得多 [3]；邱颖利用了无偏的灰色预测模型以及马尔科夫理论对预测值进行修正的方法，对陕西省苹果年产量进行预测 [4]；在以上学者的研究中都指出了灰色系统与马尔科夫理论相结合的模型预测效果较为显著。

由于水果生产是具有明显的灰度特征的，因此，虽然可以利用灰色系统模型对其进行产量的预测，但是灰色马尔科夫预测模型预测结果更加精准。本文将选用2004~2019年的贵州省水果产量统计数据，对其进行灰色马尔科夫模型预测分析。

2. 模型介绍

2.1. 灰色GM(1,1)模型

灰色预测模型就是一种灰色系统对某一数列进行预测的方法模型 [5]。

灰色预测模型步骤为：设有一原始时间数据数列如下

$X^{\left(0\right)}\left(k\right)=\left\{x^{\left(0\right)}\left(1\right),x^{\left(0\right)}\left(2\right),\cdots ,x^{\left(0\right)}\left(n\right)\right\},k=1,2,\cdots ,n$ (1)

然后对原始数列进行一次累加，生成新的数据序列：

$X^{\left(1\right)}\left(k\right)=\left\{x^{\left(1\right)}\left(1\right),x^{\left(1\right)}\left(2\right),\cdots ,x^{\left(1\right)}\left(n\right)\right\},k=1,2,\cdots ,n$ (2)

其中 $X^{\left(1\right)}\left(k\right)={\displaystyle \underset{i=1}{\overset{k}{\sum }}{x}^{\left(0\right)}\left(i\right),k=1,2,\cdots ,n}$。

由累加算子性质可知，新生成的数据序列是近似服从于指数增加的规律，紧邻均值生成的序列如下：

$Z^{\left(1\right)}=\left\{z^{\left(1\right)}\left(2\right),z^{\left(1\right)}\left(3\right),\cdots ,z^{\left(1\right)}\left(n\right)\right\}$ (3)

其中 $Z^{\left(1\right)}\left(k\right)=\frac{1}{2}\left[x^{\left(1\right)}\left(k\right)+x^{\left(1\right)}\left(k+1\right)\right],k=1,2,\cdots ,n-1$。

由此可以建立灰色微分方程：

$\frac{\text{d}{x}^{\left(1\right)}\left(t\right)}{\text{d}t}+a{x}^{\left(1\right)}\left(t\right)=b$ (4)

其中，a是发展关系数，b是灰色作用量；其白化方程也就是影子方程为：

$x^{\left(0\right)}\left(k\right)+a{z}^{\left(1\right)}\left(k\right)=b,k=2,3,\cdots ,n$ (5)

采用最小二乘法的求解：

$\hat{A}={\left(\hat{a},\hat{b}\right)}^{\text{T}}={\left(B^{\text{T}}B\right)}^{\text{T}}{B}^{\text{T}}Y$

其中， $B=\left\{\begin{array}{cc}\begin{array}{c}-Z^{\left(1\right)}\left(1\right)\\ -Z^{\left(1\right)}\left(2\right)\\ ⋮\\ -Z^{\left(1\right)}\left(n-1\right)\end{array}& \begin{array}{c}1\\ 1\\ 1\\ 1\end{array}\end{array}\right\}$， $Y=\left\{\begin{array}{c}{x}^{\left(0\right)}\left(2\right)\\ x^{\left(0\right)}\left(3\right)\\ ⋮\\ x^{\left(0\right)}\left(n\right)\end{array}\right\}$。

将求出的a、b代入(4)中就得到灰色GM(1,1)模型的时间响应函数模型如下：

${\hat{X}}^{\left(0\right)}\left(k\right)={\hat{X}}^{\left(1\right)}\left(k\right)-{\hat{X}}^{\left(1\right)}\left(k-1\right)=\left({\hat{X}}^{\left(0\right)}-\frac{\hat{a}}{\hat{b}}\right)\ast {\text{e}}^{-ab}+\frac{\hat{a}}{\hat{b}}$

然后原始数据列的预测值需经过一次累减还原得：

${\hat{X}}^{\left(0\right)}\left(k\right)={\hat{x}}^{\left(1\right)}\left(k+1\right)-{\hat{x}}^{\left(1\right)}\left(k\right),k=1,2,\cdots ,n-1$，

其中 ${\hat{X}}^{\left(0\right)}\left(k\right)=\left\{{\hat{x}}^{\left(0\right)}\left(1\right),{\hat{x}}^{\left(0\right)}\left(2\right),\cdots ,{\hat{x}}^{\left(0\right)}\left(n\right)\right\}$。

2.2. 灰色马尔科夫模型

灰色马尔科夫模型简单来说是分成了两个部分；首先利用灰色GM(1,1)模型对所需要预测的数列的粗略地判断该数列的发展趋势，再利用马尔科夫理论对灰色模型的预测值进行修正，使得预测的精准度提高。

至于利用马尔科夫理论对灰色模型的预测值修正，首先是计算出相对误差，利用相对误差进行状态划分，通过状态转移概率矩阵对灰色预测值进行修正，达到提高准确度的目的。

第一步计算相对误差，相对误差 $\epsilon \left(t\right)$ 是绝对预测误差 $\Delta \left(t\right)$ 与实际值 $X\left(t\right)$ 的比，表示如下：

$\epsilon \left(t\right)=\frac{\Delta \left(t\right)}{X\left(t\right)}$ (6)

相对误差 $\epsilon \left(t\right)$ 是一非平稳的随机序列；然后对 $\epsilon \left(t\right)$ 进行划分状态区间，划分状态区间就是将相对误差 $\epsilon \left(t\right)$ 划分为n个状态区间，状态区间可表示如下：

$E_i=\left[e_{1i},e_{2i}\right],i=1,2,\cdots ,n$

其中， $e_{1i},e_{2i}$，分别为状态E_i的上限和下限，则相对误差的状态集合为 $E=\left(E_1,E_2,\cdots ,E_n\right)$。

马尔科夫预测模型的一个很关键步骤就是对相对误差的划分，状态划分后计算状态转移概率，状态转移概率的意思是现阶段是一种状态，下一个阶段转移到另外一个状态的概率。任一状态E_i经过m步转移到状态 $E_j$ 的转移概率，表示如下：

$p_{ij}\left(m\right)=\frac{M_{ij}\left(m\right)}{M_i}$

其中， $M_{ij}\left(m\right)$ 为样本状态E_i经过m步转移到状态 $E_j$ 的次数, $M_i$ 是状态 $E_i$ 在样本中出现的次数。则可计算出m步状态转移概率矩阵可表示为：

$P\left(m\right)=\left|\begin{array}{cccc}{p}_{11}\left(m\right)& p_{12}\left(m\right)& \cdots & p_{1n}\left(m\right)\\ p_{21}\left(m\right)& p_{22}\left(m\right)& \cdots & p_{2n}\left(m\right)\\ ⋮& ⋮& & ⋮\\ p_{n1}\left(m\right)& p_{n2}\left(m\right)& \cdots & p_{nn}\left(m\right)\end{array}\right|$

而且在状态转移矩阵中需知各行元素之和等于1；在一步状态转移矩阵中每一行的最大元素不止一个时，则需要继续计算二步状态转移矩阵甚至是计算到第n步状态转移概率矩阵，这样做只是为了使状态转移矩阵的每一行的元素的最大值只有一个，计算k步状态概率转移矩阵的计算方式如下：

$P\left(k\right)={\left[P\left(1\right)\right]}^k$

其中， $P\left(1\right)$ 为一步状态转移矩阵。

最后一步就是对灰色模型的预测值进行修正，由于预测值的修正与下一个转移状态是相关的，则当所修正的值转移到下一个状态时 $E_j$，则灰色预测值的修正公式就表示如下：

$\hat{Y}\left(k\right)=\frac{\hat{X}\left(k\right)}{1\pm 0.5\left|e_{1j}+e_{2j}\right|}$

其中 $e_{1j},e_{2j}$，分别为状态 $E_j$ 的上限和下限，当预测值比实际值高时则取“+”号进行修正，预测值比实际值低时取“−”号进行修正。

3. 实证分析

本文将基于2004~2019年的贵州省水果产量数据(数据来源：2020年贵州省年鉴)进行实证分析，数据如表1所示：

Table 1. Guizhou 2004~2019 fruit yield data

表1. 贵州2004~2019年水果产量数据

3.1. 灰色模型建立

利用上述数据以及模型进行演算，得到灰色模型GM(1,1)的估计参数发展系数−a以及灰色做用量b分别为： $-a=0.1312277,b=42.345$，因此灰色模型的预测公式为：

${\hat{x}}^{\left(0\right)}\left(k+1\right)=42.345{\text{e}}^{0.13128k},k=1,2,\cdots ,n$ (8)

得到灰色模型预测实际情况如表2所示：

Table 2. Grey model predicted values

表2. 灰色模型预测值

经过计算，方差比C = 0.23 < 0.35，小残差概率P = 1，关联度为r = 0.705 > 0.6，可以看到预测效果不错。但精度为86.33%，有空间可以提高。拟合效果图如下：

Figure 1. The fitting diagram of predicted value (▽) and actual value (×) of grey model

图1. 灰色模型预测值(▽)与实际值(×)拟合图

从上图1中可观察到在前期和后期的预测中，灰色模型预测值比实际值要偏低，中间的拟合效果较好。

3.2. 灰色马尔科夫模型建立

状态空间的划分是很重要的一步，常见的划分方法有样本空间等距法，聚类分析法、均值法以及常数划分法，由于状态划分的不同会使得状态转移矩阵会存在不同，状态划分不同会使得预测结果出现差异，状态划分不准确也会使得预测结果不准确。

通过GM(1,1)模型得到的灰色预测值，然后得到的相对残差，对相对残差进行排序。得到序列 $\epsilon \left(t\right)$ = (0, 0.0005664747, 0.00920529, 0.01567989, 0.02947272, 0.031259770, 0.03216678, 0.03687669, 0.05495186, 0.142346, 0.1534074, 0.1835321, 0.2320934, 0.3279857, 0.3998361, 0.4010769) (后面都将保留四位小数)，建立相对误差范围，若是将 $\epsilon \left(t\right)$ 是等距的区间，对于本文的数据会使得状态分布不均匀，本文利用了聚类的方法，对相对误差进行聚类如图2所示：

在图示中“single”、“complete”、“median”、“average”分别表示最短距离法、最长距离法、中间距离法、类平均法这四类聚类的计算方法，在这四张图下的聚类，都有显示，状态区间分为三组比较合适；因此，本文将状态区间划分如表3所示：

Figure 2. Clustering demonstration chart

图2. 聚类演示图

Table 3. Markov Model division

表3. 马尔科夫模型划分

因此原始数据的状态具体情况如表4所示：

Table 4. Raw data specific status display

表4. 原始数据具体状态展示

状态之间的转移状况如下表5所示：

Table 5. Statistics of numbers of state transitions

表5. 状态转移个数统计

由上述状态可以得到以下状态转移概率：

$P\left(1\right)=\left|\begin{array}{ccc}\frac{2}{3}& \frac{1}{3}& 0\\ 0& \frac{3}{4}& \frac{1}{4}\\ 0& \frac{1}{8}& \frac{7}{8}\end{array}\right|$

$P\left(2\right)={\left[P\left(1\right)\right]}^2=\left|\begin{array}{ccc}0.444& 0.472& 0.084\\ 0& 0.594& 0.406\\ 0& 0.203& 0.797\end{array}\right|$

$P\left(3\right)={\left[P\left(1\right)\right]}^3=\left|\begin{array}{ccc}0.296& 0.513& 0.191\\ 0& 0.496& 0.503\\ 0& 0.252& 0.748\end{array}\right|$

$P\left(4\right)={\left[P\left(1\right)\right]}^4=\left|\begin{array}{ccc}0.198& 0.507& 0.295\\ 0& 0.435& 0.565\\ 0& 0.282& 0.718\end{array}\right|$

3.3. 预测值修正

利用公式(7)对灰色模型的预测值进行修正，例如2005年的修正值为：

$\frac{57.47}{1-0.5\left|0.3280+0.4011\right|}=90.4398$

其它年份的具体修正情况如下表6所示：

Table 6. Comparison of the two models

表6. 两种模型的对比情况

从表6中可以看到修正后的数据更接近真实数据。将修正后的数据、实际数据以及灰色模型的预测值做拟合图如图3所示，在图3中可以直观的看到修正后的预测值更加贴合实际数据。

Figure 3. Predicted after Markov model correction (●), actual (×) and grey model predicted (▽) values after Markov model correction

图3. 马尔科夫模型修正后的预测值(●)、实际值(×)和灰色模型预测值(▽)

4. 数据预测

本文将再对2020~2023年的水果产量先利用灰色模型GM(1,1)进行预测，再利用公式对其进行修正。

对2020年的水果产量预测需要考虑2020年处于哪一个状态，见表7，由于2019年的处于状态2，2020年的水果产量处于状态1的概率为0，处于状态2的概率为3/4，处于状态3的概率为1/4；利用加权平均法取修正预测值：

Table 7. Guizhou Province fruit yield probability forecast table in 2020

表7. 2020年贵州省水果产量概率预测表

在GM(1,1)预测2020年的水果产量为411.4653 (万吨)，再通过公式对其修正，得到以下结果：

$\left(411.4653/\left(1-0.5\times 0.3744\right)\right)\times 0.75+\left(411.4653/\left(1-0.5\times 0.0609\right)\right)\times 0.25=485.7709(\; 万\; 吨\; )$

因此得到以下表8：

Table 8. Two models prediction results

表8. 两种模型预测结果

5. 结论

研究水果产量波动规律对省内农业安全问题具有一定的实际意义，本文利用2004~2019年的贵州省水果产量的数据(贵州省2020年统计年鉴)，先通过灰色模型对水果产量进行预测，然后利用马尔科夫模型对预测值修正，在修正中本文没有使用等距的方法划分状态区间，而是利用了聚类的方法对相对误差分成了三类，再继续计算得到最终结果。在本文中可知灰色模型的平均相对误差值为13.66971%，精度为86.33%；而灰色马尔科夫模型的平均相对误差值为2.87%，误差大幅度地减小了，相对精度97.13%，也比灰色模型预测的精度提高了。说明灰色马尔科夫模型能更好地预测短期水果产量。