加权马尔可夫链模型对上海地区年降水量的预测

Prediction of Annual Precipitation in Shanghai by Weighted Markov Chain Model

摘要: 采用均值–标准差分级法，以上海地区1970~2020年的年降水量数据为样本降水序列，根据上海地区降水量特点，确定了样本降水序列的分级标准和状态。根据马尔可夫理论和统计学原理，验证了样本降水序列满足马尔可夫性(马氏性)，进而以规范化的各阶自相关系数为权重，建立了适用于该地区降水量的加权马尔可夫链的预测模型。以此模型预测了上海地区2021年和2022年的年降雨量，预测结果比较精确。

Abstract: Using mean-standard deviation classification method, taking annual precipitation data of Shanghai from 1970 to 2020 as sample precipitation series, according to the characteristics of precipitation in Shanghai, the classification standard and status of sample precipitation series were determined. According to Markov theory and statistical principle, it is verified that the sample precipitation series satisfies Markov property. Then, the weighted Markov chain prediction model suitable for the precipitation in this region is established by taking the normalized autocorrelation coefficients of each order as the weight. This model predicts the annual rainfall of Shanghai in 2021 and 2022, and the prediction results are relatively accurate.

1. 引言

降水量对农业生产及重大旱涝灾害的产生起着十分重要的作用。但气象条件具有复杂性、多样性和不可预知性，导致降水量在不同地区和不同时间的不均衡性。因此精确的降水量的预报对农业生产和旱涝灾害的控制非常重要。自从1999年冯耀龙等 [1] 首次在我国提出加权马尔可夫链预测的概念以来，国内很多学者利用加权马尔可夫链预测模型成功的对不同地区降雨量进行了预测 [2] - [10] 。

加权马尔可夫链预测模型是一个步骤清楚、计算方便且客观可靠的方法，值得深入研究和推广。本文从上海气象局官网选取了上海地区1970~2020年51年的年降水量数据为样本降水序列，采用均值-标准差分级法，根据上海地区降水量特点，建立了适用于上海地区年降雨量的分级标准，并用加权马尔可夫链对上海地区的2021和2022年的年降雨量状态进行了预测。预测结果所在的状态跟真实的降水量所在的状态一致。再结合模糊集理论的级别特征值，对降水量做了具体预测。结果显示跟真实值相比，平均误差为4.5%，预测精度较高。

2. 加权马尔可夫模型

马尔可夫过程是一类随机过程，是研究事物的状态和状态转移规律的理论。由俄国数学家A. A. 马尔可夫于1907年提出。它通过不同状态的初始概率和状态之间的转移概率来确定状态的变化趋势，从而达到预测的目的。马尔可夫链具有无后效性特征，这也被后人称为马尔可夫性(“马氏性”)。马尔科夫链预测模型的基本原理。先分别依据前面若干年降水量所对应的状态对某时段的状态做预测，然后按照前面各时段与该时段相关关系的强弱对转移概率进行加权求和。

3. 计算步骤

3.1. 样本分组(均值–均方差法)

上海地区1970~2020年的年平均降水量数据 $\left\{x_1,x_2\text{,}\cdots \text{,}{x}_n\right\}$ 的样本均值和样本均方差分别为 $\bar{x}=1253.7451$ ， $s=1253.7451$ 。以样本均值 $\bar{x}$ 为中心，将数据序列分成如下五组(本均值–均方差(标准差)分组法)

$\left(-\infty ,\bar{x}-{\alpha }_{1}s\right),\left[\bar{x}-{\alpha }_{1}s,\bar{x}-{\alpha }_{2}s\right],\left[\bar{x}-{\alpha }_{2}s,\bar{x}+{\alpha }_{3}s\right],\left[\bar{x}+{\alpha }_{3}s,\bar{x}+{\alpha }_{4}s\right],\left[\bar{x}+{\alpha }_{4}s,+\infty \right).$

根据上海地区降水量特点，取 ${\alpha }_1={\alpha }_4=1.340,{\alpha }_2={\alpha }_3=0.185$ 。依据上面的样本均值–均方差分组法将1970~2020年的年降水量分为5个等级(组)，分别为丰涝年、偏丰年、正常年、偏旱年、干旱年，具体分级情况见下表1。

Table 1. Scale of annual precipitation

表1. 年降水量分级表

根据表1的分级标准，定义将1970~2020的降水量所对应状态如表2 (因表格大小限制，下表中每一行列出7年的状态，每行从左到右分别是7年对应的状态)。

Table 2. Annual precipitation status table

表2. 年降水量状态表

3.2. “马氏性”检验

采用 ${\chi }^2$ (卡方分布)来检验该序列是否具有“马氏性”。设m为状态数，即m = 5。对服从自由度为 ${\left(m-1\right)}^2$ 的 ${\chi }^2$ 分布，给定一个显著性水平 $\alpha$ ，查表得 ${\chi }_{\alpha }^2\left({\left(m-1\right)}^2\right)$ 的值。如果 ${\chi }^2>{\chi }_{\alpha }^2\left({\left(m-1\right)}^2\right)$ ，则拒绝原假设，即认为序列具备“马氏性”。反之，则序列不具备“马氏性”。 ${\chi }^2$ 统计量为：

${\chi }^2=2{\displaystyle \underset{i=1}{\overset{m}{\sum }}{\displaystyle \underset{j=1}{\overset{m}{\sum }}{f}_{{}_{ij}}\cdot \left|\ln \frac{p_{ij}}{p_{•j}}\right|}}$

其中概率 $p_{ij}\text{=}\frac{f_{ij}}{{\displaystyle \underset{i=1}{\overset{m}{\sum }}{\displaystyle \underset{j=1}{\overset{m}{\sum }}{f}_{ij}}}}$ ，边际概率 $p_{•j}=\frac{{\displaystyle \underset{i=1}{\overset{m}{\sum }}{f}_{ij}}}{{\displaystyle \underset{i=1}{\overset{m}{\sum }}{\displaystyle \underset{j=1}{\overset{m}{\sum }}{f}_{ij}}}}$ ( $f_{ij}$ 表示本年是状态i而经过一年后是状态j的年份数量)。如果 $p_{ij}=0$ ，则 $\ln \frac{p_{ij}}{p_{•j}}=-\infty$ ，且 $f_{ij}=0$ ，规定 $0\times \infty =0$ 。

3.3. 计算各阶的相关系数和权重

令 $r_k$ 表示第k阶(k为步长(滞时)， $k\in E=\left\{1,2,3,4,5\right\}$ )的自相关系数。其公式为：

$r_k=\frac{{\displaystyle \underset{i=1}{\overset{n-k}{\sum }}\left(x_i-\bar{x}\right)\left(x_{i+k}-\bar{x}\right)}}{{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(x_i-\bar{x}\right)}^2}},k\in E,n=51.$

定义权重： $w_k\triangleq \frac{\left|r_k\right|}{{\displaystyle \underset{k\in E}{\sum }\left|r_k\right|}}.$ 故 ${\displaystyle \underset{k\in E}{\sum }{w}_k=1}$ 且 $w_k\ge 0$ 。

3.4. 利用各阶转移概率矩阵预测状态

以前面若干时段的状态作为初始状态，利用其相应的各阶转移概率矩阵，即可得该时段的概率 $p_{{}_{ij}}^{(k)}=0$ ，k为步长。将同一状态的各时段概率加权和作为处于该状态的概率 $P_{\left(j\right)}={\displaystyle \underset{k=1}{\overset{5}{\sum }}{w}_k{p}_{ij}^{\left(k\right)}}$ 。这些概率中最大的 $P_{\left(i\right)}=\underset{1\le j\le 5}{\max }\left\{p_{\left(j\right)}\right\}$ 所对应的状态即为该时段的预测状态。

3.5. 利用模糊集理论中的级别特征值预测年降水量的值

根据状态概率可得到当年降雨量预测值 $X_{n+1}$ ，依据模糊集理论中的级别特征值求法，对各状态定义相应的权重 $W_i=p_i^{\eta }÷{\displaystyle \underset{i=1}{\overset{5}{\sum }}{p}_i^{\eta }},1\le i\le 5$ ， $\eta$ 为最大概率的作用系数，通常取为2或4。级别特征值 $H={\displaystyle \underset{i=1}{\overset{5}{\sum }}i\cdot W_i}$ 。如果年降雨量预测状态为i，则年降雨量预测值为 $X_{n+1}=\frac{H{D}_i}{i}$ ，期中 $D_i$ 为状态i区间值的下限。

上面公式在应用时结果并不十分理想，因此本文在大量计算的基础上对公式进行了分类改进 [9] 。

$\{\begin{array}{l}{X}_{n+1}=\frac{H{D}_i}{i-\text{1}\text{.2}}\left(i=4\right)\hfill \\ X_{n+1}=\frac{H{D}_i}{i+0.2}\left(i=3,H\ge 3.6\right)\hfill \\ X_{n+1}=\frac{H{D}_i}{i-0.\text{2}}\left(i=3,H<3.6\right)\hfill \\ X_{n+1}=\frac{H{D}_i}{i+\text{1}}\left(i=2,H>2.8\right)\hfill \\ \begin{array}{l}{X}_{n+1}=\frac{H{D}_i}{i+0.3}\left(i=2,H\le 2.8\right)\\ X_{n+1}=\frac{H{D}_i}{i+0.8}(i=1，H>2)\end{array}\hfill \end{array}$

4. 计算各步长的转移概率矩阵

4.1. 步长为1年的转移概率矩阵

取步长为1年，对1970~2020这51年进行统计总结状态变化次数如表3 (如果本次是1970年的状态，下次状态就是1971年的状态。表3中状态 $i,1\le i\le 5$ 所在的行与 $j,1\le j\le 5$ 所在的列交叉位置的数字记为 $f_{ij}$ ，表示本年是状态i而下一年是状态j的年份数量。如 $f_{11}=1,f_{25}=4,\cdots$ )。

Table 3. Statistical table of precipitation state of step size 1

表3. 步长为1年的降水量状态统计表

由表3可得步长为1年的转移概率矩阵

$P_1=\left[\begin{array}{ccccc}0.2& 0.6& 0& 0.2& 0\\ 0.1333& 0.1333& 0.0667& 0.4& 0.2667\\ 0.2857& 0.2857& 0& 0.4286& 0\\ 0& 0.3529& 0.2941& 0.2941& 0.0589\\ 0.1667& 0.3333& 0.1667& 0.1667& 0.1667\end{array}\right]$ 。

根据表3中的数据计算可得各种状态的边际概率(见表4)

Table 4. The marginal probability of each state

表4. 各种状态的边际概率

4.2. “马氏性”检验

由表3和表4及步长为1年的转移概率矩阵，计算得卡方分布 ${\chi }^2=27.9988$ 。给定显著性水平 $\alpha =0.05$ ，自由度 ${\left(m-1\right)}^2={\left(5-1\right)}^2=16$ 。查表可知 ${\chi }_{0.05}^2\left(16\right)=26.296$ ， ${\chi }^2>{\chi }_{0.05}^2\left(16\right)$ ，拒绝原假设，满足“马氏性”。

4.3. 各阶的自相关系数和权重

由3.3中公式得：(表5)

Table 5. Autocorrelation coefficient and weight of each order

表5. 各阶的自相关系数和权重

4.4. 步长为2~5年的转移概率矩阵

类似于步长为1年的情形，可得步长为2~5年的转移概率矩阵分别为：

$P_2=\left[\begin{array}{ccccc}0.25& 0.25& 0& 0.5& 0\\ 0.0667& 0.5333& 0.1333& 0.2& 0.0667\\ 0.1429& 0.4286& 0.1429& 0.2857& 0\\ 0.1765& 0.1176& 0.1765& 0.2353& 0.2941\\ 0& 0.1667& 0.1667& 0.6667& 0\end{array}\right]$ $P_3=\left[\begin{array}{ccccc}0& 0.25& 0.25& 0.25& 0.25\\ 0.1333& 0.2& 0.1333& 0.4667& 0.0667\\ 0.1667& 0.5& 0.1667& 0& 0.1667\\ 0.1765& 0.3529& 0.1176& 0.2353& 0.1176\\ 0& 0.3333& 0.1667& 0.5& 0\end{array}\right]$

$P_4=\left[\begin{array}{ccccc}0.25& 0.25& 0.25& 0.25& 0\\ 0.0667& 0.4& 0.1333& 0.3333& 0.0667\\ 0.1667& 0.3333& 0& 0.1667& 0.3333\\ 0.1875& 0.3125& 0.1875& 0.25& 0.0625\\ 0& 0.1667& 0.1667& 0.5& 0.1667\end{array}\right]$ $P_5=\left[\begin{array}{ccccc}0.25& 0& 0.25& 0.5& 0\\ 0.2143& 0.0714& 0.2143& 0.3571& 0.1429\\ 0& 0.0667& 0.1667& 0.1667& 0\\ 0.0625& 0.375& 0.125& 0.3125& 0.125\\ 0.1667& 0.6667& 0& 0& 0.1667\end{array}\right]$

5. 结果预测

5.1. 2021年降水量的预测结果

用2016~2020年的年降水量状态，对应的状态转移概率矩阵及表5对2021年的年降水量状态进行预测，结果见表6：

Table 6. Using 2016~2020 data to predict precipitation status in 2021

表6. 用2016~2020年的数据预测 2021年的降水量状态

由表6得0.2950最大，其对应的状态为2。故预测2021年降水量所对应的状态为2。而根据得到的数据中可以得知2021年的降水量为1478.5，对应的为状态2，预测成功。根据模糊集理论计算可得H = 2.69，年降水量1523.40。与实际降水量相对误差为3.0%。预测效果较好。

5.2. 2022年降水量的预测结果

同样的方法利用2017~2021年的年降水量数据，可以预测2022年降水量。

由表7可知，0.3130最大，因此最大加权预测2022年降水量所对应的状态为4。而根据得到的数据中可以得知2022年的降水量为1079.4，对应的为状态4，预测成功。根据模糊集理论计算可得H = 3.20，年降水量1014.97。与实际降水量相对误差为5.9%。预测效果较好。

Table 7. Using 2017~2021 data to predict precipitation status in 2022

表7. 用2017~2021年的数据预测2022年的降水量状态

5.3. 遍历性与平稳分布

以相依性最强的步长为2的马尔可夫链进行分析，由于降雨量的5个状态是互通的，没有周期性。因此这是一个不可约的正常返的马尔可夫链。此链具有遍历性，其极限分布即为为平稳分布。令状态为i，j的平稳分布分别为 ${\pi }_i,{\pi }_j$ ，各状态重现的周期为 $T_i$ ，则有方程组：

$\{\begin{array}{c}{\pi }_j\text{=}{\displaystyle \underset{i=1}{\overset{5}{\sum }}{\pi }_i{p}_{ij}}\\ {\displaystyle \underset{i=1}{\overset{5}{\sum }}{\pi }_j=1,{\pi }_j\ge 0}\end{array}$ 。

设各状态重现的周期为 $T_i\left(T_i=1/{\pi }_i\right)$ ，根据步长为2的状态转移矩阵，可得 ${\pi }_j$ 与 $T_j$ ，见表8：

Table 8. Limit distribution and recurrence period of each state

表8. 极限分布与各状态重现周期

由表8可知，按照本文的分级标准，在1970~2020年共51年的降水过程中，偏丰年出现的概率最大，平均每隔3.0371年出现一次；丰涝年出现的概率最小，平均每隔11.0375年才出现一次。

6. 模型拓展训练

将上面所建模型应用到不同城市(选择南方和北方共10个城市)，取上述城市最近2000~2020年的降水量数据，预测2021年的降水量。计算显示，在预测降水量状态方面，除了杭州、南京的错了以外，别的城市的降水量状态都预测准确。在用模糊集理论预测降雨量方面(见表9)，相对误差超过20%的占20%，相对误差在20%内的占80%，相对误差在10%内的占70%，总体预测效果较好。南方的城市相对来说误差比较高，有的甚至预测失败，这说明所建模型还有需要改进的地方。

Table 9. Predicted values for different cities

表9. 不同城市的预测值

7. 结语

加权马尔可夫链对降水量的预测北方地区的准确率要高于南方 [9] [10] 。上海地区属于南方，降水量大，本文根据上海的降雨量特点做了有针对性的分组，并在大量计算的基础上对模糊集理论的公式做了一些修改以使得预测结果更加准确。结果显示平均误差为4.5%，预测精度较高。将所建模型应用到不同城市，结果显示大部分城市预测结果较好，但有部分城市预测失败。这说明模型有需要提高的地方，这是后面继续研究要解决的问题。

基金项目

本文由上海工程技术大学大学生创新项目：螺纹钢期货价格影响因素分析(项目编号：CS2221003)资助。