ssm_sde: Univariate state space model with continuous SDE dynamics

In bssm: Bayesian Inference of Non-Linear and Non-Gaussian State Space Models

Univariate state space model with continuous SDE dynamics

Description

Constructs an object of class ssm_sde by defining the functions for the drift, diffusion and derivative of diffusion terms of univariate SDE, as well as the log-density of observation equation. We assume that the observations are measured at integer times (missing values are allowed).

Usage

ssm_sde(
  y,
  drift,
  diffusion,
  ddiffusion,
  obs_pdf,
  prior_pdf,
  theta,
  x0,
  positive
)

Arguments

y	Observations as univariate time series (or vector) of length n.
drift, diffusion, ddiffusion	An external pointers for the C++ functions which define the drift, diffusion and derivative of diffusion functions of SDE.
obs_pdf	An external pointer for the C++ function which computes the observational log-density given the the states and parameter vector theta.
prior_pdf	An external pointer for the C++ function which computes the prior log-density given the parameter vector theta.
theta	Parameter vector passed to all model functions.
x0	Fixed initial value for SDE at time 0.
positive	If TRUE, positivity constraint is forced by abs in Milstein scheme.

Details

As in case of ssm_nlg models, these general models need a bit more effort from the user, as you must provide the several small C++ snippets which define the model structure. See vignettes for an example and cpp_example_model.

Value

An object of class ssm_sde.

Examples

 # Takes a while on CRAN
library("sde")
set.seed(1)
# theta_0 = rho = 0.5
# theta_1 = nu = 2
# theta_2 = sigma = 0.3
x <- sde.sim(t0 = 0, T = 50, X0 = 1, N = 50,
       drift = expression(0.5 * (2 - x)),
       sigma = expression(0.3),
       sigma.x = expression(0))
y <- rpois(50, exp(x[-1]))

# source c++ snippets
pntrs <- cpp_example_model("sde_poisson_OU")

sde_model <- ssm_sde(y, pntrs$drift, pntrs$diffusion,
 pntrs$ddiffusion, pntrs$obs_density, pntrs$prior,
 c(rho = 0.5, nu = 2, sigma = 0.3), 1, positive = FALSE)

est <- particle_smoother(sde_model, L = 12, particles = 500)

ts.plot(cbind(x, est$alphahat, 
  est$alphahat - 2*sqrt(c(est$Vt)), 
  est$alphahat + 2*sqrt(c(est$Vt))), 
  col = c(2, 1, 1, 1), lty = c(1, 1, 2, 2))


# Takes time with finer mesh, parallelization with IS-MCMC helps a lot
out <- run_mcmc(sde_model, L_c = 4, L_f = 8, 
  particles = 50, iter = 2e4,
  threads = 4L)

bssm documentation built on Nov. 2, 2023, 6:25 p.m.