naginterfaces.library.tsa.multi_​kalman_​sqrt_​invar

naginterfaces.library.tsa.multi_kalman_sqrt_invar(transf, a, b, stq, q, c, r, s, tol=0.0)

multi_kalman_sqrt_invar performs a combined measurement and time update of one iteration of the time-invariant Kalman filter using a square root covariance filter.

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

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

Parameters

transfstr, length 1

Indicates whether to transform the input matrix pair $(A,C)$ to lower observer Hessenberg form. The transformation will only be required on the first call to multi_kalman_sqrt_invar.

$transf=\text{'T'}$

The matrices in arrays $a$ and $c$ are transformed to lower observer Hessenberg form and the matrices in $b$ and $s$ are transformed as described in Notes.

$transf=\text{'H'}$

The matrices in arrays $a$, $c$ and $b$ should be as returned from a previous call to multi_kalman_sqrt_invar with $transf=\text{'T'}$.

afloat, array-like, shape $(n,n)$

If $transf=\text{'T'}$, the state transition matrix, $A$.

If $transf=\text{'H'}$, the transformed matrix as returned by a previous call to multi_kalman_sqrt_invar with $transf=\text{'T'}$.

bfloat, array-like, shape $(n,l)$

If $transf=\text{'T'}$, the noise coefficient matrix $B$.

If $transf=\text{'H'}$, the transformed matrix as returned by a previous call to multi_kalman_sqrt_invar with $transf=\text{'T'}$.

stqbool

If $stq=True$, the state noise covariance matrix $Q_i$ is assumed to be the identity matrix. Otherwise the lower triangular Cholesky factor, $Q_i^{1/2}$, must be provided in $q$.

qfloat, array-like, shape $(:,:)$

Note: the required extent for this argument in dimension 1 is determined as follows: if $not\left(stq\right)$: $l$; otherwise: $1$.

Note: the required extent for this argument in dimension 2 is determined as follows: if $not\left(stq\right)$: $l$; otherwise: $0$.

If $stq=False$, $q$ must contain the lower triangular Cholesky factor of the state noise covariance matrix, $Q_i^{1/2}$. Otherwise $q$ is not referenced.

cfloat, array-like, shape $(m,n)$

If $transf=\text{'T'}$, the measurement coefficient matrix, $C$.

If $transf=\text{'H'}$, the transformed matrix as returned by a previous call to multi_kalman_sqrt_invar with $transf=\text{'T'}$.

rfloat, array-like, shape $(m,m)$

The lower triangular Cholesky factor of the measurement noise covariance matrix $R_i^{1/2}$.

sfloat, array-like, shape $(n,n)$

If $transf=\text{'T'}$ the lower triangular Cholesky factor of the state covariance matrix, $S_i$.

If $transf=\text{'H'}$ the lower triangular Cholesky factor of the covariance matrix of the transformed state vector $S_i^{\ast }$ as returned from a previous call to multi_kalman_sqrt_invar with $transf=\text{'T'}$.

tolfloat, optional

The tolerance used to test for the singularity of $H_i^{1/2}$. If $0.0\le tol<m^2\times \text{machine precision}$, then $m^2\times \text{machine precision}$ is used instead. The inverse of the condition number of $H^{1/2}$ is estimated by a call to lapacklin.dtrcon. If this estimate is less than $tol$ then $H^{1/2}$ is assumed to be singular.

Returns

afloat, ndarray, shape $(n,n)$

If $transf=\text{'T'}$, the transformed matrix, $U^{\ast }A{U}^{\ast T}$, otherwise $a$ is unchanged.

bfloat, ndarray, shape $(n,l)$

If $transf=\text{'T'}$, the transformed matrix, $U^{\ast }B$, otherwise $b$ is unchanged.

cfloat, ndarray, shape $(m,n)$

If $transf=\text{'T'}$, the transformed matrix, $C{U}^{\ast T}$, otherwise $c$ is unchanged.

sfloat, ndarray, shape $(n,n)$

The lower triangular Cholesky factor of the transformed state covariance matrix, $S_{i+1}^{\ast }$.

kfloat, ndarray, shape $(n,m)$

The Kalman gain matrix for the transformed state vector premultiplied by the state transformed transition matrix, $U^{\ast }A{K}_i$.

hfloat, ndarray, shape $(m,m)$

The lower triangular matrix $H_i^{1/2}$.

ufloat, ndarray, shape $(n,:)$

If $transf=\text{'T'}$ the $n\times n$ transformation matrix $U^{\ast }$, otherwise $u$ is not referenced.

Raises

NagValueError

(errno $1$)

On entry, $transf=⟨value⟩$.

Constraint: $transf=\text{'T'}$ or $\text{'H'}$.

(errno $1$)

On entry, $tol=⟨value⟩$.

Constraint: $tol\ge 0.0$.

(errno $1$)

On entry, $l=⟨value⟩$.

Constraint: $l\ge 1$.

(errno $1$)

On entry, $m=⟨value⟩$.

Constraint: $m\ge 1$.

(errno $1$)

On entry, $n=⟨value⟩$.

Constraint: $n\ge 1$.

(errno $2$)

The matrix $H_i^{1/2}$ is singular.

Notes

In the NAG Library the traditional C interface for this routine uses a different algorithmic base. Please contact NAG if you have any questions about compatibility.

The Kalman filter arises from the state space model given by

$$\begin{array}{c}\begin{array}{cc}{X}_{i+1}=A{X}_i+B{W}_i\text{,}& Var\left(W_i\right)=Q_i\\ & \\ Y_i=C{X}_i+V_i\text{,}& Var\left(V_i\right)=R_i\end{array}\end{array}$$

where $X_i$ is the state vector of length $n$ at time $i$, $Y_i$ is the observation vector of length $m$ at time $i$ and $W_i$ of length $l$ and $V_i$ of length $m$ are the independent state noise and measurement noise respectively. The matrices $A,B$ and $C$ are time invariant.

The estimate of $X_i$ given observations $Y_1$ to $Y_{i-1}$ is denoted by ${\hat{X}}_{i|i-1}$ with state covariance matrix $Var\left({\hat{X}}_{i|i-1}\right)=P_{i|i-1}=S_i{S}_i^T$ while the estimate of $X_i$ given observations $Y_1$ to $Y_i$ is denoted by ${\hat{X}}_{i|i}$ with covariance matrix $Var\left({\hat{X}}_{i|i}\right)=P_{i|i}$. The update of the estimate, ${\hat{X}}_{i|i-1}$, from time $i$ to time $(i+1)$ is computed in two stages. First, the measurement-update is given by

$${\hat{X}}_{i|i}={\hat{X}}_{i|i-1}+K_i[Y_i-C{\hat{X}}_{i|i-1}]$$

where $K_i=P_{i|i}{C}^T{[C{P}_{i|i}{C}^T+R_i]}^{-1}$ is the Kalman gain matrix. The second stage is the time-update for $X$, which is given by

$${\hat{X}}_{i+1|i}=A{\hat{X}}_{i|i}+D_i{U}_i$$

where $D_i{U}_i$ represents any deterministic control used.

The square root covariance filter algorithm provides a stable method for computing the Kalman gain matrix and the state covariance matrix. The algorithm can be summarised as

$$\begin{array}{c}\left(\begin{array}{ccc}{R}_i^{1/2}& 0& C{S}_i\\ & & \\ 0& B{Q}_i^{1/2}& A{S}_i\end{array}\right)U=\left(\begin{array}{ccc}{H}_i^{1/2}& 0& 0\\ & & \\ G_i& S_{i+1}& 0\end{array}\right)\end{array}$$

where $U$ is an orthogonal transformation triangularizing the left-hand pre-array to produce the right-hand post-array. The triangularization is carried out via Householder transformations exploiting the zero pattern of the pre-array. The relationship between the Kalman gain matrix $K_i$ and $G_i$ is given by

$$A{K}_i=G_i{\left(H_i^{1/2}\right)}^{-1}\text{.}$$

In order to exploit the invariant parts of the model to simplify the computation of $U$ the results for the transformed state space $U^{\ast }X$ are computed where $U^{\ast }$ is the transformation that reduces the matrix pair $(A,C)$ to lower observer Hessenberg form. That is, the matrix $U^{\ast }$ is computed such that the compound matrix

$$\begin{array}{c}\left[\begin{array}{c}C{U}^{\ast T}\\ U^{\ast }A{U}^{\ast T}\end{array}\right]\end{array}$$

is a lower trapezoidal matrix. Further the matrix $B$ is transformed to $U^{\ast }B$. These transformations need only be computed once at the start of a series, and multi_kalman_sqrt_invar will, optionally, compute them. multi_kalman_sqrt_invar returns transformed matrices $U^{\ast }A{U}^{\ast T}$, $U^{\ast }B$, $C{U}^{\ast T}$ and $U^{\ast }A{K}_i$, the Cholesky factor of the updated transformed state covariance matrix $S_{i+1}^{\ast }$ (where $U^{\ast }{P}_{i+1|i}{U}^{\ast T}=S_{i+1}^{\ast }{S}_{i+1}^{\ast T}$) and the matrix $H_i^{1/2}$, valid for both transformed and original models, which is used in the computation of the likelihood for the model. Note that the covariance matrices $Q_i$ and $R_i$ can be time-varying.

References

Vanbegin, M, van Dooren, P and Verhaegen, M H G, 1989, Algorithm 675: FORTRAN subroutines for computing the square root covariance filter and square root information filter in dense or Hessenberg forms, ACM Trans. Math. Software (15), 243–256

Verhaegen, M H G and van Dooren, P, 1986, Numerical aspects of different Kalman filter implementations, IEEE Trans. Auto. Contr. (AC-31), 907–917