集合卡尔曼滤波在河道洪水预报中的应用

1975080514 1981090600 1983100413 1989070910
1975-8-5 14:00 1981-9-6 0:00 1983-10-4 13:00 1989-7-9 10:00
1975-8-14 8:00 1981-9-14 10:00 1983-10-13 7:00 1989-7-17 23:00
4
表3
洪号 开始时间
白马寺 4 场检验洪水的 2t 时段模拟洪水预报特征值
结束时间 预报洪峰 3 流量(m ) 2093 1151 1919 1413 实测洪 峰 流量 3 (m ) 1920 1110 1880 1320 峰值 绝对误 差 3 (m ) 173 41 39 93 平均 绝对误 差 3 (m ) 18.11 9.74 19.98 15.50 确定性 系数 0.9888 0.9916 0.9937 0.9889
实测 预报
0
50
100 150 时 段 T（ 小 时 ） 洪 号 1983100413实 测 与 预 报 流 量
200
250
2000
流 量 Q(立 方 米 )
1500 1000 500 0 0 50 100 150 时 段 T（ 小 时 ） 洪 号 1989070910实 测 与 预 报 流 量 200
n n 1 X bt MX at
（1）
式中：X 为状态变量；M 为预报模式。 分析部分：在预报部分得到模型预测集合，用集合平均代替真实值，同化观测资料求得分析集 合平均解。河道洪水流量是一个随机动态过程，测量流量值是含误差的不定值。集合 Kalman 滤波参
1
数估计方法通过对观测值施加观测噪声，可以模拟现场观测值的不确定性。对于集合 Kalman 滤波而 言，加入一定的噪声还可以提高滤波器的稳定性，保证估计量的标准差不会急剧下降[6]。采用 Monte Carlo 模拟高斯白噪声（均值为 0，标准差为观测值乘以噪声比例系数） ，流量观测向量变为一组经 过扰动的观测集合，其第 i 个样本可表示为：
预报部分
分析部分 [1] 预报误差的分析
n n B n H T ( X bt X bn )( HX bt HX bn )T
变量的集合预报
n n 1 X bt MX at
n n HB n H T HX bt HX bn HX bt HX bn



T
[2] 计算增益矩阵
1975080514 1981090600 1983100413 1989070910
1975-8-5 14:00 1981-9-6 0:00 1983-10-4 13:00 1989-7-9 10:00
1975-8-14 8:00 1981-9-14 10:00 1983-10-13 7:00 1989-7-17 23:00
洪号 开始时间
白马寺 4 场检验洪水的模拟洪水预报特征值
结束时间 预报洪峰 3 流量(m ) 1539 1034 1892 1224 实测洪 峰 流量 3 (m ) 1920 1110 1880 1320 峰值 绝对误 差 3 (m ) 381 76 12 96 平均 绝对误 差 3 (m ) 76.2 36 61.81 37.35 确定性 系数 0.9327 0.9550 0.9747 0.9772
30Байду номын сангаас0
洪 号 1975080514实 测 与 预 报 流 量 实测 预报t 预 报 2t
流 量 Q(立 方 米 )
2000 1000 0
0
50
1500
100 150 时 段 T（ 小 时 ） 洪 号 1981090600实 测 与 预 报 流 量
200
250
流 量 Q(立 方 米 )
1000 500 0
X i' HX i i
式中：H 为观测算子；为观测误差。 根据实测和预测值，结合预测误差矩阵，就可以直接更新模型参数：
n n n X at X bt K ( X tn HX bt )
（2）
（3） （4） （5）
K B n H T ( HB n H T O n ) l
2 集合 Kalman 滤波同化技术
Kalman 滤波不假设背景场的误差是定常的，而是利用协方差矩阵预报方程来计算[4]。其原理是 先进的，但是在实际应用中却因为计算量和存储量太大，在可预见的将来还无法应用到像最优插值 方法现在应用到的那种规模的数值模式和观测系统。Kalman 光滑是 Kalman 滤波的一种完善，它可 以利用后一时刻的资料来改进前一时刻的估计。Kalman 光滑分为固定区间光滑和固定时滞光滑两 种，其计算量比 Kalman 滤波更大。目前都在研究 Kalman 滤波的简化形式，比如基于粗空间分辨率 的模式来构造协方差矩阵预报方程。集合 Kalman 滤波(Evensen, 1994)是针对 Kalman 滤波中的协方 差矩阵预报模式有时会出现计算不稳定等问题而提出的，其主要思想是抛弃协方差矩阵预报模式， 而是直接利用 Monte-Carol 方法多次积分的模式， 从而得到背景场的误差协方差矩阵， 这个方法的计 [5] 算存储量比 Kalman 滤波小很多 。 集合 Kalman(EnKF)是一种纯的统计 Monte-Carlo 罗方法，模式状态的集合在状态空间中随时间 演变，集合的均值为模式预报值的最优估计，集合的发散分布代表了误差方差。在测量时刻，每一 个测量值被另一个集合所取代，集合的期望值就是测量值的最优估计，集合的方差反映了测量误差。 EnKF 同化方法包含预报部分和分析部分[6-8]。 预报部分：对于一个水文计算模型，模式将分析阶段得到的分析集合作为初始场向前预报至下 一个观测资料时刻（即数据同化时刻 t） ：
实测 预报
250
1500
流 量 Q(立 方 米 )
1000 500 0
实测 预报
0
50
100 150 时 段 T（ 小 时 ）
200
250
图 3 1975 年-1989 年共 4 场检验洪水实测与预报流量过程对比
计算结果显示，应用集合 Kalman 模型预报结果优于最小二乘法的预报结果，试验效果较好。 为了进一步研究延长预见期后的情况，下面应用集合 Kalman 滤模型在逐时段t 基础上进行 2t 时段模拟洪水预报。计算预报结果见表 3 、图 4 所示。
基于集合 Kalman 滤波的河道洪水预报研究
岳延兵 1,2，李致家 1，李振兴 2
(1.河海大学水文水资源学院、江苏 南京 210098；2.山西水利职业技术学院、山西 运城 044004)
摘
要：基于集合 Kalman 滤波的数据同化技术应用逐渐成熟，研究了集合 Kalman 滤波法与最小二乘法对马斯京根
1 河道洪水预报研究简述
河道洪水预报，是依据洪水波在河道中的运动、传播规律为基础，由河道上断面的水位、流量 值预报下游断面水位、流量值的方法[1-2]。常见的河道洪水预报方法有两类：相应水位（流量）法和 洪水演算法。按照描述洪水传播规律形式的不同，洪水演算法可分为水文学方法、水力学方法和基 于系统理论的方法。基于系统理论的方法直接描述河道系统输入与输出之间的因果关系，不依赖于 河道地形资料。基于系统理论的河道洪水演算方法具体可分为线性系统模型、时间序列模型、神经 网络模型等[2]。目前利用小波分析、遗传算法、神经网络与混沌分析等进行信息融合的预测方法也 逐渐发展起来[3]。
n n B n H T ( X bt X bn )( HX bt HX bn )T
n n HB n H T HX bt HX bn HX bt HX bn



T
（6）
式中：K 为卡尔曼增益矩阵；O 为观测误差协方差矩阵；B 为背景场误差协方差矩阵。 集合 Kalman 滤波进行数据同化的计算流程如图 1 所示。
模型参数率定的河道洪水预报技术，利用长水——白马寺实测洪水进行了对比检验，讨论了集合 Kalman 滤波的多时 段预报洪水过程及其特点。试验结果表明，应用集合 Kalman 滤波技术优于最小二乘法的预报效果，可有效延长洪水 预报的时段和提高预报精度。 关键词：集合 Kalman 滤波；数据同化；信息融合；水文预测
2
Wt K X t I t 1 X t Ot Ot 1 C 0 I t 1 C1 I t C 2 Ot
（7） （8）
式中：I 为上游输入流量；O 为下游输出流量；K 为蓄量流量关系曲线的坡度；X 为流量比重系 数；C 为流量演算方程系数。 利用长水——白马寺实测洪水数据，应用最小二乘法对马斯京根模型参数率定结果：X=0.35， K=11.86 ，t=12。计算得：C0=0.134875，C1=0.740462，C2=0.124663。 对白马寺站进行了 1975 年-1989 年共 4 场模拟洪水预报，结果见表 1、图 2 所示。 表1
2000
流 量 Q(立 方 米 )
洪 号 1975080514实 测 与 预 报 流 量
流 量 Q(立 方 米 )
2000
实测 预报
洪 号 1983100413实 测 与 预 报 流 量 实测 预报
1500 1000 500 0 0 50 100 150 时 段 T（ 小 时 ） 洪 号 1981090600实 测 与 预 报 流 量 200
3
2000
洪 号 1975080514实 测 与 预 报 流 量 实测 预报
流 量 Q(立 方 米 )
1500 1000 500 0 0 50 100 150 时 段 T（ 小 时 ） 洪 号 1981090600实 测 与 预 报 流 量 200
250
1500
流 量 Q(立 方 米 )
1000 500 0
200
250
0
50
100 150 时 段 T（ 小 时 ）
200
250
图2
1975 年~1989 年共 4 场检验洪水实测与预报流量过程对比
4 集合 Kalman 滤波同化洪水预报模型与检验
对白马寺站进行了 1975 年-1989 年共 4 场逐时段预见期的模拟洪水预报，应用集合 Kalman 滤 波模型预报结果见表 2 、图 3 所示。 表2
K B n H T ( HB n H T O n ) l
[3] 对变量的分析
n n n X at X bt K ( X tn HX bt )
图 1 集合 Kalman 滤波的计算流程
3 最小二乘法率定马斯京根模型参数与检验
为了便于对比研究，首先利用最小二乘法率定马斯京根模型参数并检验预报效果。马斯京根模 型的槽蓄方程式和流量演算方程式如下[9]：
洪号 开始时间
白马寺 4 场检验洪水的逐时段模拟洪水预报特征值
结束时间 预报洪峰 3 流量(m ) 1977 1117 1913 1358 实测洪 峰 3 流量(m ) 1920 1110 1880 1320 峰值 绝对误 差 3 (m ) 57 7 33 38 平均 绝对误 差 3 (m ) 7.54 3.73 8.41 4.95 确定性 系数 0.9971 0.9982 0.9983 0.9985
实测 预报t 预 报 2t
0
50
100 150 时 段 T（ 小 时 ）
200
250
2000
洪 号 1983100413实 测 与 预 报 流 量 实测 预报t 预 报 2t
流 量 Q(立 方 米 )
1500 1000 500 0 0 50 100 150 时 段 T（ 小 时 ） 洪 号 1989070910实 测 与 预 报 流 量 200
1975080514 1981090600 1983100413 1989070910
1975-8-5 14:00 1981-9-6 0:00 1983-10-4 13:00 1989-7-9 10:00
1975-8-14 8:00 1981-9-14 10:00 1983-10-13 7:00 1989-7-17 23:00
1500 1000 500 0 0 50 100 150 时 段 T（ 小 时 ） 洪 号 1989070910实 测 与 预 报 流 量 200
250
250
1500
流 量 Q(立 方 米 )
1500
流 量 Q(立 方 米 )
1000 500 0
实测 预报
1000 500 0
实测 预报
0
50
100 150 时 段 T（ 小 时 ）