fks.SP {FKF.SP}

Fast Kalman Smoother through Sequential Processing

Description

The Kalman smoother is a backwards algorithm that is run after the Kalman filter that allows the user to refine estimates of previous states to produce "smoothed" estimates of state variables. This function performs the "Kalman smoother" algorithm using sequential processing, an approach that can substantially improve processing time over the traditional Kalman filtering/smoothing algorithms. The primary application of Kalman smoothing is in conjunction with expectation-maximization to estimate the parameters of a state space model. This function is called after running fkf.SP. fks.SP wraps the C-function fks_SP which relies upon the linear algebra subroutines of BLAS (Basic Linear Algebra Subprograms).

Usage

fks.SP(FKF.SP_obj)

Arguments

FKF.SP_obj

An S3-object of class "fkf.SP", returned by fkf.SP when using the argument verbose = TRUE.

Details

fks.SP is typically called after the fkf.SP function to calculate "smoothed" estimates of state variables and their corresponding variances. Smoothed estimates are used when utilizing expectation-maximization (EM) to efficiently estimate the parameters of a state space model.

Sequential Processing Kalman smoother solution:

The fks.SP function uses the solution to the Kalman smoother through sequential processing provided in the textbook of Durbin and Koopman (2001).

Given a state space model has been filtered through the sequential processing Kalman filter algorithm described in fkf.SP, the smoother can be reformulated for the univariate series:

\[y_t'=(y_{(1,1)},y_{(1,2)},\cdots,y_{(1,p_{1})},y_{(2,1)},\cdots,y_{(t,p_{t})})\]

The sequential processing Kalman smoother approach iterates backwards through both observations and time, i.e.: \(i=p_{t}, \cdots, 1\) and \(t=n,\cdots,1\), where \(p_{t}\) is the number of observations at time \(t\) and \(n\) is the total number of observations.

The initialisations are:

\[r_{(n,p_{n})} = 0\] \[N_{(n,p_{n})}=0\]

Then, \(r\) and \(N\) are recursively calculated through:

\[L_{t,i} = I_{m} - K_{t,i} Z_{t,i}\] \[r_{(t,i-1)} = Z_{t,i}' F_{t,i}^{-1} v_{t,i} + L_{t,i}' r_{t,i}\] \[N_{t,i-1} = Z_{t,i}' F_{t,i}^{-1} Z_{t,i} + L_{t,i}' N_{t,i} L_{t,i}\] \[r_{t-1,p_{t}} = T_{t-1}' r_{t,0}\] \[N_{t-1,p_{t}} = T_{t-1}' N_{t,0} T_{t-1}\]

for \(i=p_{t},\cdots,1\) and \(t=n,\cdots,1\)

The equations for \(r_{t-1},p_{t}\) and \(N_{t-1,p_{t}}\) do not apply for \(t=1\)

Under this formulation, the values for \(r_{t,0}\) and \(N_{t,0}\) are the same as the values for the smoothing quantities of \(r_{t-1}\) and \(N_{t-1}\) of the standard smoothing equations, respectively.

The standard smoothing equations for \(\hat{a_{t}}\) and \(V_t\) are used:

\[\hat{a_{t}} = a_{t} + P_{t} r_{t-1}\] \[V_t = P_t - P_t N_{t-1} P_t\]

Where:

\[a_{t}=a_{t,1}\] \[P_{t} = P_{t,1}\]

In the equations above, \(r_{t,i}\) is an \(m \times 1\) vector, \(I_{m}\) is an \(m \times m\) identity matrix, \(K_{t,i}\) is an \(m \times 1\) column vector, \(Z_{t,i}\) is a \(1 \times m\) row vector, and both \(F_{t,i}^{-1}\) and \(v_{t,i}\) are scalars. The reduced dimensionality of many of the variables in this formulation compared to traditional Kalman smoothing can result in increased computational efficiency.

Finally, in the formulation described above, \(a_{t}\) and \(P_{t}\) correspond to the values of att and ptt returned from the fkf.SP function, respectively.

Value

An S3-object of class fks.SP, which is a list with the following elements:

ahatt	A m * n-matrix containing the smoothed state variables, i.e. ahatt[,t] = \(a_{t|n}\)
Vt	A m * m * n-array containing the variances of ahatt, i.e. Vt[,,t] = \(P_{t|n}\)

References

Aspinall, T. W., Harris, G., Gepp, A., Kelly, S., Southam, C., and Vanstone, B. (2022). The Estimation of Commodity Pricing Models with Applications in Capital Investments. Available Online.

Durbin, James, and Siem Jan Koopman (2001). Time series analysis by state space methods. Oxford university press.

Examples

### Perform Kalman Filtering and Smoothing through sequential processing:
#Nile's annual flow:
yt <- Nile

# Incomplete Nile Data - two NA's are present:
yt[c(3, 10)] <- NA

dt <- ct <- matrix(0)
Zt <- Tt <- matrix(1)
a0 <- yt[1]   # Estimation of the first year flow
P0 <- matrix(100)       # Variance of 'a0'

# Parameter estimation - maximum likelihood estimation:
# Unknown parameters initial estimates:
GGt <- HHt <- var(yt, na.rm = TRUE) * .5
HHt = matrix(HHt)
GGt = matrix(GGt)
yt = rbind(yt)
# Filter through the Kalman filter - sequential processing:
Nile_filtered <- fkf.SP(HHt = matrix(HHt), GGt = matrix(GGt), a0 = a0, P0 = P0, dt = dt, ct = ct,
                  Zt = Zt, Tt = Tt, yt = rbind(yt), verbose = TRUE)
# Smooth filtered values through the Kalman smoother - sequential processing:
Smoothed_Estimates <- fks.SP(Nile_filtered)