【集合(ensemble)数据同化】集合卡尔曼滤波(EnKF)、集合均方根滤波(EnSRF)、集合变换卡尔曼滤波(ETKF)和集合调整卡尔曼滤波(EAKF)的关系

本文回顾集合卡尔曼滤波和均方根滤波相关知识,介绍卡尔曼滤波误差协方差演化,分析矩阵均方根与扰动矩阵关系,阐述集合卡尔曼滤波及‘确定性’卡尔曼均方根滤波,还介绍集合均方根滤波的几种变体,指出EnSRF能消除抽样误差,某些情况估计更精确。

笔者按

        {~~~~~~~}        正所谓“温故而知新”。最近,我又开启对“集合xxx滤波”的回顾和学习。
        {~~~~~~~}        其实关于集合数据同化我是有这样一个心路历程:
        {~~~~~~~}        刚开始接触集合同化,便将一切“集合xxx滤波”归为“EnKF”,觉得理所当然也无可厚非;
        {~~~~~~~}        后来随着对“EAKF”的研究,知晓了EnKF和EAKF的区别在于“观测是否扰动”以及“随机性和确定性”,但我始终声称“EAKF是一种确定性的EnKF”并乐此不疲;
        {~~~~~~~}        直到最近被审稿人指出“EAKF是一种确定性的EnKF”的提法令他困惑,我才终于意识到是时候、有必要好好掰扯掰扯这些“集合xxx滤波”了。
        {~~~~~~~}        于是在导师的帮助下,参考Tippett等(2003),在这篇分享中对卡尔曼滤波(Kalman filter)、均方根滤波(square root filter,SRF)、集合卡尔曼滤波(ensemble Kalman filter,EnKF)、集合均方根滤波(ensemble square root filter,EnSRF)、集合变换卡尔曼滤波(ensemble transform Kalman filter,ETKF)和集合调整卡尔曼滤波(ensemble adjustment Kalman filter,EAKF)的公式推导以及它们的关系作以整理。
        {~~~~~~~}        这次的分享可能对初学者不太友好,要求读者对于数据同化、集合滤波、线性代数等数理知识有较好的基础。其实每过一段时间温习一遍类似的数学公式推导(如数理方程、数值模式、同化算法等等)对于沉淀知识净化心灵很有效果(哈哈)。后面有时间我会更新一些入门级的分享。这也是我第一次发布学习笔记,有不恰当之处我会及时改正,请读者多多包涵。欢迎感兴趣的同行一起交流!

1 前言

        {~~~~~~~}        集合数据同化方法利用状态估计方法以及预测和分析误差协方差矩阵的低秩形式来同化观测。这类方法的一个关键要素是采用适当的统计手段将预测集合转换为分析集合。这种转换可以通过将观测视为随机变量来“随机性”实现(Houtekamer和Mitchell,1998Burgers等,1998Evensen,2003),这类方法完善了Evensen(1994)最初提出的集合卡尔曼滤波(EnKF);也可以要求更新后的分析扰动满足卡尔曼滤波分析误差协方差的方程来“确定性”实现。“随机性”更新分析集合的EnKF方法和后续发展的“确定性”更新分析集合的方法都可被视作卡尔曼均方根滤波(Kalman SRF)的范畴Bierman,1977Maybeck,1982Heemink等,2001Tippett等,2003)。第2部分依次介绍了卡尔曼滤波、SRF框架下的“随机性”EnKF,并提出“确定性”Kalman SRF的概念。第3部分介绍了四种不同的“确定性”Kalman SRF(即EnSRF)变体。第4部分是简单的总结。

2 卡尔曼均方根滤波(Kalman SRF)

2.1 卡尔曼滤波误差协方差的演化

        {~~~~~~~}        为方便后续的表达,首先给出一些常用符号的含义。
        {~~~~~~~}         n n n,状态变量维度;
        {~~~~~~~}         p p p,观测个数;
        {~~~~~~~}         P k f {\bf{P}}_k^f Pkf,预报误差协方差矩阵,维度为 n × n n \times n n×n
        {~~~~~~~}         P k a {\bf{P}}_k^a Pka,分析误差协方差矩阵,维度为 n × n n \times n n×n
        {~~~~~~~}         M k {\bf{M}}_k^{} Mk,切线性动力模型;
        {~~~~~~~}         Q k { {\bf{Q}}_k} Qk,模型误差协方差矩阵,维度为 n × n n \times n n×n
        {~~~~~~~}         H k {\bf{H}}_k^{} Hk,观测投影算子,维度为 p × n p \times n p×n
        {~~~~~~~}         R k {\bf{R}}_k^{} Rk,观测误差协方差矩阵,维度为 p × p p \times p p×p
        {~~~~~~~}         K k {\bf{K}}_k^{} Kk,卡尔曼增益矩阵,维度为 n × p n \times p n×p
        {~~~~~~~}        卡尔曼滤波的核心是预报误差协方差和分析误差协方差的演化,它们的演化方程以及卡尔曼增益矩阵的定义如下:
P k f = M k P k − 1 a M k T + Q k          ( 1 ) {\bf{P}}_k^f = {\bf{M}}_k^{}{\bf{P}}_{k - 1}^a{\bf{M}}_k^{\rm{T}} + { {\bf{Q}}_k}~~~~~~~~(1) Pkf=MkPk1aMkT+Qk        (1) P k a = ( I − K k H k ) P k f          ( 2 ) {\bf{P}}_k^a = \left( { {\bf{I}} - {\bf{K}}_k^{}{\bf{H}}_k^{}} \right){\bf{P}}_k^f~~~~~~~~(2) Pka=(IKkHk)Pkf        (2) K k ≡ P k f H k T ( H k P k f H k T + R k ) − 1          ( 3 ) {\bf{K}}_k^{} \equiv {\bf{P}}_k^f{\bf{H}}_k^{\rm{T}}{\left( { {\bf{H}}_k^{}{\bf{P}}_k^f{\bf{H}}_k^{\rm{T}} + {\bf{R}}_k^{}} \right)^{ - 1}}~~~~~~~~(3) KkPkfHkT(HkPkfHkT+Rk)1        (3)其中, k k k表示当前分析步。可以看出,在以模型误差作为强迫的动力作用下,分析误差的传播由式(1)给出;式(2)表示在最优的数据同化方案中同化观测得到的分析误差协方差总是小于预报误差协方差。

2.2 矩阵均方根及其与扰动矩阵的关系

        {~~~~~~~}        为了在卡尔曼SRF的框架下分析EnKF和EnSRF,有必要清楚矩阵均方根是什么及其与扰动矩阵(扰动向量集合)的关系。对于一个维度为 n × n n \times n n×n的协方差矩阵 (对称正定矩阵),如果它的秩为 m m m,那么总有一个维度为 n × m n \times m n×m的矩阵 P = Z Z T {\bf{P}} = {\bf{Z}}{ {\bf{Z}}^{\rm{T}}} P=ZZT这里允许 m ≪ n m \ll n mn。不难知道, P k f {\bf{P}}_k^f Pkf P k a {\bf{P}}_k^a Pka都是对称正定矩阵,因此可以分解为 P k f = Z k f Z k f T {\bf{P}}_k^f = {\bf{Z}}_k^f{\bf{Z}}_k^{f{\rm{T}}} Pkf=ZkfZkfT P k a = Z k a Z k a T {\bf{P}}_k^a = {\bf{Z}}_k^a{\bf{Z}}_k^{a{\rm{T}}} Pka=ZkaZkaT其中 Z k f {\bf{Z}}_k^f Zkf Z k a {\bf{Z}}_k^a Zka分别是 P k f {\bf{P}}_k^f Pkf P k a {\bf{P}}_k^a Pka的矩阵均方根。根据集合卡尔曼滤波(EnKF)的定义, P k a {\bf{P}}_k^a Pka可以利用维度为 n × m n \times m n×m的扰动矩阵(扰动向量集合) S k a {\bf{S}}_k^a Ska计算,即 P k a = ( m − 1 ) − 1 S k a S k a T {\bf{P}}_k^a = {\left( {m - 1} \right)^{ - 1}}{\bf{S}}_k^a{\bf{S}}_k^{a{\rm{T}}} Pka=(m1)1SkaSkaT这里 m m m表示集合数。于是,容易得出 Z k a {\bf{Z}}_k^a Zka S k a {\bf{S}}_k^a Ska的关系是 Z k a = ( m − 1 ) − 1 2 S k a {\bf{Z}}_k^a = {\left( {m - 1} \right)^{ - \frac{1}{2}}}{\bf{S}}_k^a Zka=(m1)21Ska二者仅相差一个系数。因此,个人理解,无妨认为 Z k a {\bf{Z}}_k^a Zka是广义的分析向量的集合;同理, Z k f {\bf{Z}}_k^f Zkf是广义的预报集合。

2.3 集合卡尔曼滤波(EnKF)

        {~~~~~~~}        下面利用 Z k f {\bf{Z}}_k^f Zkf Z k a {\bf{Z}}_k^a Zka重新计算式(1)和式(2)。如果忽略 Q k { {\bf{Q}}_k} Qk,则有 Z k f Z k f T = M k Z k − 1 a Z k − 1 a T M k T = M k Z k − 1 a ( M k Z k − 1 a ) T {\bf{Z}}_k^f{\bf{Z}}_k^{f{\rm{T}}} = {\bf{M}}_k^{}{\bf{Z}}_{k - 1}^a{\bf{Z}}_{k - 1}^{a{\rm{T}}}{\bf{M}}_k^{\rm{T}} = {\bf{M}}_k^{}{\bf{Z}}_{k - 1}^a{\left( { {\bf{M}}_k^{}{\bf{Z}}_{k - 1}^a} \right)^{\rm{T}}} ZkfZkfT=MkZk1aZk1aTMkT=MkZk1a(MkZk1a)T于是容易得到预报集合的演化如下 Z k f = M k Z k − 1 a          ( 4 ) {\bf{Z}}_k^f = {\bf{M}}_k^{}{\bf{Z}}_{k - 1}^a~~~~~~~~(4) Zkf=MkZk1a        (4)从集合同化的角度看,式(4)意味着 Z k f {\bf{Z}}_k^f Zkf的每个列向量是 M k {\bf{M}}_k^{} Mk作用到 Z k − 1 a {\bf{Z}}_{k-1}^a Zk1a的每个列向量的结果,因此无妨认为 Z k f {\bf{Z}}_k^f Zkf Z k a {\bf{Z}}_k^a Zka的不同列就代表着不同的集合成员。事实上,模式误差的影响是不可以忽略的。如果有维度为 n × q n \times q n×q Z k d {\bf{Z}}_k^d Zkd

评论 6
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值