astsa (version 2.1)

ffbs: Forward Filtering Backward Sampling

Description

FFBS algorithm for state space models

Usage

ffbs(y, A, mu0, Sigma0, Phi, sQ, sR, Ups = NULL, Gam = NULL, input = NULL)

Value

Xs

An array of sampled states

X0n

The sampled initial state (because R is 1-based)

Arguments

y

Data matrix, vector or time series.

A

Observation matrix. Can be constant or an array with dim=c(q,p,n) if time varying.

mu0

Initial state mean.

Sigma0

Initial state covariance matrix.

Phi

State transition matrix.

sQ

State error covariance matrix is Q = sQ%*%t(sQ) -- see details below. In the univariate case, it is the standard deviation.

sR

Observation error covariance matrix is R = sR%*%t(sR) -- see details below. In the univariate case, it is the standard deviation.

Ups

State input matrix.

Gam

Observation input matrix.

input

matrix or vector of inputs having the same row dimension as y.

Author

D.S. Stoffer

Details

For a linear state space model, the FFBS algorithm provides a way to sample a state sequence \(x_{0:n}\) from the posterior \(\pi(x_{0:n} \mid \Theta, y_{1:n})\) with parameters \(\Theta\) and data \(y_{1:n}\).

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.

\(x_t\) is p-dimensional, \(y_t\) is q-dimensional, and \(u_t\) is r-dimensional. Note that \(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.

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/.

Examples

if (FALSE) {

## -- this is just one pass --##
# 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 one pass of FFBS, plot data, states and sampled states  
run = ffbs(y, A=1, mu0=0, Sigma0=10, Phi=1, sQ=1, sR=3)
tsplot(cbind(y,run$Xs), spaghetti=TRUE, type='o', col=c(8,4), pch=c(1,NA))
legend('topleft', legend=c("y(t)","xs(t)"), lty=1, col=c(8,4), bty="n", pch=c(1,NA))
}