Simulation of a Gaussian State Space Model
Function simulateSMM
simulates states, signals, disturbances or missing observations of
the Gaussian state space model either conditional on the data (simulation smoother) or
unconditionally.
simulateSSM( object, type = c("states", "signals", "disturbances", "observations", "epsilon", "eta"), filtered = FALSE, nsim = 1, antithetics = FALSE, conditional = TRUE )
object |
Gaussian state space object of class |
type |
What to simulate. |
filtered |
Simulate from p(α[t]|y[t-1],...,y[1]) instead of p(α|y). |
nsim |
Number of independent samples. Default is 1. |
antithetics |
Use antithetic variables in simulation. Default is |
conditional |
Simulations are conditional to data. If |
Simulation smoother algorithm is based on article by J. Durbin and S.J. Koopman (2002).
The simulation filter (filtered = TRUE
) is a straightforward modification
of the simulations smoother, where only filtering steps are performed.
Function can use two antithetic variables, one for location and other for scale, so output contains four blocks of simulated values which correlate which each other (ith block correlates negatively with (i+1)th block, and positively with (i+2)th block etc.).
Note that KFAS versions 1.2.0 and older, for unconditional simulation the initial
distribution of states was fixed so that a1
was set to the smoothed estimates
of the first state and the initial variance was set to zero. Now original
a1
and P1
are used, and P1inf
is ignored (i.e. diffuse states are
fixed to corresponding elements of a1
).
An n x k x nsim array containing the simulated series, where k is number of observations, signals, states or disturbances.
Durbin J. and Koopman, S.J. (2002). A simple and efficient simulation smoother for state space time series analysis, Biometrika, Volume 89, Issue 3
set.seed(123) # simulate new observations from the "fitted" model model <- SSModel(Nile ~ SSMtrend(1, Q = 1469), H = 15099) # signal conditional on the data i.e. samples from p(theta | y) # unconditional simulation is not reasonable as the model is nonstationary signal_sim <- simulateSSM(model, type = "signals", nsim = 10) # and add unconditional noise term i.e samples from p(epsilon) epsilon_sim <- simulateSSM(model, type = "epsilon", nsim = 10, conditional = FALSE) observation_sim <- signal_sim + epsilon_sim ts.plot(observation_sim[,1,], Nile, col = c(rep(2, 10), 1), lty = c(rep(2, 10), 1), lwd = c(rep(1, 10), 2)) # fully unconditional simulation: observation_sim2 <- simulateSSM(model, type = "observations", nsim = 10, conditional = FALSE) ts.plot(observation_sim[,1,], observation_sim2[,1,], Nile, col = c(rep(2:3, each = 10), 1), lty = c(rep(2, 20), 1), lwd = c(rep(1, 20), 2)) # illustrating use of antithetics model <- SSModel(matrix(NA, 100, 1) ~ SSMtrend(1, 1, P1inf = 0), H = 1) set.seed(123) sim <- simulateSSM(model, "obs", nsim = 2, antithetics = TRUE) # first time points sim[1,,] # correlation structure between simulations with two antithetics cor(sim[,1,]) out_NA <- KFS(model, filtering = "none", smoothing = "state") model["y"] <- sim[, 1, 1] out_obs <- KFS(model, filtering = "none", smoothing = "state") set.seed(40216) # simulate states from the p(alpha | y) sim_conditional <- simulateSSM(model, nsim = 10, antithetics = TRUE) # mean of the simulated states is exactly correct due to antithetic variables mean(sim_conditional[2, 1, ]) out_obs$alpha[2] # for variances more simulations are needed var(sim_conditional[2, 1, ]) out_obs$V[2] set.seed(40216) # no data, simulations from p(alpha) sim_unconditional <- simulateSSM(model, nsim = 10, antithetics = TRUE, conditional = FALSE) mean(sim_unconditional[2, 1, ]) out_NA$alpha[2] var(sim_unconditional[2, 1, ]) out_NA$V[2] ts.plot(cbind(sim_conditional[,1,1:5], sim_unconditional[,1,1:5]), col = rep(c(2,4), each = 5)) lines(out_obs$alpha, lwd=2)
Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.