Methods of Spatial Data Processing Based on Bayesian Approach for Environmental Monitoring

Abstract

One of the important tasks of monitoring of environment is the problem of obtaining the values of environmental parameters on a regular grid. At present, such problems are solved using all available observational data as well as the mathematical model of the process of interest to us. The mathematical formulation of the problem is included in the set of tasks of the so-called inverse modelling. If the probabilistic formulation of the problem is considered, the Bayesian approach is applied. This approach is used in popular algorithms, such as the ensemble Kalman filter, the ensemble Kalman smoothing, the particle method. The report provides a brief overview of modern methods, as well as approaches to their practical implementation.

Keywords

• data assimilation
• ensemble Kalman filter
• satellite data

1. Introduction

2. The inverse modeling problem

.

Depending on which time period is needed for state estimation, it is possible to define three estimation problems: p(xl | yk,1), k  l - forecast, p(xk | yk,1) - filtration, p(xk,0 | yk,1) - smoothing, wheгe xk,0  {xk , xk1, , x0}, yk,1  {yk , , y1} . Notation and definitions are taken from the review [1].

In the linear Gaussian case, the solution of the filtering problem is the Kalman filter, the solution of the smoothing problem is the Kalman smoothing. In [3], it was proposed to use the Monte Carlo method for solving filtering and smoothing problems, the so-called ensemble filtering and smoothing algorithms. In the ensemble Kalman filter in the nonlinear case, the prediction of the covariance matrix is performed using a nonlinear model, and the condition of Gaussianity is violated. Also in this case, the estimate at the analysis step is an approximation of the estimate of the minimum variance. To solve the smoothing problem in the nonlinear case, iterative smoothing algorithms using ensembles are currently being considered [4]. As noted in [4], iterative methods provide a good approximation for weakly non-linear models. In the nonlinear non-Gaussian case, the particle method is used, which is also based on the Bayesian approach [1, 2].

4. The ensemble Kalman filter

The ensemble Kalman filter was proposed by Evensen G. [3]. Consider a nonlinear dynamic system as a process equation

xtk  f (xtk1)  ηtk1 and an observation equation

yk  h(xtk )  εtk , where h is the observation operator, generally speaking, non-linear, transforming the forecast values into the observable variable, ηtk 1 is the ‘model noise’ vector, εtk is the vector of observation errors, xtk is the vector of estimated vEa[rεiatkblεetksT ] atRtk ,tEhe[ηtk1m ηotmk1enTt]  Qoftk1 . Wtimeewill taks,sumeεtkxtk toanbde a ‘tηrutk e1’ valauree. Gaussian random variables: The stochastic ensemble Kalman filter consists of an ensemble of forecasts {xkf ,n , n  1, , N}

xkf ,n  f (xak,n1)  ηnk1 xka,n  xkf ,n  Kk (ynk  εnk  h(xkf ,n )) . and an ensemble of analysis {xka,n , n  1, , N}

The ensembles (1) and (2) provide a sample of ‘true’ values. Here the sample mean will be the optimal estimate, and the deviations from the mean will be the ensembles of the analysis and forecast errors, respectively. To implement the ensemble version of the Kalman filter algorithm, an ensemble of ob1seNrvation errors {εkn , n  1, , N} and ennosiseem{bηlekn1o,fnfor1e,cas,tNe}rr,oErs[ηkn{d1xηkf ,knn1Tx]kf ,n Qkxkfa,nre, nspe1c,ifie,dN.}K, kwihsearemaxtkfr,inxofNthne1fxorkf m,n , and an ensemble of model

Kk  Pkf HkT (Hk Pkf HkT  Rk )1 , where Pkf and Rk are matrices that are estimated by the ensemble averaging as

Pkf N11 nN1 dxkf ,n dxkf ,n T , Rk N11 nN1 εkn εkn T , And Hk is the linearized operator of h(xkf ,n ) with respect to x f ,n : h(xk )  h(xkf ,nk)  H ε k kf .

Formulas (1) - (2) are a stochastic ensemble Kalman filter. The deterministic ensemble Kalman filter (analysis step) consists of an equation for the mean

xka,n  xkf ,n  Kk (ynk  h(xkf ,n )) and an estimate of the ensemble of analysis errors such that the corresponding covariance matrix satisfies the Kalman filter equation Pka  (I  Kk Hk )Pkf [ 5 ].

5. The problem of optimal smoothing. Ensemble smoothing (EnKS - ensemble Kaman smoother)

In the classical optimal estimation theory, the problem of optimal filtering is the problem of obtaining an optimal estimate at the end of the considered time interval. Also, there is an optimal smoothing problem, which is the problem of obtaining the optimal estimate for the given time interval. We suppouse, that an optimal estimate minimizes the trace of the covariance matrix of estimation errors [6].

Ensemble algorithms are also applied for the optimal smoothing problem. In [3] was shown that the problem of ensemble smoothing can be solved sequentially in time by using the available data at each time step to estimate values over the given time interval. In this case, the formulas for the optimal estimate are similar to those of the ensemble Kalman filter. With this variant of ensemble smoothing (EnKS), calculations computing in the opposite direction of time are not required. It should be noted that the EnKS algorithm is equivalent to the variational data assimilation algorithm with the use of the corresponding weight matrices (4DVAR) [1, 2]. (1) (2)

6. The parameters estimation in the data assimilation algorithm

Consider the equation of process as observations are xtk  f (xtk1,αtk1)  ηtk1,

yk  h(xtk ,αtk )  εtk , where αtk is the vector of parameters. We assume that the parameter does not change with time: αtk1  αtk . Consider the generalized estimation problem for vector z  [x,α]T . Omitting the intermediate calculations, we write down the result of the estimation procedure in general: xa  x f  PxxhxT (hxPxxhxT  R)1[y  h(x f ,α f )]  Px hT (hxPxxhxT  R)1[y  h(x f ,α f )] , αa  α f  P xhxT (hxPxxhxT  R)1[y  h(x f ,α f )]  P hT (hxPxxhxT  R)1[y  h(x f ,α f )] .

In these formulas, the index ‘k’ is omitted. The analysis step is considered, the index ‘a’ means analysis, the index ‘f’ means the prediction, P x - the cross-covariance of errors x and α , P - the covariance matrix of errors α , hx and h - linearized operators with respect to x and α . If h does not depend on α , the estimation xa is carried out using the same formula as in the classical Kalman filter [6]. In modern works on data assimilation, such an approach is applied to the estimation of greenhouse gas fluxes [7].

7. Variational approach to the data assimilation problem

The variational approach to the data assimilation problem is very popular and is used in the operational data assimilation systems in the leading prognostic centers [1, 2]. The data analysis algorithm called 3DVAR (3Dimensional VARiational) consists in finding a vector xa that minimizes the functional

2J(x)  (x  xb )T B1(x  xb ) [y0  h(x)]T R1[y0  h(x)].

In this formula, B and R are the covariance matrices of forecast and observation errors, respectively. The operator h is, generally speaking, nonlinear, transforms prognostic value to the observations (observed variables).

Minimum of functional is the solution of the equation

J(xa )  0.

If the gradient of the functional is written in the following form:

J(x)  B1(x  xb )  HTR1H(x  xb )  HTR1{y0  h(xb )}, then the solution of the problem is

xa  xb  (B1  HTΡ1H){y0  h(xb )}.

Where H is the linearized operator h with respect to xb . This formula is equivalent to the analysis step of the Kalman filter algorithm [1, 2, 6].

The algorithm 4DVAR is a generalization of 3DVAR to the space-time case. 4DVAR allows you to use observations from a certain time interval (tn  t0 ) . It is assumed that observations {yi0,i  0, , N} are specified at this interval. It is required to find a minimum of functional

1 1 N J[x(t0 )]  [x(t0 )  xb (t0 )]T B01[x(t0 )  xb (t0 )]  [h(xi )  yi0]T R1[h(xi )  yi0]

2 2 i0 provided that

x(tn )  Mn[x(t0 )].

In this optimization problem, the control variable x(t0 ) is the value at the initial time. That is, the algorithm 4DVAR finds an initial condition such that a forecast with this initial value best approximates observations in a given time interval. The minimum functional in 4DVAR is searched by iterative methods (for example, the quasi-Newton method). To implement the iterative process, we need to estimate the gradient of the functional J(x). The gradient of functional is calculated using the conjugate linearized model. At present, there are data assimilation algorithms based on 4DVAR, with use the ensemble approach to estimate changes in covariances over time [1, 2].

8. Problems of implementation of ensemble algorithms

Tasks of data assimilation, the estimation of the parameters in the environment modelling has a very large dimension, requires huge computing resources. The use of ensemble algorithms partially solves this problem, but still the problem remains extremely time-consuming. To carry out the analysis step of ensemble Kalman filter, algorithms of transformation of the ensemble of forecasts are used to obtain the ensemble of the analysis. In this case, the analysis is performed only for the ensemble mean value, and then the ensemble of analyses is calculated. There are deterministic and stochastic variants of such algorithms for realization of both ensemble Kalman filter and ensemble Kalman smoother [ 5 ].

A variant of a stochastic Kalman filter with ensemble transformations is the ensemble π-algorithm [8-10]. This algorithm performs operations with matrices of the order of the ensemble dimension. In addition, this algorithm can be used to implement ensemble smoothing (EnKS). An important property of the algorithm is its locality. The estimate the desired value can be implemented independently for the given subdomains. This property can be used for estimation of parameters locally in the given region.

When using ensembles in data assimilation, there are a number of problems associated with a small sample size (ensemble size). In particular, it is a problem of high correlation values of forecast errors at long distances. In this case, the so-called localization may be used, namely, the element-by-element multiplication of covariance matrices by the function decreasing with the distance [ 1,5 ].

The small sample size is also one of the reasons for the divergence of the algorithm. The algorithms of the ensemble Kalman filter (EnKF) and ensemble Kalman smoothing (EnKS) use the ‘inflation factor’. In this case the elements of the covariance matrix are multiplied by a certain multiplier (multiplicative inflation), or an additional perturbation is added to the ensemble of perturbations (additive inflation). It should be noted that these techniques are also used in the iterative smoothing algorithms using ensembles, as well as in the particle method using ensembles [ 1, 2, 5 ].

Consider the features of the assimilation of satellite data [11]. Satellite data have a number of distinctive features compared to ground-based and upper-air observations. These features include the following:

(a) The satellite measures radiation information, so it is necessary to convert model variables into observable variables using the radiation transfer equation. Thus, the observation operator h is nonlinear.

b) Observations are received continuously over time.

C) Observation errors correlate, the matrix of observational errors has large dimension and it is non-diagonal. d) There is a of systematic error in observations (bias). This makes it difficult to use the formulas of the above algorithms.

Satellite data can be used in two ways: directly, in the form of radiation data, or it is required to obtain meteorological values from satellite data in advance and then use them in the data assimilation procedure. The first option is more preferable for a number of reasons: (a) The transformation of radiation data uses forecast information, so the errors of the recovered data and forecast will be correlated. b) There is the problem of estimating the covariance matrix of observation errors.

Satellite data is a huge amount of information. However, as noted in a number of papers, a large number of observations do not always give the best result, since these data are an array of correlating random variables and therefore they are not too informative [11].

9. The data assimilation in the environment modelling

The modern environmental monitoring is based on the mathematical modeling. In particular, models of the propagation of greenhouse gases, active gas constituents, and aerosols in the atmosphere are considered.

Information on all these substances is measured using both ground-based measuring instruments and satellites. A large number of investigations are currently devoted to the environmental data assimilation problem, using the basic mathematical apparatus developed in the meteorological data assimilation [7].

An important task in monitoring of the environment is to assess the fluxes from the earth's surface of greenhouse gases with the help of the data assimilation system [7].

The problems of the environment modeling are currently using models of the atmosphere with a "chemical" block. As an example, the popular model WRF-Chem [12] may be considered. Modern databases of meteorological fields (ERA-interim reanalysis [13], for example) are used for numerical experiments. ECMWF has created a MACC database (Monitoring Atmospheric Composition and Climate), which includes active gas components, aerosols and greenhouse gases [14]. In the regional modelling, an important problem is the determination of boundary values. Reanalysis data may be used for this purpose. Also, in the literature, methods of data assimilation for determination values at the boundary of the region are considered [7].

10.

Conclusion

The report provides a brief overview of current approaches to the problem of data assimilation in the environment modelling. The main aspects of practical realization of this task are considered. [2] [3] [4]

Carrassi A., Bocquet M., Bertino L., Evensen G. Data assimilation in the geosciences: An overview of methods, issuers and perspectives // Wiley interdisciplinary reviews: Climate Change. 2018. Vol. 131. Issue5, e535.

Nakamura G., Potthast R. Inverse Modeling. IOP Publishing. 2015.

Evensen, G. Data assimilation. The ensemble Kalman filter.Berlin Heideberg: Spriger-Verlag, 2009.

Evensen G. Analysis of iterative ensemble smoother for solving inverse problems //Computational Geosciences. 2018. Vol. 22. P.885-908.

Houtekamer, H.L. Zhang, F. Review of the ensemble Kalman filter for atmospheric data assimilation //

References

[5]