An Efficient Parallel Implementation of the Ensemble Kalman Filter Based on Shrinkage Covariance Matrix Estimation

Groundwater contamination assessment can be useful in taking proper actions during environmental emergency. A three-dimensional deterministic model was taken into consideration to simulate the advective-diffusive transport of non-conservative contaminant in groundwater. Multiple stochastic data assimilation techniques, such as, ensemble Kalman filter (EnKF), local ensemble transform Kalman filter (LETKF), and the global form of the LETKF, denoted as GETKF were applied to the model. The results show that data assimilation improved contaminant concentration prediction. The EnKF method reduced the root-mean-square-error (RMSE) of the contaminant prediction from 12.5 mg/L to 1.31 mg/L, whereas the LETKF reduced that to 0.46 mg/L and GETKF reduced that to 0.38 mg/L. Approximately 89.48%, 96.30% and 96.82% improvement were made by EnKF, LETKF, and GETKF, respectively. The sensitivity analysis suggest that these data assimilation techniques are very sensitive to the observation noise, process noise, and ensemble size.

A particle filter (PF) is an ensemble data assimilation method that does not assume Gaussian error distributions. Recent studies proposed local PFs (LPFs), which use localization as in the ensemble Kalman filter, to apply the PF for high-dimensional dynamics efficiently. Among others, Penny and Miyoshi developed an LPF in the form of the ensemble transform matrix of the Local Ensemble Transform Kalman Filter (LETKF). The LETKF has been widely accepted for various geophysical systems including numerical weather prediction (NWP) models. Therefore, implementing consistently with an existing LETKF code is useful. This study developed a software platform for the LPF and its Gaussian mixture extension (LPFGM) by making slight modifications to the LETKF code with a simplified global climate model known as Simplified Parameterizations, Primitive Equation Dynamics (SPEEDY). A series of idealized twin experiments were accomplished under the ideal model assumption. With large inflation by the relaxation to prior spread, the LPF showed stable filter performance with dense observations but became unstable with sparse observations. The LPFGM showed more accurate and stable performances than the LPF with both dense and sparse observations. In addition to the relaxation parameter, regulating the resampling frequency and the amplitude of Gaussian kernels was important for the LPFGM. With a spatially inhomogeneous observing network, the LPFGM was superior to the LETKF in sparsely observed regions where the background ensemble spread and non-Gaussianity are larger. The SPEEDY-based LETKF, LPF, and LPFGM systems were available as open-source software on Github (https://github.com/skotsuki/speedy-lpf) and can be adapted to various models relatively easily like the LETKF.

A particle filter (PF) is an ensemble data assimilation method that does not assume Gaussian error distributions. Recent studies proposed local PFs (LPFs), which use localization as in the ensemble Kalman filter, to apply the PF for high-dimensional dynamics efficiently. Among others, Penny and Miyoshi developed an LPF in the form of the ensemble transform matrix of the Local Ensemble Transform Kalman Filter (LETKF). The LETKF has been widely accepted for various geophysical systems including numerical weather prediction (NWP) models. Therefore, implementing consistently with an existing LETKF code is useful. This study developed a software platform for the LPF and its Gaussian mixture extension (LPFGM) by making slight modifications to the LETKF code with a simplified global climate model known as Simplified Parameterizations, Primitive Equation Dynamics (SPEEDY). A series of idealized twin experiments were accomplished under the ideal model assumption. With large inflation by the relaxation to prior spread, the LPF showed stable filter performance with dense observations but became unstable with sparse observations. The LPFGM showed more accurate and stable performances than the LPF with both dense and sparse observations. In addition to the relaxation parameter, regulating the resampling frequency and the amplitude of Gaussian kernels was important for the LPFGM. With a spatially inhomogeneous observing network, the LPFGM was superior to the LETKF in sparsely observed regions where the background ensemble spread and non-Gaussianity are larger. The SPEEDY-based LETKF, LPF, and LPFGM systems were available as open-source software on Github (https://github.com/skotsuki/speedy-lpf) and can be adapted to various models relatively easily like the LETKF.

Ensemble methods present a practical framework for parameter estimation, performance prediction, and uncertainty quantification in subsurface flow and transport modeling. In particular, the ensemble Kalman filter (EnKF) has received significant attention for its promising performance in calibrating heterogeneous subsurface flow models. Since an ensemble of model realizations is used to compute the statistical moments needed to perform the EnKF updates, large ensemble sizes are needed to provide accurate updates and uncertainty assessment. However, for realistic problems that involve large-scale models with computationally demanding flow simulation runs, the EnKF implementation is limited to small-sized ensembles. As a result, spurious numerical correlations can develop and lead to inaccurate EnKF updates, which tend to underestimate or even eliminate the ensemble spread. Ad hoc practical remedies, such as localization, local analysis, and covariance inflation schemes, have been developed and applied to reduce the effect of sampling errors due to small ensemble sizes. In this paper, a fast linear approximate forecast method is proposed as an alternative approach to enable the use of large ensemble sizes in operational settings to obtain more improved sample statistics and EnKF updates. The proposed method first clusters a large number of initial geologic model realizations into a small number of groups. A representative member from each group is used to run a full forward flow simulation. The flow predictions for the remaining realizations in each group are approximated by a linearization around the full simulation results of the representative model (centroid) of the respective cluster. The linearization can be performed using either adjoint-based or ensemble-based gradients. Results from several numerical experiments with two-phase and three-phase flow systems in this paper suggest that the proposed method can be applied to improve the EnKF performance in large-scale problems where the number of full simulation is constrained.

[1] In ensemble Kalman filter, space localization is used to reduce the impact of long-distance sampling errors in the ensemble estimation of the forecast error covariance. When two variables are not physically correlated, their error covariance is still estimated by the ensemble and, therefore, it is dominated by sampling errors. We introduce a “variable localization” method, zeroing out such covariances between unrelated variables to the problem of assimilating carbon dioxide concentrations into a dynamical model using the local ensemble transform Kalman filter (LETKF) in an observing system simulation experiments (OSSE) framework. A system where meteorological and carbon variables are simultaneously assimilated is used to estimate surface carbon fluxes that are not directly observed. A range of covariance structures are explored for the LETKF, with emphasis on configurations allowing nonzero error covariance between carbon variables and the wind field, which affects transport of atmospheric CO2, but not between CO2 and the other meteorological variables. Such variable localization scheme zeroes out the background error covariance among prognostic variables that are not physically related, thus reducing sampling errors. Results from the identical twin experiments show that the performance in the estimation of surface carbon fluxes obtained using variable localization is much better than that using a standard full covariance approach. The relative improvement increases when the surface fluxes change with time and model error becomes significant.

<p>Observations of the Venus Climate Orbiter “Akatsuki” provide us with horizontal distributions of the horizontal winds derived from cloud tracking of the Ultraviolet Imager (UVI) and of temperature observed by the Longwave Infrared Camera (LIR). However, these observations are limited in altitude, local time (day or night side), and frequency. Then it is difficult to elucidate the general circulation of the Venus atmosphere, including various temporal and spatial scales, only from observations. In this study, we produced a Venus dataset (analysis) that has high temporal and spatial resolutions by assimilating horizontal winds derived by the Akatsuki observations. At the top of the cloud layer of Venus, there are planetary-scale atmospheric waves that are excited by the solar heating and move with the sun, called the thermal tides. In this presentation, we focused on thermal tides to verify the analysis.</p> <p>We use the Venus atmospheric data assimilation system “ALEDAS-V" (Sugimoto et al., 2017) [1] for assimilation and the Venus atmospheric general circulation model “AFES-Venus" (Sugimoto et al., 2014) [2] for ensemble forecasts. AFES-Venus is a full nonlinear dynamical GCM on the assumption of hydrostatic balance, designed for the Venus atmosphere. ALEDAS-V uses the Local Ensemble Transform Kalman Filter, and is the first data assimilation system for the Venus atmosphere. We assimilated the cloud top (~70km) zonal and meridional winds obtained by tracking morphology, using Akatsuki UVI data (Horinouchi et al., 2021) [3] from September 1st to December 31st, 2018. The assimilation data (analysis) from October 1st to November 30th, 2018, is analyzed, because the root-mean-square-deviations (RMSD) from FR (free run; the case without data assimilation) are stable.</p> <p>Figures (a) and (d) show the observed zonal and meridional winds, respectively. The zonal wind has a local minimum near 11 LT (local time) around the equator (Figure a). The meridional wind is the weakest at the equator and increases with latitude, and the amplitude is maximum around noon (Figure d) in the local time direction. Note that these winds obtained from observations exist only the dayside equatorward of 50° latitudes (Figure a and d).</p> <p>Figures (b) and (e) show the deviations from the zonal means of zonal and meridional winds at an altitude of 70 km in the FR, respectively. For zonal wind, diurnal (zonal wavenumber 1) and semidiurnal (zonal wavenumber 2) tides are dominant at latitudes poleward and equatorward of 30, respectively (Figure b). The zonal wind deviation has a local minimum at 14-15LT, which is ~ 2 hours behind the observation (Figures a and b). The meridional wind deviation is polar and equatorial on the dayside and nightside, respectively (Figure e), and this distribution is consistent with Akatsuki's observation (Figure d).</p> <p>Figures (c) and (f) show the zonal and meridional winds as a result of assimilation, respectively. The zonal wind in the equatorial region have a local minimum near 11 LT. The assimilation improved the semidiurnal tide closer to the observations (Figures a and c). The meridional wind is not so different from FR. This is probably because FR was originally very similar to observations (Figures d and f). These results are consistent with a previous study by Sugimoto et al. (2019) [4]. In addition, while the observed winds exist only on the dayside, the results of assimilation show that the horizontal winds field is modified significantly even on the nightside. It is suggested that spatially limited data assimilation can improve the general circulation of GCM.</p> <p>In the future work, we are planning to release the assimilation dataset as the “objective analysis data” of Venus for the first time in the world.</p> <p><img src="" alt="" width="743" height="482" /></p> <p>[1] Sugimoto, N., et al. Development of an ensemble Kalman filter data assimilation system for the Venusian atmosphere. Scientific Reports 7(1), 9321 (2017).</p> <p>[2] Sugimoto, N., et al. Baroclinic instability in the Venus atmosphere simulated by GCM. J. Geophys. Res. Planets 119, 1950–1968 (2014).</p> <p>[3] Horinouchi, T., et al. Venus Climate Orbiter Akatsuki Cloud Motion Vector Data Set v1.0, JAXA Data Archives and Transmission System (2021).</p> <p>[4] Sugimoto, N., et al. Impact of data assimilation on thermal tides in the case of Venus Express wind observation. Geophys. Res. Lett. 46, 4573–4580 (2019).</p>

The ensemble Kalman filter (EnKF) is a popular estimation technique in the geosciences. It is used as a numerical tool for state vector prognosis and parameter estimation. The EnKF can, for example, help to evaluate the geothermal potential of an aquifer. In such applications, the EnKF is often used with small or medium ensemble sizes. It is therefore of interest to characterize the EnKF behavior for these ensemble sizes. For seven ensemble sizes (50, 70, 100, 250, 500, 1,000, and 2,000) and seven EnKF variants (damped, iterative, local, hybrid, dual, normal score, and classical EnKF), we computed 1,000 synthetic parameter estimation experiments for two setups: a 2‐D tracer transport problem and a 2‐D flow problem with one injection well. For each model, the only difference among synthetic experiments was the generated set of random permeability fields. The 1,000 synthetic experiments allow to calculate the probability density function of the root‐mean‐square error (RMSE) of the characterization of the permeability field. Comparing mean RMSEs for different EnKF variants, ensemble sizes and flow/transport setups suggests that multiple synthetic experiments are needed for a solid performance comparison. In this work, 10 synthetic experiments were needed to correctly distinguish RMSE differences between EnKF variants smaller than 10%. For detecting RMSE differences smaller than 2%, 100 synthetic experiments were needed for ensemble sizes 50, 70, 100, and 250. The overall ranking of the EnKF variants is strongly dependent on the physical model setup and the ensemble size.

Recent numerical weather prediction systems have significantly improved medium-range forecasts by implementing hybrid background error covariance, for which climatological (static) and ensemble-based (flow-dependent) error covariance are combined. While the hybrid approach has been investigated mainly in variational systems, this study aims at exploring methods for implementing the hybrid approach for the local ensemble transform Kalman filter (LETKF). Following Kretschmer et al., the present study constructed hybrid background error covariance by adding collections of climatological perturbations to the forecast ensemble. In addition, this study proposes a new localization method that attenuates the ensemble perturbation (Z-localization) instead of inflating observation error variance (R-localization). A series of experiments with a simplified global atmospheric model revealed that the hybrid LETKF resulted in smaller forecast errors than the LETKF, especially in sparsely observed regions. Due to the larger ensemble enabled by the hybrid approach, optimal localization length scales for the hybrid LETKF were larger than those for the LETKF. With the LETKF, the Z-localization resulted in similar forecast errors as the R-localization. However, Z-localization has an advantage in enabling us to apply different localization scales for flow-dependent perturbation and climatological static perturbations with the hybrid LETKF. The optimal localization for climatological perturbations was slightly larger than that for flow-dependent perturbations. This study also proposes optimal eigendecomposition (OED) ETKF formulation to reduce computational costs. The computational expense of the OED ETKF formulation became significantly smaller than that of standard ETKF formulations as the number of climatological perturbations was increased beyond a few hundred.

Abstract. Data assimilation is the process of combining real-world observations with a modelled geophysical field. The increasing abundance of satellite retrievals of atmospheric trace gases makes chemical data assimilation an increasingly viable method for deriving more accurate analysed fields and initial conditions for air quality forecasts. We implemented a three-dimensional optimal interpolation (OI) scheme to assimilate retrievals of NO2 tropospheric columns from the Ozone Monitoring Instrument into the Danish Eulerian Hemispheric Model (DEHM, version V2009.0), a three-dimensional, regional-scale, offline chemistry-transport model. The background error covariance matrix, B, was estimated based on differences in the NO2 concentration field between paired simulations using different meteorological inputs. Background error correlations were modelled as non-separable, horizontally homogeneous and isotropic. Parameters were estimated for each month and for each hour to allow for seasonal and diurnal patterns in NO2 concentrations. Three experiments were run to compare the effects of observation thinning and the choice of observation errors. Model performance was assessed by comparing the analysed fields to an independent set of observations: ground-based measurements from European air-quality monitoring stations. The analysed NO2 and O3 concentrations were more accurate than those from a reference simulation without assimilation, with increased temporal correlation for both species. Thinning of satellite data and the use of constant observation errors yielded a better balance between the observed increments and the prescribed error covariances, with no appreciable degradation in the surface concentrations due to the observation thinning. Forecasts were also considered and these showed rather limited influence from the initial conditions once the effects of the diurnal cycle are accounted for. The simple OI scheme was effective and computationally feasible in this context, where only a single species was assimilated, adjusting the three-dimensional field for this compound. Limitations of the assimilation scheme are discussed.

The assimilation of surface‐sensitive microwave radiance observations over land can improve numerical weather prediction accuracy at Environment and Climate Change Canada. However, the benefits of these observations are limited by large errors in surface emissivity and skin temperature, along with challenges accounting for these errors in the data assimilation system. Typically, the contribution of surface emissivity error is treated as a deficiency in the observation operator and represented by inflating a fixed observation error covariance matrix (R). However, this study demonstrates the error contributions of surface emissivity in observation space for each radiance observation profile can vary substantially due to variations in the surface emissivity Jacobians, even if the surface emissivity error covariances are fixed. An approach is introduced to account for this variability by representing the surface emissivity errors in the background error covariance matrix (B) and including surface emissivity as an additional analysis variable. Using an idealized one‐dimensional variational assimilation framework, the approach is compared with a reference experiment that uses a fixed R to represent the average error contribution of surface emissivity. The proposed approach results in up to 20% and 35% greater error reduction in the mid–lower tropospheric air temperature and skin temperature analyses, respectively, when compared to the reference experiment. Both methods would be mathematically equivalent in estimating air and skin temperatures if the variations in the error contribution of surface emissivity were fully represented within the R in the reference approach, which may require significant technical changes to the assimilation system. However, additional benefits of the proposed approach could potentially be realized in a 4D‐EnVar assimilation system by propagating the updated surface emissivity and within the outer loop of an incremental formulation. Additionally, it is demonstrated that well‐represented background error correlations between air and skin temperatures in B can further improve the error reductions.

Data assimilation is a methodology for estimating accurately the state of a time-evolving complex system like the atmosphere from observational data and a numerical model of the system. It has become an indispensable tool for meteorological researches as well as for numerical weather prediction, as represented by extensive use of reanalysis datasets for research purposes. New advances of data assimilation methods emerged from the 1980s. This review paper presents the theoretical background and implementation of two advanced data assimilation methods: four-dimensional variational assimilation (4D-Var) and ensemble Kalman filtering (EnKF), which currently draw much attention in the meteorological community. Recent research results in Japan on those methods are reviewed, especially on mesoscale applications of 4D-Var and tests of the local ensemble transform Kalman filter (LETKF). Comparison of 4D-Var and EnKF is also briefly discussed. An outline of the mesoscale 4D-Var system of the Japan Meteorological Agency, which is the first operational 4D-Var for a mesoscale model, is given in Appendix.

Bayesian estimation theory provides a general approach for the state estimate of linear or nonlinear and Gaussian or non-Gaussian systems. In this study, we first explore two Bayesian-based methods: ensemble adjustment Kalman filter (EAKF) and sequential importance resampling particle filter (SIR-PF), using a well-known nonlinear and non-Gaussian model (Lorenz '63 model). The EAKF, which is a deterministic scheme of the ensemble Kalman filter (EnKF), performs better than the classical (stochastic) EnKF in a general framework. Comparison between the SIR-PF and the EAKF reveals that the former outperforms the latter if ensemble size is so large that can avoid the filter degeneracy, and vice versa. The impact of the probability density functions and effective ensemble sizes on assimilation performances are also explored. On the basis of comparisons between the SIR-PF and the EAKF, a mixture filter, called ensemble adjustment Kalman particle filter (EAKPF), is proposed to combine their both merits. Similar to the ensemble Kalman particle filter, which combines the stochastic EnKF and SIR-PF analysis schemes with a tuning parameter, the new mixture filter essentially provides a continuous interpolation between the EAKF and SIR-PF. The same Lorenz '63 model is used as a testbed, showing that the EAKPF is able to overcome filter degeneracy while maintaining the non-Gaussian nature, and performs better than the EAKF given limited ensemble size.

One of the main limitations in current wave data assimilation systems is the lack of an accurate representation of the structure of the background errors. In this work, models for the observational error variance, background error variance and background error correlations are developed based on the results of previous studies. These are tested in a global wave data assimilation system and the resulting wave forecasts are verified against independent observations from buoys. Forecasts of significant wave height show substantial improvement over the Australian Bureau of Meteorology's current operational wave forecasting system. However, forecasts of peak period are not similarly improved. The regional impacts of the new assimilation scheme are found to vary on a seasonal basis. Overall, it is shown that the inclusion of a latitudinally dependent background error, and improved specification of the background and observational error variances can reduce the root-mean-square error of 24-hour forecast Significant Wave Height by almost 10%.

Previous research suggests that an ensemble Kalman filter (EnKF) data assimilation and modeling system can produce accurate atmospheric analyses and forecasts at 30–50-km grid spacing. This study examines the ability of a mesoscale EnKF system using multiscale (36/12 km) Weather Research and Forecasting (WRF) model simulations to produce high-resolution, accurate, regional surface analyses, and 6-h forecasts. This study takes place over the complex terrain of the Pacific Northwest, where the small-scale features of the near-surface flow field make the region particularly attractive for testing an EnKF and its flow-dependent background error covariances. A variety of EnKF experiments are performed over a 5-week period to test the impact of decreasing the grid spacing from 36 to 12 km and to evaluate new approaches for dealing with representativeness error, lack of surface background variance, and low-level bias. All verification in this study is performed with independent, unassimilated observations. Significant surface analysis and 6-h forecast improvements are found when EnKF grid spacing is reduced from 36 to 12 km. Forecast improvements appear to be a consequence of increased resolution during model integration, whereas analysis improvements also benefit from high-resolution ensemble covariances during data assimilation. On the 12-km domain, additional analysis improvements are found by reducing observation error variance in order to address representativeness error. Removing model surface biases prior to assimilation significantly enhances the analysis. Inflating surface wind and temperature background error variance has large impacts on analyses, but only produces small improvements in analysis RMS errors. Both surface and upper-air 6-h forecasts are nearly unchanged in the 12-km experiments. Last, 12-km WRF EnKF surface analyses and 6-h forecasts are shown to generally outperform those of the Global Forecast System (GFS), North American Model (NAM), and the Rapid Update Cycle (RUC) by about 10%–30%, although these improvements do not extend above the surface. Based on these results, future improvements in multiscale EnKF are suggested.

Abstract. Prediction of spatiotemporal chaotic systems is important in various fields, such as numerical weather prediction (NWP). While data assimilation methods have been applied in NWP, machine learning techniques, such as reservoir computing (RC), have recently been recognized as promising tools to predict spatiotemporal chaotic systems. However, the sensitivity of the skill of the machine-learning-based prediction to the imperfectness of observations is unclear. In this study, we evaluate the skill of RC with noisy and sparsely distributed observations. We intensively compare the performances of RC and local ensemble transform Kalman filter (LETKF) by applying them to the prediction of the Lorenz 96 system. In order to increase the scalability to larger systems, we applied a parallelized RC framework. Although RC can successfully predict the Lorenz 96 system if the system is perfectly observed, we find that RC is vulnerable to observation sparsity compared with LETKF. To overcome this limitation of RC, we propose to combine LETKF and RC. In our proposed method, the system is predicted by RC that learned the analysis time series estimated by LETKF. Our proposed method can successfully predict the Lorenz 96 system using noisy and sparsely distributed observations. Most importantly, our method can predict better than LETKF when the process-based model is imperfect.