On the Efficiency of Covariance Localisation of the Ensemble Kalman Filter Using Augmented Ensembles
- 1School of Earth Sciences and Centre of Excellence for Climate Extremes, University of Melbourne, Parkville, VIC, Australia
- 2NOAA Earth System Research Lab, Boulder, CO, United States
- 3School of Atmospheric Sciences, Nanjing University, Nanjing, China
A Commentary on
On the Efficiency of Covariance Localisation of the Ensemble Kalman Filter Using Augmented Ensembles
by Farchi, A., and Bocquet, M. (2019). On the Efficiency of Covariance Localisation of the Ensemble Kalman Filter Using Augmented Ensembles. Front. Appl. Math. Stat. 5:3. doi: 10.3389/fams.2019.00003
Introduction
In discussing Equation (39) (Equation (25) of Bocquet [1]), Farchi and Bouquet [2] state that “This perturbation update has been rediscovered by Bishop et al. [3] and included in their gain ETKF (GETKF) algorithm. However, the update formula used in the GETKF is prone to numerical cancellation errors as opposed to Equation (39)”. Here, we note:
(i) The predecessor of the GETKF eigenvalue form of the modified gain matrix equation appeared in Posselt and Bishop [4, 5]—before Bocquet [1].
(ii) The spectral shift theorem reduces the differences in the numerical cancellation errors referred to by Farchi and Bouquet.
(iii) The eigenvalue form enables Wang et al.'s [6] corrections for ensemble rank deficiency.
(iv) A proof of the equivalence of the eigenvalue form and Bouquet's form.
On page 12, Farchi and Bouquet [2] also state that “Such an extension had been discussed by Bishop et al. [3] but without numerical illustration.” This is incorrect. Lei et al. [7] used the GETKF to show that model space ensemble covariance localization provided satellite data assimilation (DA) performance comparable to 3DEnsVar.
To be specific about the forms of the modified gain matrix, let
where K is the total number of ensemble members in the ensemble forecast and where the n-vector is the ith member of the prior ensemble forecast and where the p-vector is the ith member of the prior ensemble forecast of the p-vector y of p observations. When p<K, the numerical cost of the pxp eigen decomposition
is less than the K × K eigen decomposition
In (2), E is a pxp eigenvector matrix for which , Γpxp is a pxp diagonal matrix of eigenvalues. In (3), C is a K × K orthonormal matrix of eigenvectors and ΓKxK is a K × K diagonal matrix of eigenvalues. At least K-p of the eigenvalues in ΓK × K will be equal to zero in the case of K>p. Equation's (2) and (3) are directly connected to the verbose singular value decomposition where where 0px(K−p) is a px(K−p) matrix of zeros. However, since the columns of C associated with zero eigenvalues cannot contribute to products of the matrix with other vectors, it is more efficient to work with the concise svd given by where lists the p columns of having non-zero eigenvalues. Posselt and Bishop [4, 5] note that is given by
and hence can be computed without performing an eigen decomposition of the larger K × K matrix in (3). Posselt and Bishop [4, 5] prove that for a linear observation operator , if
(see [5], Equation A10) then
The analysis perturbations are given by , hence,
is the perturbation update equation implied by Posselt and Bishop [4, 5].
In the above notation, and when propagation of small amplitude ensemble perturbations by the non-linear model is replaced by the propagation of raw ensemble perturbations by the non-linear model (i.e., no tangent linear model approximation is made), Bocquet's Equation (25) [1] for the ensemble perturbation update takes the form,
A fundamental difference between (7) and (8) is that while Bouquet multiplies Zf by the KxK matrix Posselt and Bishop multiply it by the Kxp matrix . When K>p, Posselt and Bishop's form only requires the eigenvector decomposition of a pxp matrix, whereas Bouquet's form requires the inversion of a larger K × K matrix. However, when p>K, the eigen decomposition (3) is cheaper than (2), becomes identical to the K × K matrix and becomes the concise svd of . In this case, EpxK is efficiently given by and (7) becomes
Dividing Equation (9) by recovers Equation (24) of Bishop et al. [3].
The above shows that Bishop et al.'s Equation (24) [3] was not “rediscovered” from Bocquet's [1] form as implied by Farchi and Bocquet [2]. It is an extension of Posselt and Bishop's [4, 5] eigenvalue form to the case of K>p. Equation (9) is just an eigenvalue form of the modified gain matrix of Whitaker and Hamill's [8] Ensemble Square Root Filter.
Equivalence of (9) and (8)
Bocquet's Equation (25) [1] can be derived from (9) with the following steps:
(i) Drop the dimension subscripts and manipulate [I − (Γ + I)−1/2]Γ−1 as follows
(ii) Use (10) in (9) to give
(iii) But
(iv) Using (12) and (3) in (11) gives
Equation (13) is equivalent to (8) and Bocquet's Equation (25) [1].
Numerical Issues, Condition Numbers, and Understanding
Numerical cancellation errors increase when the condition number of the matrix increases. Let us define the scalars and to, respectively, denote the maximum and minimum of the eigenvalues listed in the eigenvalue matrix Γpxp. The condition number of is . Because can be zero, can be infinite. In contrast, is bounded above by . However, note that the matrix has the eigenvalue decomposition
and hence has which is bounded above by where α is a positive scalar. Hence, α can be chosen to create a matrix that is better conditioned than . Once the eigen decomposition CΛCT of (14) has been obtained, one obtains the eigenvalues required by the GETKF or ETKF using Γ = Λ − αI. Thus, condition number differences between the Bouquet and eigenvalue form are easily eliminated.
The eigenvalue form lends understanding to the performance of DA schemes in much the same way that Empirical Orthogonal Functions lend understanding to climate variability. Wang et al. [6] used this understanding to correct gross aspects of the eigenvalue overestimation that occurs when the size of the ensemble is much smaller than the rank of the true observation space forecast error covariance matrix.
Discussion
Bocquet [1] and Farchi and Bocquet [2] may have overlooked Posselt and Bishop's [4, 5] work because:
(i) It is difficult to find all relevant literature to one's own work.
(ii) Bishop et al. [3] did not cite Posselt and Bishop [4, 5].
(iii) The equivalence of Posselt and Bishop's [4, 5] form and Bocquet's [1] form is not obvious.
Similarly, Bishop et al. [3] overlooked Bocquet's [1] work because of (i) and (iii). This note serves to clarify the origins and uses of modified gain matrices used in ensemble DA.
Author Contributions
CB led the study and wrote the text. JW and LL carefully reviewed the text.
Funding
CB acknowledges support from the Australian Research Council's Centre of Excellence in Climate Extremes (CE170100023). LL acknowledges joint sponsorship by the National Key R&D Program of China through grant 2017YFC1501603 and the National Natural Science Foundation of China through grant 41675052. JW acknowledges support through the NOAA High-Impact Weather Prediction Project (HIWPP) under award NA14OAR4830123, and the NOAA/NWS Next-Generation Global Prediction System (NGGPS) project.
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
References
1. Bocquet M. Localization and the iterative ensemble Kalman smoother. Q J R Meteorol Soc. (2016) 142:1075–89. doi: 10.1002/qj.2711
2. Farchi A, Bocquet M. On the efficiency of covariance localisation of the ensemble Kalman filter using augmented ensembles. Front Appl Math Stat. (2019) 5:3. doi: 10.3389/fams.2019.00003
3. Bishop CH, Whitaker JS, Lei L. Gain form of the ensemble transform Kalman filter and its relevance to satellite data assimilation with model space ensemble covariance localization. Mon Wea Rev. (2017) 145:4575–92. doi: 10.1175/MWR-D-17-0102.1
4. Posselt DJ, Bishop CH. Nonlinear parameter estimation: comparison of an ensemble Kalman smoother with a Markov chain Monte Carlo algorithm. Mon Wea Rev. (2012) 140:1957–74. doi: 10.1175/MWR-D-11-00242.1
6. Wang X, Hamill TM, Whitaker JS, Bishop CH. A comparison of hybrid ensemble transform Kalman filter–optimum interpolation and ensemble square root filter analysis schemes. Mon Wea Rev. (2007) 135:1055–76. doi: 10.1175/MWR3307.1
7. Lei L, Whitaker JS, Bishop C. Improving assimilation of radiance observations by implementing model space localization in an ensemble Kalman filter. J Adv Model Earth Syst. (2018) 10:3221–32. doi: 10.1029/2018MS001468
Keywords: ensemble Kalman filter, modified Kalman gain, eigenvalue form of Kalman gain, GETKF, ETKF
Citation: Bishop CH, Whitaker JS and Lei L (2020) Commentary: On the Efficiency of Covariance Localisation of the Ensemble Kalman Filter Using Augmented Ensembles. Front. Appl. Math. Stat. 6:2. doi: 10.3389/fams.2020.00002
Received: 26 November 2019; Accepted: 29 January 2020;
Published: 17 March 2020.
Edited by:
Raluca Eftimie, University of Dundee, United KingdomReviewed by:
Xin Tong, National University of Singapore, SingaporeCopyright © 2020 Bishop, Whitaker and Lei. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Craig H. Bishop, craig.bishop@unimelb.edu.au