naginterfaces.library.tsa.kalman_​unscented_​state

naginterfaces.library.tsa.kalman_unscented_state(y, lx, ly, f, h, x, st, data=None, spiked_sorder='C')

kalman_unscented_state applies the Unscented Kalman Filter (UKF) to a nonlinear state space model, with additive noise.

kalman_unscented_state uses direct communication for evaluating the nonlinear functionals of the state space model.

For full information please refer to the NAG Library document for g13ek

https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/g13/g13ekf.html

Parameters

yfloat, array-like, shape $\left(\text{my}\right)$

$y_t$, the observed data at the current time point.

lxfloat, array-like, shape $(\text{mx},\text{mx})$

$L_x$, such that $L_x{L}_x^T={\Sigma }_x$, i.e., the lower triangular part of a Cholesky decomposition of the process noise covariance structure. Only the lower triangular part of $lx$ is referenced.

If ${\Sigma }_x$ is time dependent, the value supplied should be for time $t$.

lyfloat, array-like, shape $(\text{my},\text{my})$

$L_y$, such that $L_y{L}_y^T={\Sigma }_y$, i.e., the lower triangular part of a Cholesky decomposition of the observation noise covariance structure. Only the lower triangular part of $ly$ is referenced.

If ${\Sigma }_y$ is time dependent, the value supplied should be for time $t$.

fcallable fxt = f(xt, data=None)

The state function, $F$.

Parameters

xtfloat, ndarray, shape $(\text{mx},n)$

$X_t$, the sigma points generated in (a). For the $j$th sigma point, the value for the $i$th parameter is held in $xt[\text{i}-1,\text{j}-1]$, for $\text{j}=1,2,\dots ,n$, for $\text{i}=1,2,\dots ,m_x$.

dataarbitrary, optional, modifiable in place

User-communication data for callback functions.

Returns

fxtfloat, array-like, shape $(\text{mx},n)$

$F\left(X_t\right)$.

For the $j$th sigma point the value for the $i$th parameter should be held in $fxt[\text{i}-1,\text{j}-1]$, for $\text{j}=1,2,\dots ,n$, for $\text{i}=1,2,\dots ,m_x$.

hcallable hyt = h(my, yt, data=None)

The measurement function, $H$ as described in (e).

Parameters

myint

$m_y$, the number of observed variables.

ytfloat, ndarray, shape $(\text{mx},n)$

$Y_t$, the sigma points generated in (d). For the $j$th sigma point, the value for the $i$th parameter is held in $yt[\text{i}-1,\text{j}-1]$, for $\text{j}=1,2,\dots ,n$, for $\text{i}=1,2,\dots ,m_x$, where $m_x$ is the number of state variables and $n$ is the number of sigma points.

dataarbitrary, optional, modifiable in place

User-communication data for callback functions.

Returns

hytfloat, array-like, shape $(my,n)$

$H\left(Y_t\right)$.

For the $j$th sigma point the value for the $i$th parameter should be held in $hyt[\text{i}-1,\text{j}-1]$, for $\text{j}=1,2,\dots ,n$, for $\text{i}=1,2,\dots ,m_y$.

xfloat, array-like, shape $\left(\text{mx}\right)$

${\hat{x}}_{t-1}$ the state vector for the previous time point.

stfloat, array-like, shape $(\text{mx},\text{mx})$

$S_t$, such that $S_{t-1}{S}_{t-1}^T=P_{t-1}$, i.e., the lower triangular part of a Cholesky decomposition of the state covariance matrix at the previous time point. Only the lower triangular part of $st$ is referenced.

dataarbitrary, optional

User-communication data for callback functions.

spiked_sorderstr, optional

If $lx$ and $st$ are spiked (i.e., have unit extent in all but one dimension, or have size $1$), $spiked\_sorder$ selects the storage order to associate with them in the NAG Engine:

spiked_sorder = $\text{'C'}$

row-major storage will be used;

spiked_sorder = $\text{'F'}$

column-major storage will be used.

Two-dimensional arrays returned from callback functions in this routine must then use the same storage order.

Returns

xfloat, ndarray, shape $\left(\text{mx}\right)$

${\hat{x}}_t$ the updated state vector.

stfloat, ndarray, shape $(\text{mx},\text{mx})$

$S_t$, the lower triangular part of a Cholesky factorization of the updated state covariance matrix.

Raises

NagValueError

(errno $11$)

On entry, $\text{mx}=⟨value⟩$.

Constraint: $\text{mx}\ge 1$.

(errno $21$)

On entry, $\text{my}=⟨value⟩$.

Constraint: $\text{my}\ge 1$.

(errno $301$)

A weight was negative and it was not possible to downdate the Cholesky factorization.

(errno $302$)

Unable to calculate the Kalman gain matrix.

(errno $303$)

Unable to calculate the Cholesky factorization of the updated state covariance matrix.

Warns

NagCallbackTerminateWarning

(errno $61$)

User requested termination in $f$.

(errno $71$)

User requested termination in $h$.

Notes

kalman_unscented_state applies the Unscented Kalman Filter (UKF), as described in Julier and Uhlmann (1997b) to a nonlinear state space model, with additive noise, which, at time $t$, can be described by:

$$\begin{array}{c}\begin{array}{cc}{x}_{t+1}& =F\left(x_t\right)+v_t\\ y_t& =H\left(x_t\right)+u_t\end{array}\end{array}$$

where $x_t$ represents the unobserved state vector of length $m_x$ and $y_t$ the observed measurement vector of length $m_y$. The process noise is denoted $v_t$, which is assumed to have mean zero and covariance structure ${\Sigma }_x$, and the measurement noise by $u_t$, which is assumed to have mean zero and covariance structure ${\Sigma }_y$.

Unscented Kalman Filter Algorithm

Given ${\hat{x}}_0$, an initial estimate of the state and $P_0$ and initial estimate of the state covariance matrix, the UKF can be described as follows:

1. Generate a set of sigma points (see Sigma Points):

   $$\begin{array}{cc}{X}_t& =[{\hat{x}}_{t-1}{\hat{x}}_{t-1}+\gamma \sqrt{P_{t-1}}{\hat{x}}_{t-1}-\gamma \sqrt{P_{t-1}}]\end{array}$$

2. Evaluate the known model function $F$:

   $$\begin{array}{cc}{F}_t& =F\left(X_t\right)\end{array}$$

   The function $F$ is assumed to accept the $m_x\times n$ matrix, $X_t$ and return an $m_x\times n$ matrix, $F_t$. The columns of both $X_t$ and $F_t$ correspond to different possible states. The notation $F_{t,i}$ is used to denote the $i$th column of $F_t$, hence the result of applying $F$ to the $i$th possible state.

3. Time Update:

   $$\begin{array}{c}\begin{array}{cc}{\hat{x}}_{\overline{t}}& ={\sum }_{i=1}^n{W}_i^m{F}_{t,i}\\ P_{\overline{t}}& ={\sum }_{i=1}^n{W}_i^c(F_{t,i}-{\hat{x}}_{\hat{t}}){(F_{t,i}-{\hat{x}}_{\hat{t}})}^T+{\Sigma }_x\end{array}\end{array}$$

4. Redraw another set of sigma points (see Sigma Points):

   $$\begin{array}{cc}{Y}_t& =[{\hat{x}}_{\hat{t}}{\hat{x}}_{\hat{t}}+\gamma \sqrt{P_{\overline{t}}}{\hat{x}}_{\hat{t}}-\gamma \sqrt{P_{\overline{t}}}]\end{array}$$

5. Evaluate the known model function $H$:

   $$\begin{array}{cc}{H}_t& =H\left(Y_t\right)\end{array}$$

   The function $H$ is assumed to accept the $m_x\times n$ matrix, $Y_t$ and return an $m_y\times n$ matrix, $H_t$. The columns of both $Y_t$ and $H_t$ correspond to different possible states. As above $H_{t,i}$ is used to denote the $i$th column of $H_t$.

6. Measurement Update:

   $$\begin{array}{c}\begin{array}{cc}{\hat{y}}_t& ={\sum }_{i=1}^n{W}_i^m{H}_{t,i}\\ P_{y{y}_t}& ={\sum }_{i=1}^n{W}_i^c(H_{t,i}-{\hat{y}}_t){(H_{t,i}-{\hat{y}}_t)}^T+{\Sigma }_y\\ P_{x{y}_t}& ={\sum }_{i=1}^n{W}_i^c(F_{t,i}-{\hat{x}}_{\hat{t}}){(H_{t,i}-{\hat{y}}_t)}^T\\ K_t& =P_{x{y}_t}{P}_{y{y}_t}^{-1}\\ {\hat{x}}_t& ={\hat{x}}_{\hat{t}}+K_t(y_t-{\hat{y}}_t)\\ P_t& =P_{\overline{t}}-K_t{P}_{y{y}_t}{K}_t^T\end{array}\end{array}$$

Here $K_t$ is the Kalman gain matrix, ${\hat{x}}_t$ is the estimated state vector at time $t$ and $P_t$ the corresponding covariance matrix. Rather than implementing the standard UKF as stated above kalman_unscented_state uses the square-root form described in the Haykin (2001).

Sigma Points

A nonlinear state space model involves propagating a vector of random variables through a nonlinear system and we are interested in what happens to the mean and covariance matrix of those variables. Rather than trying to directly propagate the mean and covariance matrix, the UKF uses a set of carefully chosen sample points, referred to as sigma points, and propagates these through the system of interest. An estimate of the propagated mean and covariance matrix is then obtained via the weighted sample mean and covariance matrix.

For a vector of $m$ random variables, $x$, with mean $\mu$ and covariance matrix $\Sigma$, the sigma points are usually constructed as:

$$X_t=[\mu \mu +\gamma \sqrt{\Sigma }\mu -\gamma \sqrt{\Sigma }]$$

When calculating the weighted sample mean and covariance matrix two sets of weights are required, one used when calculating the weighted sample mean, denoted $W^m$ and one used when calculating the weighted sample covariance matrix, denoted $W^c$. The weights and multiplier, $\gamma$, are constructed as follows:

$$\begin{array}{c}\begin{array}{cc}\lambda & ={\alpha }^2(L+\kappa )-L\\ \gamma & =\sqrt{L+\lambda }\\ W_i^m& =\{\begin{array}{cc}\frac{\lambda }{L+\lambda }& i=1\\ \frac{1}{2(L+\lambda )}& i=2,3,\dots ,2L+1\end{array}\\ W_i^c& =\{\begin{array}{cc}\frac{\lambda }{L+\lambda }+1-{\alpha }^2+\beta & i=1\\ \frac{1}{2(L+\lambda )}& i=2,3,\dots ,2L+1\end{array}\end{array}\end{array}$$

where, usually $L=m$ and $\alpha ,\beta$ and $\kappa$ are constants. The total number of sigma points, $n$, is given by $2L+1$. The constant $\alpha$ is usually set to somewhere in the range ${10}^{-4}\le \alpha \le 1$ and for a Gaussian distribution, the optimal values of $\kappa$ and $\beta$ are $3-L$ and $2$ respectively.

The constants, $\kappa$, $\alpha$ and $\beta$ are given by $\kappa =3-m_x$, $\alpha =1.0$ and $\beta =2$. If more control is required over the construction of the sigma points then the reverse communication function, kalman_unscented_state_revcom(), can be used instead.

References

Haykin, S, 2001, Kalman Filtering and Neural Networks, John Wiley and Sons

Julier, S J, 2002, The scaled unscented transformation, Proceedings of the 2002 American Control Conference (Volume 6), 4555–4559

Julier, S J and Uhlmann, J K, 1997, A consistent, debiased method for converting between polar and Cartesian coordinate systems, Proceedings of AeroSense ‘97, International Society for Optics and Phonotonics, 110–121

Julier, S J and Uhlmann, J K, 1997, A new extension of the Kalman Filter to nonlinear systems, International Symposium for Aerospace/Defense, Sensing, Simulation and Controls (Volume 3) (26)

See also

naginterfaces.library.examples.tsa.kalman_unscented_state_ex.main()