Learn R Programming

KFAS (version 1.2.5)

importanceSSM: Importance Sampling of Exponential Family State Space Model

Description

Function importanceSSM simulates states or signals of the exponential family state space model conditioned with the observations, returning the simulated samples of the states/signals with the corresponding importance weights.

Usage

importanceSSM(model, type = c("states", "signals"), filtered = FALSE, nsim = 1000, save.model = FALSE, theta, antithetics = FALSE, maxiter = 50)

Arguments

model
Exponential family state space model of class SSModel.
type
What to simulate, "states" or "signals". Default is "states"
filtered
Simulate from $p(\alpha_t|y_{t-1},...,y_1)$ instead of $p(\alpha|y)$. Note that for large models this can be very slow. Default is FALSE.
nsim
Number of independent samples. Default is 1000.
save.model
Return the original model with the samples. Default is FALSE.
theta
Initial values for the conditional mode theta.
antithetics
Logical. If TRUE, two antithetic variables are used in simulations, one for location and another for scale. Default is FALSE.
maxiter
Maximum number of iterations used in linearisation. Default is 50.

Value

A list containing elements A list containing elements

Details

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

Examples

Run this code
data("sexratio")
model <- SSModel(Male ~ SSMtrend(1, Q = list(NA)), u = sexratio[,"Total"], data = sexratio,
                distribution = "binomial")
fit <- fitSSM(model, inits = -15, method = "BFGS")
fit$model$Q #1.107652e-06
# Computing confidence intervals for sex ratio
# Uses importance sampling on response scale (1000 samples with antithetics)
set.seed(1)
imp <- importanceSSM(fit$model, nsim = 250, antithetics = TRUE)
sexratio.smooth <- numeric(length(model$y))
sexratio.ci <- matrix(0, length(model$y), 2)
w <- imp$w/sum(imp$w)
for(i in 1:length(model$y)){
  sexr <- exp(imp$sample[i,1,])
  sexratio.smooth[i]<-sum(sexr*w)
  oo <- order(sexr)
  sexratio.ci[i,] <- c(sexr[oo][which.min(abs(cumsum(w[oo]) - 0.05))],
                   sexr[oo][which.min(abs(cumsum(w[oo]) - 0.95))])
}

## Not run: 
# # Filtered estimates
# impf <- importanceSSM(fit$model, nsim = 250, antithetics = TRUE,filtered=TRUE)
# sexratio.filter <- rep(NA,length(model$y))
# sexratio.fci <- matrix(NA, length(model$y), 2)
# w <- impf$w/rowSums(impf$w)
# for(i in 2:length(model$y)){
#   sexr <- exp(impf$sample[i,1,])
#   sexratio.filter[i] <- sum(sexr*w[i,])
#   oo<-order(sexr)
#   sexratio.fci[i,] <- c(sexr[oo][which.min(abs(cumsum(w[i,oo]) - 0.05))],
#                     sexr[oo][which.min(abs(cumsum(w[i,oo]) - 0.95))])
# }
# 
# ts.plot(cbind(sexratio.smooth,sexratio.ci,sexratio.filter,sexratio.fci),
#         col=c(1,1,1,2,2,2),lty=c(1,2,2,1,2,2))
# ## End(Not run)

Run the code above in your browser using DataLab