一种降水时空预测方法及电子设备

1.本发明涉及降水时空预测领域，特别是涉及一种降水时空预测方法及电子设备。


背景技术：

2.降水的空间分布特征是影响径流模拟及一系列其他水文问题的主要控制因素。其中，径流过程模拟的不确定性主要来源于不确定性的降水输入，而洪水预测最敏感的因素就是降水的空间变化。然而，离散的雨量站点观测难以充分认识降水的时空分布特征，需要稳健的空间插值技术利用已知的数据来获取分布式的降水数据，支撑区域水文、水资源分析以及水资源的高效利用与管理。同时，通过插值预测的高分辨率空间降水场可以帮助识别发生洪水和山洪的局部高强度降水，有利于及时发布预警信息，并通过落实应急举措规避洪涝灾害的风险。此外，空间化的降水信息对于生态环境保护与治理也具有十分重要的意义。
3.目前，空间插值方法分类有很多种，可以概括为确定性和地理统计的插值方法，而每种插值方法都有相应的适用条件。确定性方法中最具代表性的是反距离权重法和泰森多边形法，这类模型方法对降水的空间变化通常采用均化处理方式，其在一定条件下对局部区域内降水插值具有良好的效果，但由于忽略了降水的随机性和空间相关性，难以客观反映降水的空间分布不均匀性。地理统计的插值方法则不仅考虑对象的空间相关性和依赖性，还考虑对象的统计特性。
4.但是，日和小时尺度的降水量存在“零值”较多的现象，使得小时和日降水的空间插值具有较大的困难，静态的插值方法难以反映降水的时序动态特征和依赖关系。


技术实现要素：

5.本技术提供一种降水时空预测方法及电子设备，既能有效表征降水时间动态变化，又能反映降水相关参数空间变异过程，获得可靠的降水预测结果。
6.第一方面，本发明实施例提供了一种降水时空预测方法，包括：
7.获取降水观测数据、获取降水空间预报因子、获取时空变异参数的先验分布；
8.根据所述降水观测数据和所述降水空间预报因子，建立广义线性回归模型；
9.构建高斯过程模型；
10.根据所述广义线性回归模型、所述高斯过程模型以及所述时空变异参数的先验分布，构建贝叶斯层次模型，确定空间点降水量在观测条件下的后验概率分布形式；
11.根据所述后验概率分布形式，确定空间点降水量预测结果。
12.在一些实施例中，所述方法还包括：处理降水观测数据，得到处理后的降水观测数据，具体包括：
13.使用box-cox降水转换方法处理所述降水观测数据，用最大似然法估计最优变换参数λ，得到处理后的降水观测数据，所述box-cox降水转换方法，包括：
[0014][0015]
其中，si为空间位置，t为时刻，λ为变换参数，z(si,t)为t时刻在空间位置si的降水观测数据，t∈[1,t],1≤i≤n，n为降水观测站点数目。
[0016]
在一些实施例中，所述方法还包括：筛选出具有显著相关性的降水空间预报因子，具体包括：
[0017]
通过假设检验的方法，计算所述降水空间预报因子与所述降水观测数据具有显著相关性的概率值；
[0018]
设定预设值，比较概率值和预设值的大小：
[0019]
若概率值大于预设值，则所述降水空间预报因子与所述降水观测数据不具有显著相关性；
[0020]
若概率值小于或等于预设值，则所述降水空间预报因子与所述降水观测数据具有显著相关性。
[0021]
在一些实施例中，所述获取时空变异参数的先验分布，包括：
[0022]
初始化时空变异参数，将所述时空变异参数设定在预设范围，其中，所述时空变异参数符合逆伽马分布；
[0023]
指定所述时空变异参数的逆伽马分布的超参数，得到所述时空变异参数的先验分布。
[0024]
在一些实施例中，所述根据所述降水观测数据和所述降水空间预报因子，建立广义线性回归模型，包括：
[0025][0026]yt
为所述降水观测数据矩阵的转置,o
t
为平均动态过程，ε
t
为白噪声，ε
t
～n(0，σ
ε2
)，服从以0为期望值，以σ
ε2
为标准差的正态分布，σ
ε2
为时空变异第一参数，α为截距，βj(si)＝(βj(s1)，βj(s2)，...，βj(sn))'为对应协变量预报因子的空间变化参数，x
jt
为对角矩阵，对角线上的元素代表某空间位置上降水空间预报因子xj的值，1≤i≤n，n为降水观测站点数目，m为所述筛选出的具有显著相关性的降水空间预报因子的数目，η
t
为空间相关的随机误差项。
[0027]
在一些实施例中，所述构建高斯过程模型，包括：
[0028]
所述广义线性回归模型的参数η
t
和βj(si)的变异服从高斯过程，获得以高斯过程表示的空间变化过程如下：
[0029][0030]
η
t
为空间相关的随机误差项，βj(si)为对应协变量预报因子的空间变化参数，1≤i≤n，n为降水观测站点数目，gp表示高斯过程，为时空变异第二参数，为时空变异第三参数，κ(
·
)公式如下其中，si和sj为空间中任意两个点的空间位置，||s
i-sj||为两个空间位置之间的距离；参数φ为时空变异第四参数，参数ν为时空变异第五参数，当距离变化时，参数φ为相关性的衰减程度，参数ν表示随机场的平滑程度；г(ν)为标准的伽马函数，k为ν阶的第二类贝塞尔函数。
[0031]
在一些实施例中，所述根据所述后验概率分布形式，确定空间点降水量预测结果，包括：
[0032]
通过马尔可夫链蒙特卡罗mcmc抽样方法对所述后验概率分布形式进行抽样，以确定空间点降水量预测结果，其中，所述马尔可夫链蒙特卡罗mcmc抽样方法，包括：
[0033]
抽样获得各空间位置降水预测的后验分布，计算所述后验分布的均值统计量，将所述均值统计量确定为所述空间点降水量预测结果。
[0034]
在一些实施例中，所述方法还包括：通过克里金法修正所述空间点降水量预测结果，具体包括：
[0035]
获取待预测降水的研究区域的空间信息；
[0036]
通过克里金法，计算在待预测降水的空间位置第一范围内的降水观测站点对所述待预测降水的空间位置的降水的贡献权重，并根据所述在待预测降水的空间位置第一范围内的降水观测站点是否发生降水判断所述待预测降水的空间位置降水发生的可能性，可能性取值为1和0，1代表降水发生，0代表降水不发生，进而获得所述空间点降水量预测结果的修正值，得到修正后的空间点降水量预测结果，其中，所述克里金法，包括：
[0037][0038]
其中，λi为在待预测降水的空间位置第一范围内的降水观测站点对所述待预测降水的空间位置s0的降水的贡献权重，p(si)为观测站点是否发生降水的概率，取值为0或1，z为在待预测降水的空间位置s0第一范围内的观测站点的数目，w为判断降水是否发生的阈值，0≤w≤1。
[0039]
在一些实施例中，所述方法还包括：校验降水空间预测模型的预测效果，所述降水空间适配模型包括所述广义线性回归模型、所述高斯过程模型、所述时空变异参数的先验分布以及所述贝叶斯层次模型，具体包括：
[0040]
根据所述修正后的空间点降水量预测结果，采用交叉验证方法，通过去除任一参考观测站点的降水数据，利用剩余站点的观测数据运行模拟过程；
[0041]
最后计算pmcc,所述pmcc能反映所述降水空间预测模型预测效果，所述pmcc越小，所述降水空间预测模型的预测效果越好，所述计算pmcc，包括：
[0042][0043]
其中，n为观测站点数目，r为时间序列，t为时刻，
[0044]
e(y
l
(si，t)
rep-y
l
(si，t))2为数学期望，var y
l
(si，t)
rep
为方差，y
l
(si，t)
rep
为t时刻第i个降水观测站点的降水量预测结果，y
l
(si，t)为t时刻第i个降水观测站点的降水观测数据，pmcc为所述降水空间预测模型的预测效果。
[0045]
第二方面，本发明实施例还提供了一种计算机设备，包括：
[0046]
至少一个处理器；以及，
[0047]
与至少一个处理器通信连接的存储器；其中，
[0048]
存储器存储有可被至少一个处理器执行的指令，指令被至少一个处理器执行，以使至少一个处理器能够执行第一方面的方法。
[0049]
第三方面，本发明实施例还提供了一种非易失性计算机可读存储介质，计算机可读存储介质存储有计算机可执行指令，当计算机可执行指令被处理器执行时，使处理器执行第一方面的方法。
[0050]
本发明实施例的有益效果在于：通过提供一种降水时空预测方法，所述方法，包括：获取降水观测数据、降水空间预报因子、时空变异参数的先验分布；根据所述降水观测数据和所述降水空间预报因子，建立广义线性回归模型；构建高斯过程模型；根据所述广义线性回归模型、所述高斯过程模型以及所述时空变异参数的先验分布，构建贝叶斯层次模型，确定空间点降水量在观测条件下的后验概率分布形式；根据所述后验概率分布形式，确定空间点降水量预测结果。
[0051]
本发明实施例通过引入所述降水空间预报因子，让所述广义线性回归模型能反映所述研究区域的空间变化过程，能够提高预测精度，有效表征降水时间动态变化，反映降水相关参数空间变异过程，获得可靠的降水预测结果。
附图说明
[0052]
一个或多个实施例通过与之对应的附图中的图片进行示例性说明，这些示例性说明并不构成对实施例的限定，附图中具有相同参考数字标号的元件表示为类似的元件，除非有特别申明，附图中的图不构成比例限制。
[0053]
图1是本发明实施例提供的一种降水时空预测方法的应用场景示意图；
[0054]
图2是本发明实施例提供的一种降水时空预测方法的流程图示意图；
[0055]
图3是图2中的步骤s20的细化流程图；
[0056]
图4是本发明实施例提供的一种处理降水观测数据的流程图；
[0057]
图5是本发明实施例提供的一种筛选降水空间预报因子的流程图；
[0058]
图6是本发明实施例提供的一种修正空间点降水量预测结果的流程图；
[0059]
图7是本发明实施例提供的一种降水时空预测方法的整体流程图；
[0060]
图8是本发明实施例提供的一种降水时空预测方法的处理设备的硬件结构示意图。
具体实施方式
[0061]
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0062]
本发明实施例提供的降水时空预测方法可以应用于如图1所示的应用场景，在图1所示的应用场景中，包括研究区域10，所述研究区域10多个降水观测站点11，在本发明实施例中，图1所示的研究区域可以是北京市，该北京市内可以具有8个降雨观测站点。
[0063]
其中，所述研究区域10只是示意性的，其中，所述研究区域指的是待预测降水的一定地域范围内的区域，例如一个区、一个城市、或者地域上接近的若干个城市组成的地理区域，类似长三角地区、华北地区等，或者指的是一个自治区、省等地域概念，所述研究区域的大小、形状不做限制；所述降水观测站点11可以是任何的，能获取当地降水数据的站点或者物体，其中降水观测站点11的面积大小不做限制，但是所述降水观测站点具有一定的的观测范围，所述降水观测站点的分布也是任意的。
[0064]
在本发明实施例中，在进行降水预测前，需要先获得研究区域空间位置信息，包括：以研究区域的空间范围作为边界，设定空间分辨率，分辨率越高，获得的空间分布越精确；应用gis工具提取对应空间分辨率的栅格中心点的空间位置，即地理坐标，获得空间位置信息，所述空间位置信息为所述空间位置的网格数据，每个所述空间位置都有地形和降雨信息，所述网格数据包括每个所述空间位置的地形和降雨信息。
[0065]
在本发明实施例中，请参阅图2,图2是本发明实施例提供的一种降水时空预测方法的流程示意图；
[0066]
如图2所示，该降水时空预测方法，包括：
[0067]
步骤s10：获取降水观测数据、获取降水空间预报因子、获取时空变异参数的先验分布；
[0068]
在本发明实施例中，所述降水观测数据可以通过降水观测站点获得，可以理解的是，所述降水观测站点获取的降水观测数据仅是在所述降水观测站点的观测范围内的降水观测数据。
[0069]
在本发明实施例中，所述降水空间预报因子包括协变量，协变量是对降水有影响，能反映空间变化的变量，如高程、坡度、坡向、站点距离海洋的距离、气象雷达的栅格数据等，协变量可以通过多种途径获得，如公开的遥感数据，地图数据等。
[0070]
在本发明实施例中，所述时空变异参数包括时空变异第一参数σ
ε2
、时空变异第二参数σ
β2
、时空变异第三参数σ
η2
、时空变异第四参数φ和时空变异第五参数v，所述时空变异第一参数σ
ε2
、时空变异第二参数σ
β2
、时空变异第三参数σ
η2
、时空变异第四参数φ和时空变异第五参数v均符合逆伽马分布，所述获取时空变异参数的先验分布指的是获取所述时空变异第一参数σ
ε2
、时空变异第二参数σ
β2
、时空变异第三参数σ
η2
、时空变异第四参数φ和时
空变异第五参数v的第一逆伽马函数形状、第二逆伽马函数形状、第三逆伽马函数形状、第四逆伽马函数形状、第五逆伽马分布形状。
[0071]
在本发明实施例中，步骤s10还包括获取所述研究区域的空间位置信息；
[0072]
所述获取研究区域的空间位置信息，包括：
[0073]
以研究区域的空间范围作为边界，设定空间分辨率，分辨率越高，获得的空间分布越精确，应用gis工具提取对应空间分辨率的栅格中心点的空间位置，获得空间位置信息。
[0074]
具体的，请再参阅图3，图3是图2中的步骤s10的细化流程图；
[0075]
如图3所述，该步骤s10获取降水观测数据、获取降水空间预报因子、获取时空变异参数的先验分布，包括：
[0076]
步骤s111：初始化时空变异参数，将所述时空变异参数设定在预设范围，其中，所述时空变异参数符合逆伽马分布；
[0077]
具体的，所述时空变异参数包括时空变异第一参数σ
ε2
、时空变异第二参数σ
β2
、时空变异第三参数σ
η2
、时空变异第四参数φ和时空变异第五参数v，将所述时空变异第一参数、时空变异第二参数、时空变异第三参数、时空变异第四参数、时空变异第五参数分别设定在第一范围、第二范围、第三范围、第四范围、第五范围；
[0078]
分别指定时空变异第一参数σ
ε2
、时空变异第二参数σ
β2
、时空变异第三参数σ
η2
、时空变异第四参数φ和时空变异第五参数v的逆伽马分布的第一超参数、第二超参数、第三超参数、第四超参数、第五超参数，
[0079]
在本发明实施例中，所述时空变异第一参数、时空变异第二参数、时空变异第三参数、时空变异第四参数可以分别设定在σ
ε2
∈[0.01,0.1]，σ
β2
∈[0.01,0.1]，σ
η2
∈[0.1,10]和φ∈[0.001,0.1]范围内，时空变异第五参数v的取值范围可以根据经验选择。
[0080]
步骤s112：指定所述时空变异参数的逆伽马分布的超参数，得到所述时空变异参数的先验分布。
[0081]
具体的，分别指定时空变异第一参数σ
ε2
、时空变异第二参数σ
β2
、时空变异第三参数σ
η2
、时空变异第四参数φ和时空变异第五参数v的逆伽马分布的第一超参数、第二超参数、第三超参数、第四超参数、第五超参数，得到所述时空变异第一参数、时空变异第二参数、时空变异第三参数、时空变异第四参数、时空变异第五参数的第一先验分布、第二先验分布、第三先验分布、第四先验分布、第五先验分布。
[0082]
在本发明实施例中，所述时空变异第一参数、时空变异第二参数、时空变异第三参数、时空变异第四参数、时空变异第五参数均符合逆伽马分布，其中逆伽马函数的形状由超参数决定，通过设定超参数，可以得到不同的逆伽马函数形状，根据所述不同的逆伽马函数形状，确定所述时空变异第一参数、时空变异第二参数、时空变异第三参数、时空变异第四参数、时空变异第五参数的第一先验分布、第二先验分布、第三先验分布、第四先验分布、第五先验分布，超参数的取值根据经验决定，可以根据实验结果调参或者参考公开的论文或报道进行取值。
[0083]
在本发明实施例中，时空变异第一参数σ
ε2
为白噪声方差，也为块金方差，反映的是最小抽样尺度以下变量的变异性及测量误差，时空变异第二参数σ
β2
为βj(si)的超参数，βj(si)为对应协变量预报因子的空间变化参数，时空变异第三参数σ
η2
为空间随机过程的位点不变方差，时空变异第四参数φ空间平滑参数，时空变异第五参数v为空间衰减参数。
[0084]
通过获取所述时空变异参数的先验分布，本发明能够利用所述时空变异参数的先验分布构建广义线性回归模型和高斯过程模型。
[0085]
在获取降水观测数据之后，所述方法还包括：对所述降水观测数据进行处理。
[0086]
具体的，请再参阅图4，图4是本发明实施例提供的一种处理降水观测数据的流程图；
[0087]
步骤s121：使用box-cox转换方法处理所述降水观测数据，具体包括：
[0088]
使用box-cox转换方法处理降水观测数据z(si,t)，降水分布具有明显的偏态，需要通过降水数据转换以得到更加可靠的小时和日降水数据的模拟结果，在本发明实施例中，所述处理方法为box-cox变换方法，变换公式如下：
[0089][0090]
其中，si为空间位置，t为时刻，λ为变换参数，z(si,t)为t时刻在空间位置si的降水观测数据，t∈[1,t],t为降水观测数据的时间，例如，如果所述降水观测数据为所述研究区域内的某一观测站点的近20天内的降水观测数据，则t为20，1≤i≤n，n为降水观测站点数目，可以理解的是，n和i为正整数。
[0091]
步骤s122：用最大似然法估计最优变换参数λ
*
，得到处理后的降水观测数据，具体包括：
[0092]
根据现有降水观测数据样本，求样本的联合概率密度函数最大时的变换参数的参数值；假定样本是对数正态分布，将所有样本带入概率密度函数，该概率密度函数的乘积是关于λ的一个函数形式，确定乘积最大时的该概率密度函数，结合该概率密度函数对λ求偏导，偏导数为0时代表该概率密度函数取得极值，对应λ就是最优变换参数λ
*
，最后根据最优变换参数λ
*
求得处理后的降水观测数据。
[0093]
降水分布具有明显的偏态，通过使用box-cox转换方法处理所述降水观测数据，本发明能得到更可靠的降水观测数据的模拟结果。
[0094]
在获取降水空间预报因子之后，所述方法还包括：对所述降水空间预报因子进行筛选；
[0095]
具体的，请再参阅图5，图5是本发明实施例提供的一种筛选降水空间预报因子的流程图；
[0096]
如图5所示，所述筛选降水空间预报因子，包括：
[0097]
步骤s131：通过假设检验的方法，计算所述降水空间预报因子与所述降水观测数据具有显著相关性的概率值，具体包括：
[0098]
将所述降水空间预报因子与所述降水观测数据进行假设检验，获得p值，p值即所述降水空间预报因子与所述降水观测数据具有显著相关性的概率值，p值表示所述协变量与所述处理后的降水观测数据的统计学差异，p值越小，所述协变量与所述处理后的降水观测数据具有越显著的统计学差异，后续处理得到的空间点降水量预测结果越精准。
[0099]
步骤s132：设定预设值，所述预设值用于与概率值比较大小，在本发明实施例中，
所述预设值可以为0.05。
[0100]
步骤s133：比较预设值与概率值的大小，所述比较预设值与概率值的大小可以确定所述降水空间预报因子与所述降水观测数据是否具有显著相关性。
[0101]
步骤s134：若概率值大于预设值，则所述降水空间预报因子与所述降水观测数据不具有显著相关性。
[0102]
步骤s135：若概率值小于或等于预设值，则所述降水空间预报因子与所述降水观测数据具有显著相关性。
[0103]
通过筛选降水空间预报因子，本发明能获得与所述降水观测数据具有显著相关性的降水空间预报因子，所述具有显著相关性的降水空间预报因子对所述研究区域的降水具有显著影响。
[0104]
步骤s20：根据所述降水观测数据和所述降水空间预报因子，建立广义线性回归模型；
[0105]
所述广义线性回归模型如下：
[0106][0107]yt
为所述降水观测数据矩阵的转置,o
t
为平均动态过程，ε
t
为白噪声，ε
t
～n(0，σ
ε2
)，服从以0为期望值，以σ
ε2
为标准差的正态分布，σ
ε2
为时空变异第一参数，α为截距，βj(si)＝(βj(s1)，βj(s2)，...，βj(sn))
′
为对应协变量预报因子的空间变化参数，x
jt
为对角矩阵，对角线上的元素代表某空间位置上降水空间预报因子xj的值，1≤i≤n，n为降水观测站点数目，m为所述筛选出的具有显著相关性的降水空间预报因子的数目，η
t
为空间相关的随机误差项。
[0108]
其中,其中,,o
t
＝(o(s1,t),o(s2,t),
…
,o(sn,t))’，η
t
＝(η(s1,t),η(s2,t),
…
,η(sn,t))’是服从高斯分布的空间相关的随机误差项，即η(s,t)～n(0,∑
η
)，服从以0为期望值，以∑
η
为标准差的正态分布，其中，∑
η
＝σ
η2
s，σ
η2
为时空变异第三参数，s＝κ(φ,ν)表示空间相关性的系数，符合马特恩协方差模型。
[0109]
其中，κ(
·
)公式如下：
[0110][0111]
其中，si和sj为空间中任意两个点的空间位置，||s
i-sj||为两个空间位置之间的距离；参数φ为时空变异第四参数，参数ν为时空变异第五参数，当距离变化时，参数φ为相关性的衰减程度，参数ν表示随机场的平滑程度；γ(ν)为标准的伽马函数，k为ν阶的第二类贝塞尔函数。
[0112]
通过将所述具有显著相关性的降水空间预报因子融入广义线性回归模型中，本发
明建立的广义线性回归模型能反映降水的时空分布特点。
[0113]
步骤s30：构建高斯过程模型；
[0114]
假设协变量的空间变化参数βj(si)，空间相关的随机误差项η
t
服从高斯过程，则我们可以构建高斯过程模型如下：
[0115][0116]
η
t
为空间相关的随机误差项，βj(si)为对应协变量预报因子的空间变化参数，1≤i≤n，n为降水观测站点数目，gp表示高斯过程，为时空变异第二参数，为时空变异第三参数，κ(
·
)公式与步骤s30中的相同。
[0117]
通过构建高斯过程模型，本发明可以描述空间降水潜在的高斯过程，处理降水时间序列中高频率的无降水观测数据的删失问题
[0118]
步骤s40：根据所述广义线性回归模型、所述高斯过程模型以及所述时空变异参数的先验分布，构建贝叶斯层次模型，确定空间点降水量在观测条件下的后验概率分布形式；
[0119]
具体的，根据所述广义线性回归模型、所述高斯过程模型以及所述时空变异参数的先验分布,确定数据模型：过程模型：0～gp(xj，βj(si)，σ
η2
κ(
·
))，参数模型：θ先验分布：
[0120]
其中，y服从以o为数学期望，以为标准差的正态分布，i为单位矩阵，为时空变异第一参数，o服从以x
jt
βj(si)为数学期望，以σ
η2
κ(
·
)为标准差的正态分布，x
jt
为对角矩阵，对角线上的元素代表某空间位置上降水空间预报因子xj的值，σ
η2
为时空变异第三参数，βj(si)＝(βj(s1)，βj(s2)，...，βj(sn))'为对应协变量预报因子的空间变化参数，θ包括空间时空变异第一参数σ
ε2
、时空变异第二参数σ
β2
、时空变异第三参数σ
η2
和时空变异第四参数φ和时空变异第五参数v，κ(
·
)公式如下：
[0121][0122]
基于所述数据模型、所述过程模型和所述时空变异参数的先验分布，构建贝叶斯层次模型，获得降水观测条件下的空间点降水量后验概率分布形式如下：
[0123][0124]
其中，n＝n*t，n为降水观测站点数目，t为降水观测数据时刻的最大值，y为降水观
测数据，而y*为站点之外空间位置没有观测的数据集，即要预测的降水数据集，同时，θ包括空间时空变异第一参数σ
ε2
、时空变异第二参数σ
β2
、时空变异第三参数σ
η2
和时空变异第四参数φ和时空变异第五参数v，x
jt
为对角矩阵，对角线上的元素代表某空间位置上降水空间预报因子xj的值，βj(si)＝(βj(s1)，βj(s2)，...，βj(sn))'为对应协变量预报因子的空间变化参数，1≤i≤n，n为降水观测站点数目；π(σ
η2
)、π(σ
ε2
)、π(σ
β2
)、π(φ)和π(v)表示时空变异参数的先验分布，βj(si)是回归系数，其先验分布同样服从伽马的逆函数分布，假设其超参数分别为a和b，对应的平均值和方差分别为a/b和a/b2，通过指定合适的超参数，可确定回归系数βj(si)的分布形式，s为n
×
n的向量，表示所有观测站点之间的协方差矩阵，第i*j个元素为κ(si,sj；φ.ν)，表示观测站点i和观测站点j之间的协方差，当距离变化时，参数φ为相关性的衰减程度，而参数ν为随机场的平滑程度。
[0125]
下面，进一步获得降水观测条件下的空间点降水量后验概率分布形式：
[0126]
对于空间任意一点s0，在时刻t预测的降水分布为:y(s0,t)～n(o(s0,t),σ
ε2
)，σ
ε2
为时空变异第一参数，o(s0,t)的序列如下：
[0127][0128]
根据所述y(s0,t)～n(o(s0,t),σ
ε2
)和所述o(s0,t)的序列，进一步获得降水观测条件下的空间点降水量后验概率分布形式如下：
[0129][0130]
其中，βj(s0)|θ～n(s
12
s-1
βj(s),σ
β2
(1-s
12
s-1s21
))表示βj(s0)服从以s
12
s-1
βj(s)为期望，以σ
β2
(1-s
12
s-1s21
)为标准差的正态分布并与参数θ有关，此外，预测的o(s0,t)条件分布形式为o(s0,t)|βj(s0),均值和方差的计算如下：
[0131][0132]
其中，s
12
表示预测点与所有观测站点的协方差矩阵，s
12
为1
×
n的向量，s
21
为向量s
12
的转置，s
12
＝s
21’，n为观测站点的数量，第i个元素为κ(si,s0；φ.ν)，si为第i个观测站点的空间位置(1≤i≤n)，s0为预测点的空间位置。
[0133]
通过构建贝叶斯层次模型，本发明可以获得空间点降水量在观测条件下所述后验概率分布形式，根据所述后验概率分布形式可以确定空间点降水量预测结果。
[0134]
进一步地，所述方法还包括，根据时空变异参数的先验分布，用mcmc抽样方法估计所述时空变异参数的后验分布，所述mcmc抽样方法具体包括：
[0135]
通过构建一个马尔可夫链，基于该马尔可夫链进行随机游走，产生样本的序列，获得平稳分布的样本，之后使用该平稳分布的样本进行近似的数值计算，mcmc抽样方法在序列接近目标分布之前需要进行多次迭代，需要设定预热期，在本发明实施例中，预热期可以设定为4000次，最大迭代次数可以设定为10000次。
[0136]
步骤s50：根据所述后验概率分布形式，确定空间点降水量预测结果，具体包括：
[0137]
通过借助马尔可夫链蒙特卡罗mcmc抽样的方法抽样估计各空间位置降水预测的后验分布π(z(s0,t)|z)，后验分布的均值统计量确定为空间点降水量预测结果所述mcmc抽样方法具体包括：
[0138]
通过构建一个马尔可夫链，基于该马尔可夫链进行随机游走，产生样本的序列，获得平稳分布的样本，之后使用该平稳分布的样本进行近似的数值计算，mcmc抽样方法在序列接近目标分布之前需要进行多次迭代，需要设定预热期，在本发明实施例中，可选的，预热期可以设定为4000次，最大迭代次数可以设定为10000次。
[0139]
请再参阅图6，图6是本发明实施例提供的一种修正空间点降水量预测结果的流程图；
[0140]
步骤s60：通过克里金法修正所述空间点降水量预测结果，具体包括：
[0141]
通过克里金法计算各降水观测站点对空间位置s0降水的贡献权重λi，根据所述各降水观测站点是否发生降水，判断空间位置s0处降水发生的概率值p(s0,t)(1：降水发生；0：降水不发生)，进而获得降水修正后的空间点降水量预测结果，其中，所述克里金插值方法，如下式所示：
[0142][0143]
其中，z为在待预测降水的空间位置s0第一范围内的观测站点的数目，p(si)为第i个降水观测站点的降水观测数据，0≤w≤1为判断降水是否发生的阈值，在本发明实施例中，w的推荐范围为[0.2,0.5]。
[0144]
最终的空间降水的预测值为所述修正后的空间点降水量预测结果是具有一定空间分辨率的降水分布时间序列。
[0145]
在本发明实施例中，所述计算权重λi的函数可以是半方差函数，并进一步的，可以根据所述降水空间预测结果的修正值，调整表示随机场的平滑程度的参数v的参数值。
[0146]
通过克里金法修正空间点降水量预测结果，本发明克服了在日和小时降水序列中，降水的“零值”较多的问题，完善了对降水不发生的情况的预测结果的合理性，实现了对降水预测值的修正。
[0147]
在本发明实施例中，通过克里金法修正空间点降水量预测结果前，还需要用box-cox逆变换转换空间点降水量预测结果，得到经过box-cox逆变换后的空间点降水量预测结果，再通过克里金法修正经过box-cox逆变换后的空间点降水量预测结果，得到最终的空间点降水量预测结果，其中，box-cox逆变换可以通过python、matlab等工具实现。
[0148]
步骤s70：校验降水空间预测模型的预测效果；
[0149]
在本发明实施例中，所述降水空间适配模型包括所述广义线性回归模型、所述高
斯过程模型、所述时空变异参数的先验分布以及所述贝叶斯层次模型；
[0150]
在本发明实施例中，所述校验方法可以是交叉验证方法，通过去除任一降水观测站点的降水观测数据，利用剩余降水观测站点的降水观测数据运行模拟过程，每次去除一个降水观测站点的降水观测数据，不重复去除，共n个降水观测站点，故共运行n次模拟过程，得到n个模拟结果；
[0151]
得到模拟结果以后，计算pmcc,所述pmcc能反映所述降水空间预测模型预测效果，所述pmcc公式如下：
[0152][0153]
其中，n为观测站点数目，r为时间序列，t为时刻，e(y
l
(si，t)
rep-y
l
(si，t))2为数学期望，var y
l
(si，t)
rep
为方差，y
l
(si，t)
rep
为t时刻第i个降水观测站点的降水量预测结果，y
l
(si，t)为t时刻第i个降水观测站点的降水观测数据，pmcc为所述降水空间预测模型的预测效果。
[0154]
在本发明实施例中，可以进一步获取均方根差、平均绝对误差降水误报率和确定系数等相关评价指标来校验模型预测效果。
[0155]
通过校验方法校验所述降水空间预测模型，本发明可以确定模型的预测效果，进一步的，可以根据模型预测效果修改模型参数值。
[0156]
请再参阅图7，图7是本发明实施例提供的一种降水时空预测方法的整体流程图；
[0157]
如图7所示，所述降水时空预测方法的整体流程，包括：
[0158]
步骤s811：获取降水观测数据。
[0159]
步骤s812：降水观测数据box-cox转换；
[0160]
其中，降水观测数据box-cox转换的目的是为了解决降水分布具有明显偏态的问题，通过降水数据转换以得到更加可靠的小时和日降水数据的模拟结果。
[0161]
步骤s813:获得并筛选具有相关性的降水空间预报因子；
[0162]
具体的，筛选具有显著相关性的降水空间预报因子，筛选方法为假设检验方法。
[0163]
步骤s814：构建广义线性回归模型和高斯过程模型；
[0164]
具体的，构建广义线性回归模型和高斯过程模型，将影响因子的作用融入降水空间预测中。
[0165]
步骤s815：初始化时空变异参数；
[0166]
具体的，初始化时空变异参数，包括：该时空变异参数的函数形式符合逆伽马分布，通过指定该时空变异参数的逆伽马分布函数的超参数，得到时空变异参数的先验分布。
[0167]
步骤s816：通过马尔可夫蒙特卡罗模拟，获得时空变异参数的后验分布；
[0168]
具体的，通过马尔可夫蒙特卡罗模拟，获得时空变异参数的后验分布，该后验分布用来得到空间点降水在观测条件下的后验分布。
[0169]
步骤s817：设定空间采样分辨率；
[0170]
具体的，设定空间采样率分辨率，包括：以研究区域的空间范围作为边界，设定空间分辨率，应用gis工具提取对应空间分辨率的栅格中心点的空间位置，即地理坐标。
[0171]
步骤s818：获得空间点位置信息；
[0172]
具体的，获得空间点位置信息，该空间点位置信息由gis工具提取对应空间分辨率的栅格中心点的空间位置得到，将生成的空间位置信息，输入降水空间预测适配模型中，计算获取研究区域内降水空间分布的时间序列图谱，进而可以对降水时空模式进行深入分析和研究。
[0173]
步骤s819：获得点降水o(s0,t)的预测序列；
[0174]
具体的，获得点降水o(s0,t)的预测序列，所述o(s0,t)的预测序列与空间点位置信息有关，步骤s40已经介绍，在此不再赘述。
[0175]
步骤s820：构建贝叶斯层次模型，包括：数据模型、过程模型、参数模型。
[0176]
步骤s821：获得降水空间预测适配模型；
[0177]
具体的，获得降水空间预测适配模型，该降水空间预测适配模型包括：广义线性回归模型、高斯过程模型、时空变异参数的先验分布、贝叶斯层次模型。
[0178]
步骤s822：获得空间点降水量在观测条件下的后验分布；
[0179]
具体的，获得空间点降水量在观测条件下的后验分布，该后验分布通过降水空间适配模型确定，运用了层次贝叶斯的基本原理，通过一层层的条件分布推求预测值。
[0180]
步骤s823：通过马尔可夫蒙特卡罗模拟，获得空间点降水量预测结果；
[0181]
具体的，通过马尔可夫蒙特卡罗模拟，获得空间点降水量预测结果，该步骤是用马尔可夫蒙特卡罗模拟方法估计空间点降水量。
[0182]
步骤s824：通过box-cox逆变换转换空间点降水量预测结果；
[0183]
具体的，box-cox逆变换可以通过python、matlab等工具实现。
[0184]
步骤s825：通过克里金法修正空间点降水量预测结果；
[0185]
具体的，通过克里金法修正经box-cox逆变换转换后的空间点降水量预测结果，处理得到更准确的预测结果。
[0186]
步骤s826：获得降水空间分布的时间序列图谱；
[0187]
具体的，获得降水空间分布的时间序列图谱，该图谱为具有一定空间分辨率的降水分布时间序列。
[0188]
步骤s827：模型校验和预报效果评估(pmcc、rmse、r2等)，校验方法包括交叉验证方法，评估方法包括计算pmcc。
[0189]
本发明实施例提供的降水时空预测方法，建立了时空变异参数模型描述空间降水潜在的高斯过程，处理降水时间序列中高频率的无降水观测数据的删失问题，引入广义线性回归模型以利于引入降水空间预报因子提高降水的预测精度，结合贝叶斯层次推理和蒙特卡洛-马尔可夫抽样算法求解复杂的高斯预测过程和概率计算，获取高精度的空间降水场的时间序列。同时，本发明提出的模型方法可进一步融入时间动态参数，并可以在模型层次结构中考虑动态参数的迟滞效应。
[0190]
请参阅图8，图8是本发明的电子设备示意图。
[0191]
如图8所示，该电子设备80包括：一个或多个处理器81以及存储器82，图8中以一个处理器81为例。
[0192]
处理器81和存储器82可以通过总线或者其他方式连接，图8中以通过总线连接为例。
[0193]
处理器81，用于获取降水观测数据、获取降水空间预报因子、获取时空变异参数的
先验分布；根据所述降水观测数据和所述降水空间预报因子，建立广义线性回归模型；构建高斯过程模型；根据所述广义线性回归模型、所述高斯过程模型以及所述时空变异参数的先验分布，构建贝叶斯层次模型，确定空间点降水量在观测条件下的后验概率分布形式；根据所述后验概率分布形式，确定空间点降水量预测结果。
[0194]
本发明实施例通过引入所述降水空间预报因子，让所述广义线性回归模型能反映所述研究区域的空间变化过程，能够提高预测精度，有效表征降水时间动态变化，反映降水相关参数空间变异过程，获得可靠的降水预测结果。
[0195]
存储器82作为一种非易失性计算机可读存储介质，可用于存储非易失性软件程序、非易失性计算机可执行程序以及模块，如本技术实施例中的降水时空预测方法的程序指令/模块。处理器82通过运行存储在存储器82中的非易失性软件程序、指令以及模块，从而执行控制器的各种功能应用以及数据处理，即实现上述方法实施例的降水时空预测方法。
[0196]
存储器82可以包括存储程序区和存储数据区，其中，存储程序区可存储操作系统、至少一个功能所需要的应用程序；存储数据区可存储根据控制器的使用所创建的数据等。此外，存储器82可以包括高速随机存取存储器，还可以包括非易失性存储器，例如至少一个磁盘存储器件、闪存器件、或其他非易失性固态存储器件。在一些实施例中，存储器82可选包括相对于处理器81远程设置的存储器，这些远程存储器可以通过网络连接至控制器。上述网络的实例包括但不限于互联网、企业内部网、局域网、移动通信网及其组合。
[0197]
所述一个或者多个模块存储在所述存储器82中，当被所述一个或者多个处理器81执行时，执行上述任意方法实施例中的降水时空预测方法，例如，执行以上描述的图2中的方法步骤s10至步骤s60。
[0198]
需要说明的是，上述产品可执行本技术实施例所提供的方法，具备执行方法相应的功能模块和有益效果。未在本实施例中详尽描述的技术细节，可参见本技术实施例所提供的方法。
[0199]
本技术实施例提供了一种非易失性计算机可读存储介质，所述计算机可读存储介质存储有计算机可执行指令，该计算机可执行指令被一个或多个处理器执行，例如图8中的一个处理器81，可使得上述一个或多个处理器可执行上述任意方法实施例中的降水时空预测方法，执行以上描述的图2中的方法步骤s10至步骤s60。
[0200]
通过以上的实施例的描述，本领域普通技术人员可以清楚地了解到各实施例可借助软件加通用硬件平台的方式来实现，当然也可以通过硬件。本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程是可以通过计算机程序来指令相关的硬件来完成，所述的程序可存储于一计算机可读取存储介质中，该程序在执行时，可包括如上述各方法的实施例的流程。其中，所述的存储介质可为磁碟、光盘、只读存储记忆体(read-only memory,rom)或随机存储记忆体(randomaccessmemory,ram)等。
[0201]
最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；在本发明的思路下，以上实施例或者不同实施例中的技术特征之间也可以进行组合，步骤可以以任意顺序实现，并存在如上所述的本发明的不同方面的许多其它变化，为了简明，它们没有在细节中提供；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技
术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的范围。