11  Ensemble Kalman Filter

The idea behind the EnKF is

  1. permit fully nonlinear process and observation models;
  2. avoid gradient calculations for linearization, as in the EKF,
  3. by replacing mean and covariance with empirical, ensemble averaging.

The EnKF is obtained by - replacing the exact covariance \(P\) by the ensemble sample covariance, and - adding noise to the data in order to avoid a shrinking of the ensemble spread and to obtain the correct filtering covariance in the limit

Here are the steps for the so-called stochastic EnKF, where we add an artificial, random perturbation to the observations, and we assume we have a linear observation operator, \(H.\)

11.1 Stochastic EnKF - linear observation operator

11.1.1 Prediction/Forecast

\[\begin{align} \hat{v}_{k+1}^{n} &= \Psi ( {v}_{k}^{n} ) + \xi_{k}^{n} , \quad n=1, \ldots, N, \\ \hat{m}_{k+1} & = \frac{1}{N} \sum_{i=1}^{N} \hat{v}_{k+1}^{n} , \\ \hat{C}_{k+1} & = \frac{1}{N-1} \sum_{i=1}^{N} \left( \hat{v}_{k+1}^{n} - \hat{m}_{k+1} \right) \left( \hat{v}_{k+1}^{n} - \hat{m}_{k+1} \right)^{\mathrm{T}}. \end{align}\]

11.1.2 Correction/Analysis

\[\begin{align} {S}_{k+1} &= H \hat{C}_{k+1} H^{\mathrm{T}} + \Gamma ,\\ {K}_{k+1} &= \hat{C}_{k+1} H^{\mathrm{T}} {S}_{k+1}^{-1} , \\ {y}_{k+1}^{n} &= {y}_{k+1} + \eta^n_{k+1} , \quad n=1, \ldots, N, \\ {v}_{k+1}^{n} &= (I - K_{k+1}H )\hat{v}_{k+1}^{n} + K_{k+1} y^n_{k+1} , \quad n=1, \ldots, N. \end{align}\]

Alternatively, defining the innovation \(d = {y}_{k+1}^{n} - H \hat{v}_{k+1}^{n},\) we can write the state update more simply as

\[ {v}_{k+1}^{n} = \hat{v}_{k+1}^{n} + K_{k+1} d . \]

In words:

  1. For a given \(N_e \in \mathbb{N}\) generate i.i.d. ensemble of states random variables from the distribution of \(X(0).\)
  2. For \(t \in \mathbb{N}\) recursively repeat the following steps:
  • Advance each ensemble member in time, using the nonlinear state equation with independently generated random state noise
  • Compute the forecast sample mean and the forecast sample covariance
  • Compute the sample Kalman gain
  • Add additional perturbation to the observation vector \(Y\) using independently generated random variables \(\eta(t)\)
  • Update each forecast ensemble member

Burgers et al. [1998] shows that without the data perturbation, the covariance of the ensemble would go to the zero matrix as \(t\) goes to infinity. The data perturbation also guarantees that the relation between the forecast sample covariance and the analysis sample covariance is analogous to the relation between the forecast and analysis covariances in the standard KF.

11.2 Full nonlinear formulation of the ensemble Kalman filter

There are many ways to formulate the EnKF. Following Vetra-Carvalho, et al (Tellus A, 2018), we express the filter in terms of the anomalies of state and observations. This is indispensable for fully nonlinear state and measurement models,

\[\begin{align} x_{k+1}^{n} &= \Psi ( {x}_{k}^{n} ) + w_{k}^{n} , \quad n=1, \ldots, N_e, \\ {y}_{k+1} &= \mathcal{H} (x_{k+1}) + v_{k+1}. \end{align}\]

To fix notation:

  • state forecast \(X^{\mathrm{f}},\) dimension \((N_t \times N_x)\)
  • ensemble state forecast \(\mathbf{X}^{\mathrm{f}},\) dimension \((N_t \times N_x \times N_e)\)
  • observation, \(Y,\) dimension \((N_t \times N_y)\)
  • ensemble state anomaly, \[\mathbf{X}' = \frac{1}{\sqrt{N_e - 1}} \left(\mathbf{X} - \overline{{X}}\right),\] dimension \((N_t \times N_x \times N_e)\) with \(\overline{{X}} = (1/N_e) \sum_{e=1}^{N_e} {X}_e\)
  • ensemble observation anomaly,
    \[\mathbf{Y}' = \frac{1}{\sqrt{N_e - 1}} \left( \mathcal{H}(\mathbf{X}) - \overline{ \mathcal{H}(\mathbf{X} ) } \right),\] dimension \((N_t \times N_y \times N_e)\) with \(\overline{ \mathcal{H}( \mathbf{X} )} = (1/N_e)\sum_{e=1}^{N_e} \mathcal{H}({X}_e).\)

Then, the Kalman analysis update is

\[ \mathbf{X}^\mathrm{a} = \mathbf{X}^\mathrm{f} + \mathbf{X}'(\mathbf{Y}')^\mathrm{T} S^{-1} D, \]

with

\[\begin{align} S &= \mathbf{Y}'(\mathbf{Y}')^\mathrm{T} + R & \quad \text{(observation covariance)}, \\ D &= ( \mathbf{Y} + \mathbf{y}) - \mathcal{H}(\mathbf{X}) & \quad \text{(innovation)}, \end{align}\]

where \(y \sim \mathcal{N}(0,R)\) is the stochastic perturbation, and \(R\) is the measurement noise covariance matrix.

Or, defining the Kalman gain matrix as

\[ K = \mathbf{X}'(\mathbf{Y}')^\mathrm{T} S^{-1}, \]

we obtain the classical KF update,

\[ \mathbf{X}^\mathrm{a} = \mathbf{X}^\mathrm{f} + K D. \]

11.3 Summary of EnKF properties

  • EnKF represents error statistics by ensembles of (nonlinear) model and (nonlinear) measurement realizations.
  • EnKF performs sequential DA that processes measurements recursively in time.
  • EnKF is suitable for weather-prediction and any other complex, chaotiic dynamic systems.
  • Error propagation is nonlinear (see point 1).
  • Filter update is linear and computed in the low rank, ensemble subspace.
  • EnKF does not require any gradients, adjoints, linearizations.