Learn R Programming

bssm (version 2.0.2)

ar1_ng: Non-Gaussian model with AR(1) latent process

Description

Constructs a simple non-Gaussian model where the state dynamics follow an AR(1) process.

Usage

ar1_ng(y, rho, sigma, mu, distribution, phi, u, beta, xreg = NULL)

Value

An object of class ar1_ng.

Arguments

y

A vector or a ts object of observations.

rho

A prior for autoregressive coefficient. Should be an object of class bssm_prior.

sigma

A prior for the standard deviation of noise of the AR-process. Should be an object of class bssm_prior

mu

A fixed value or a prior for the stationary mean of the latent AR(1) process. Should be an object of class bssm_prior or scalar value defining a fixed mean such as 0.

distribution

Distribution of the observed time series. Possible choices are "poisson", "binomial", "gamma", and "negative binomial".

phi

Additional parameter relating to the non-Gaussian distribution. For negative binomial distribution this is the dispersion term, for gamma distribution this is the shape parameter, and for other distributions this is ignored. Should an object of class bssm_prior or a positive scalar.

u

A vector of positive constants for non-Gaussian models. For Poisson, gamma, and negative binomial distribution, this corresponds to the offset term. For binomial, this is the number of trials.

beta

A prior for the regression coefficients. Should be an object of class bssm_prior or bssm_prior_list (in case of multiple coefficients) or missing in case of no covariates.

xreg

A matrix containing covariates with number of rows matching the length of y. Can also be ts, mts or similar object convertible to matrix.

Examples

Run this code
model <- ar1_ng(discoveries, rho = uniform(0.5,-1,1), 
  sigma = halfnormal(0.1, 1), mu = normal(0, 0, 1), 
  distribution = "poisson")
out <- run_mcmc(model, iter = 1e4, mcmc_type = "approx",
  output_type = "summary")
  
ts.plot(cbind(discoveries, exp(out$alphahat)), col = 1:2)

set.seed(1)
n <- 30
phi <- 2
rho <- 0.9
sigma <- 0.1
beta <- 0.5
u <- rexp(n, 0.1)
x <- rnorm(n)
z <- y <- numeric(n)
z[1] <- rnorm(1, 0, sigma / sqrt(1 - rho^2))
y[1] <- rnbinom(1, mu = u * exp(beta * x[1] + z[1]), size = phi)
for(i in 2:n) {
  z[i] <- rnorm(1, rho * z[i - 1], sigma)
  y[i] <- rnbinom(1, mu = u * exp(beta * x[i] + z[i]), size = phi)
}

model <- ar1_ng(y, rho = uniform_prior(0.9, 0, 1), 
  sigma = gamma_prior(0.1, 2, 10), mu = 0., 
  phi = gamma_prior(2, 2, 1), distribution = "negative binomial",
  xreg = x, beta = normal_prior(0.5, 0, 1), u = u)

Run the code above in your browser using DataLab