Ensemble Rauch-Tung-Striebel Smoothing

The EnRTS algorithm provides a smoothing counterpart to the EnKF, based on the RTS smoother. It makes an empirical approximation of the RTS smoother gain, which is applied to the ensemble in a backwards pass to obtain a smoothing distribution. Note that, based on the outputs of the EnKF, the EnRTS step is entirely deterministic.

`cuthbertlib.ensemble_kalman.smoothing`

Implements the Ensemble Rauch-Tung-Striebel (EnRTS) smoother update.

Cf. Raanes (2016).

`update(filtered_ensemble, predicted_ensemble, next_smoothed_ensemble)`

Applies one EnRTS smoother update step.

Parameters:

Name	Type	Description	Default
`filtered_ensemble`	`Array`	Filtered ensemble at time t, shape (N, x_dim).	required
`predicted_ensemble`	`Array`	Paired forecast ensemble at time t + 1, shape (N, next_x_dim).	required
`next_smoothed_ensemble`	`Array`	Smoothed ensemble at time t + 1, shape (N, next_x_dim).	required

Returns:

Type	Description
`tuple[Array, Array]`	Tuple of the smoothed ensemble at time t and the EnRTS gain.

Source code in cuthbertlib/ensemble_kalman/smoothing.py

def update(
    filtered_ensemble: Array,
    predicted_ensemble: Array,
    next_smoothed_ensemble: Array,
) -> tuple[Array, Array]:
    """Applies one EnRTS smoother update step.

    Args:
        filtered_ensemble: Filtered ensemble at time t, shape (N, x_dim).
        predicted_ensemble: Paired forecast ensemble at time t + 1,
            shape (N, next_x_dim).
        next_smoothed_ensemble: Smoothed ensemble at time t + 1,
            shape (N, next_x_dim).

    Returns:
        Tuple of the smoothed ensemble at time t and the EnRTS gain.
    """
    filtered_dev = filtered_ensemble - jnp.mean(filtered_ensemble, axis=0)
    predicted_dev = predicted_ensemble - jnp.mean(predicted_ensemble, axis=0)

    gain = filtered_dev.T @ jnp.linalg.pinv(predicted_dev.T)

    smoothed_ensemble = (
        filtered_ensemble + (next_smoothed_ensemble - predicted_ensemble) @ gain.T
    )
    return smoothed_ensemble, gain