Ksmooth: Quick Kalman Smoother

In astsa: Applied Statistical Time Series Analysis

Quick Kalman Smoother

Description

Returns the smoother values for various linear state space models. The predicted and filtered values and the likelihood at the given parameter values are also returned (via Kfilter). This script replaces Ksmooth0, Ksmooth1, and Ksmooth2.

Usage

Ksmooth(y, A, mu0, Sigma0, Phi, sQ, sR, Ups = NULL, Gam = NULL, 
         input = NULL, S = NULL, version = 1)

Arguments

y	data matrix (n x q), vector or time series, n = number of observations. Use NA or zero (0) for missing data.
A	can be constant or an array with dimension dim=c(q,p,n) if time varying (see details). Use NA or zero (0) for missing data.
mu0	initial state mean vector (p x 1)
Sigma0	initial state covariance matrix (p x p)
Phi	state transition matrix (p x p)
sQ	state error pre-matrix (see details)
sR	observation error pre-matrix (see details)
Ups	state input matrix (p x r); leave as NULL (default) if not needed
Gam	observation input matrix (q x r); leave as NULL (default) if not needed
input	NULL (default) if not needed or a matrix (n x r) of inputs having the same row dimension (n) as y
S	covariance matrix between state and observation errors; not necessary to specify if not needed and only used if version=2; see details
version	either 1 (default) or 2; version 2 allows for correlated errors

Details

This script replaces Ksmooth0, Ksmooth1, and Ksmooth2 by combining all cases. The major difference is how to specify the covariance matrices; in particular, sQ = t(cQ) and sR = t(cR) where cQ and cR were used in Kfilter0-1-2 scripts.

The states x_t are p-dimensional, the data y_t are q-dimensional, and the inputs u_t are r-dimensional for t=1, \dots, n. The initial state is x_0 \sim N(\mu_0, \Sigma_0).

The measurement matrices A_t can be constant or time varying. If time varying, they should be entered as an array of dimension dim = c(q,p,n). Otherwise, just enter the constant value making sure it has the appropriate q \times p dimension.

Version 1 (default): The general model is

x_t = \Phi x_{t-1} + \Upsilon u_{t} + sQ\, w_t \quad w_t \sim iid\ N(0,I)

y_t = A_t x_{t-1} + \Gamma u_{t} + sR\, v_t \quad v_t \sim iid\ N(0,I)

where w_t \perp v_t. Consequently the state noise covariance matrix is Q = sQ\, sQ' and the observation noise covariance matrix is R = sR\, sR' and sQ, sR do not have to be square as long as everything is conformable. Notice the specification of the state and observation covariances has changed from the original scripts.

NOTE: If it is easier to model in terms of Q and R, simply input the square root matrices sQ = Q %^% .5 and sR = R %^% .5.

Version 2 (correlated errors): The general model is

x_{t+1} = \Phi x_{t} + \Upsilon u_{t+1} + sQ\, w_t \quad w_t \sim iid\ N(0,I)

y_t = A_t x_{t-1} + \Gamma u_{t} + sR\, v_t \quad v_t \sim iid\ N(0,I)

where S = {\rm Cov}(w_t, v_t), and NOT {\rm Cov}(sQ\, w_t, sR\, v_t).

NOTE: If it is easier to model in terms of Q and R, simply input the square root matrices sQ = Q %^% .5 and sR = R %^% .5.

Note that in either version, sQ\, w_t has to be p-dimensional, but w_t does not, and sR\, v_t has to be q-dimensional, but v_t does not.

Value

Time varying values are returned as arrays.

Xs	state smoothers
Ps	smoother mean square error
X0n	initial mean smoother
P0n	initial smoother covariance
J0	initial value of the J matrix
J	the J matrices
Xp	state predictors
Pp	mean square prediction error
Xf	state filters
Pf	mean square filter error
like	negative of the log likelihood
innov	innovation series
sig	innovation covariances
Kn	the value of the last Gain

Note

Note that Ksmooth is similar to Ksmooth-0-1-2 except that only the essential values need to be entered (and come first in the statement); the optional values such as input are set to NULL by default if they are not needed. This version is faster than the older versions. The biggest change was to how the covarainces are specified. For example, if you have code that used Ksmooth1, just use sQ = t(cQ) and sR = t(cR) here.

Author(s)

D.S. Stoffer

References

You can find demonstrations of astsa capabilities at FUN WITH ASTSA.

The most recent version of the package can be found at https://github.com/nickpoison/astsa/.

In addition, the News and ChangeLog files are at https://github.com/nickpoison/astsa/blob/master/NEWS.md.

The webpages for the texts and some help on using R for time series analysis can be found at https://nickpoison.github.io/.

See Also

Kfilter

Examples

# generate some data
 set.seed(1)
 sQ  = 1; sR = 3; n = 100  
 mu0 = 0; Sigma0 = 10; x0 = rnorm(1,mu0,Sigma0)
 w = rnorm(n); v = rnorm(n)
 x = c(x0 + sQ*w[1]);  y = c(x[1] + sR*v[1])   # initialize
for (t in 2:n){
  x[t] = x[t-1] + sQ*w[t]
  y[t] = x[t] + sR*v[t]   
  }
# run and plot the smoother  
run = Ksmooth(y, A=1, mu0, Sigma0, Phi=1, sQ, sR)
tsplot(cbind(y,run$Xs), spaghetti=TRUE, type='o', col=c(4,6), pch=c(1,NA), margins=1)
# CRAN tests need extra white space :( so margins=1 above is not necessary otherwise
legend('topleft', legend=c("y(t)","Xs(t)"), lty=1, col=c(4,6), bty="n", pch=c(1,NA))

astsa documentation built on May 29, 2024, 10:29 a.m.