importanceSSM: Importance Sampling of Exponential Family State Space Model

In helske/KFAS: Kalman Filter and Smoother for Exponential Family State Space Models

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,
  expected = FALSE,
  H_tol = 1e+15
)

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.
expected	Logical value defining the approximation of H_t in case of Gamma and negative binomial distribution. Default is FALSE which matches the algorithm of Durbin & Koopman (1997), whereas TRUE uses the expected value of observations in the equations, leading to results which match with glm (where applicable). The latter case was the default behaviour of KFAS before version 1.3.8. Essentially this is the difference between observed and expected information.
H_tol	Tolerance parameter for check max(H) > H_tol, which suggests that the approximation converged to degenerate case with near zero signal-to-noise ratio. Default is very generous 1e15.

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

Value

A list containing elements

samples	Simulated samples.
weights	Importance weights.
model	Original model in case of save.model==TRUE.

Examples

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)

helske/KFAS documentation built on Sept. 9, 2023, 8:12 a.m.