Skip to content

Commit

Permalink
updated documentation
Browse files Browse the repository at this point in the history
  • Loading branch information
helske committed May 10, 2017
1 parent ade067d commit 63af9d3
Show file tree
Hide file tree
Showing 25 changed files with 261 additions and 215 deletions.
8 changes: 2 additions & 6 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -36,11 +36,6 @@ S3method(logLik,ng_bsm)
S3method(logLik,ngssm)
S3method(logLik,nlg_ssm)
S3method(logLik,svm)
S3method(particle_simulate,bsm)
S3method(particle_simulate,gssm)
S3method(particle_simulate,ng_bsm)
S3method(particle_simulate,ngssm)
S3method(particle_simulate,svm)
S3method(particle_smoother,bsm)
S3method(particle_smoother,gssm)
S3method(particle_smoother,ng_bsm)
Expand Down Expand Up @@ -84,7 +79,6 @@ export(ng_bsm)
export(ngssm)
export(nlg_ssm)
export(normal)
export(particle_simulate)
export(particle_smoother)
export(run_mcmc)
export(sim_smoother)
Expand Down Expand Up @@ -125,4 +119,6 @@ importFrom(stats,time)
importFrom(stats,ts)
importFrom(stats,ts.union)
importFrom(stats,tsp)
importFrom(utils,head)
importFrom(utils,tail)
useDynLib(bssm)
164 changes: 82 additions & 82 deletions R/backwards_sampling.R
Original file line number Diff line number Diff line change
@@ -1,82 +1,82 @@
#' Backwards sampling
#'
#' Function \code{particle_simulate} performs backwards sampling from
#' the conditional distribution of the states $p(alpha|y, theta)$.
#'
#' @param object Model.
#' @param nsim Number of samples used in particle filtering.
#' @param nsim_store Number of samples to store.
#' @param seed Seed for RNG.
#' @param ... Ignored.
#' @rdname particle_simulate
#' @export
particle_simulate <- function(object, nsim, nsim_store, seed, ...) {
UseMethod("particle_simulate", object)
}
#' @method particle_simulate gssm
#' @rdname particle_simulate
#' @export
particle_simulate.gssm <- function(object, nsim, nsim_store = 1,
seed = sample(.Machine$integer.max, size = 1), ...) {

out <- gssm_backward_simulate(object, nsim, seed, nsim_store)

rownames(out$alpha) <- names(object$a1)
out$alpha <- aperm(out$alpha, c(2, 1, 3))
out
}

#' @method particle_simulate bsm
#' @rdname particle_simulate
#' @export
particle_simulate.bsm <- function(object, nsim, nsim_store = 1,
seed = sample(.Machine$integer.max, size = 1), ...) {

out <- bsm_backward_simulate(object, nsim, seed, nsim_store)

rownames(out$alpha) <- names(object$a1)
out$alpha <- aperm(out$alpha, c(2, 1, 3))
out
}

#' @method particle_simulate ngssm
#' @rdname particle_simulate
#' @export
particle_simulate.ngssm <- function(object, nsim, nsim_store = 1,
seed = sample(.Machine$integer.max, size = 1), ...) {

object$distribution <- pmatch(object$distribution,
c("poisson", "binomial", "negative binomial"))
out <- ngssm_backward_simulate(object, nsim, seed, nsim_store)

rownames(out$alpha) <- names(object$a1)
out$alpha <- aperm(out$alpha, c(2, 1, 3))
out
}
#' @method particle_simulate ng_bsm
#' @rdname particle_simulate
#' @export
particle_simulate.ng_bsm <- function(object, nsim, nsim_store = 1,
seed = sample(.Machine$integer.max, size = 1), ...) {

object$distribution <- pmatch(object$distribution,
c("poisson", "binomial", "negative binomial"))
out <- ng_bsm_backward_simulate(object, nsim, seed, nsim_store)

rownames(out$alpha) <- names(object$a1)
out$alpha <- aperm(out$alpha, c(2, 1, 3))
out
}

#' @method particle_simulate svm
#' @rdname particle_simulate
#' @export
particle_simulate.svm <- function(object, nsim, nsim_store = 1,
seed = sample(.Machine$integer.max, size = 1), ...) {

out <- svm_backward_simulate(object, nsim, seed, nsim_store)

rownames(out$alpha) <- names(object$a1)
out$alpha <- aperm(out$alpha, c(2, 1, 3))
out
}
#' #' Backwards sampling
#' #'
#' #' Function \code{particle_simulate} performs backwards sampling from
#' #' the conditional distribution of the states $p(alpha|y, theta)$.
#' #'
#' #' @param object Model.
#' #' @param nsim Number of samples used in particle filtering.
#' #' @param nsim_store Number of samples to store.
#' #' @param seed Seed for RNG.
#' #' @param ... Ignored.
#' #' @rdname particle_simulate
#' #' @export
#' particle_simulate <- function(object, nsim, nsim_store, seed, ...) {
#' UseMethod("particle_simulate", object)
#' }
#' #' @method particle_simulate gssm
#' #' @rdname particle_simulate
#' #' @export
#' particle_simulate.gssm <- function(object, nsim, nsim_store = 1,
#' seed = sample(.Machine$integer.max, size = 1), ...) {
#'
#' out <- gssm_backward_simulate(object, nsim, seed, nsim_store)
#'
#' rownames(out$alpha) <- names(object$a1)
#' out$alpha <- aperm(out$alpha, c(2, 1, 3))
#' out
#' }
#'
#' #' @method particle_simulate bsm
#' #' @rdname particle_simulate
#' #' @export
#' particle_simulate.bsm <- function(object, nsim, nsim_store = 1,
#' seed = sample(.Machine$integer.max, size = 1), ...) {
#'
#' out <- bsm_backward_simulate(object, nsim, seed, nsim_store)
#'
#' rownames(out$alpha) <- names(object$a1)
#' out$alpha <- aperm(out$alpha, c(2, 1, 3))
#' out
#' }
#'
#' #' @method particle_simulate ngssm
#' #' @rdname particle_simulate
#' #' @export
#' particle_simulate.ngssm <- function(object, nsim, nsim_store = 1,
#' seed = sample(.Machine$integer.max, size = 1), ...) {
#'
#' object$distribution <- pmatch(object$distribution,
#' c("poisson", "binomial", "negative binomial"))
#' out <- ngssm_backward_simulate(object, nsim, seed, nsim_store)
#'
#' rownames(out$alpha) <- names(object$a1)
#' out$alpha <- aperm(out$alpha, c(2, 1, 3))
#' out
#' }
#' #' @method particle_simulate ng_bsm
#' #' @rdname particle_simulate
#' #' @export
#' particle_simulate.ng_bsm <- function(object, nsim, nsim_store = 1,
#' seed = sample(.Machine$integer.max, size = 1), ...) {
#'
#' object$distribution <- pmatch(object$distribution,
#' c("poisson", "binomial", "negative binomial"))
#' out <- ng_bsm_backward_simulate(object, nsim, seed, nsim_store)
#'
#' rownames(out$alpha) <- names(object$a1)
#' out$alpha <- aperm(out$alpha, c(2, 1, 3))
#' out
#' }
#'
#' #' @method particle_simulate svm
#' #' @rdname particle_simulate
#' #' @export
#' particle_simulate.svm <- function(object, nsim, nsim_store = 1,
#' seed = sample(.Machine$integer.max, size = 1), ...) {
#'
#' out <- svm_backward_simulate(object, nsim, seed, nsim_store)
#'
#' rownames(out$alpha) <- names(object$a1)
#' out$alpha <- aperm(out$alpha, c(2, 1, 3))
#' out
#' }
7 changes: 5 additions & 2 deletions R/kfilter.R
Original file line number Diff line number Diff line change
Expand Up @@ -60,14 +60,16 @@ kfilter.svm <- function(object, ...) {
kfilter(gaussian_approx(object))
}

#' Extended Kalman Filtering
#' (Iterated) Extended Kalman Filtering
#'
#' Function \code{ekf} runs the extended Kalman filter for the given
#' Function \code{ekf} runs the (iterated) extended Kalman filter for the given
#' non-linear Gaussian model of class \code{nlg_ssm},
#' and returns the filtered estimates and one-step-ahead predictions of the
#' states \eqn{\alpha_t} given the data up to time \eqn{t}.
#'
#' @param object Model object
#' @param iekf_iter If \code{iekf_iter > 0}, iterated extended Kalman filter
#' is used with \code{iekf_iter} iterations.
#' @return List containing the log-likelihood,
#' one-step-ahead predictions \code{at} and filtered
#' estimates \code{att} of states, and the corresponding variances \code{Pt} and
Expand Down Expand Up @@ -95,6 +97,7 @@ ekf <- function(object, iekf_iter = 0) {
#' states \eqn{\alpha_t} given the data up to time \eqn{t}.
#'
#' @param object Model object
#' @param alpha,beta,kappa Tuning parameters for the UKF.
#' @return List containing the log-likelihood,
#' one-step-ahead predictions \code{at} and filtered
#' estimates \code{att} of states, and the corresponding variances \code{Pt} and
Expand Down
2 changes: 1 addition & 1 deletion R/loglik.R
Original file line number Diff line number Diff line change
Expand Up @@ -80,6 +80,6 @@ logLik.nlg_ssm <- function(object, nsim_states, method = "bsf", seed = 1,
object$R, object$Z_gn, object$T_gn, object$a1, object$P1,
object$theta, object$log_prior_pdf, object$known_params,
object$known_tv_params, object$n_states, object$n_etas,
as.integer(object$time_varying), as.integer(object$state_varying), nsim, seed,
as.integer(object$time_varying), as.integer(object$state_varying), nsim_states, seed,
max_iter, conv_tol, iekf_iter, pmatch(method, c("psi", "bsf", "apf", "ekf", "psi_df")))
}
13 changes: 11 additions & 2 deletions R/models.R
Original file line number Diff line number Diff line change
Expand Up @@ -1228,15 +1228,24 @@ mv_gssm <- function(y, Z, H, T, R, a1, P1, xreg = NULL, beta, state_names,
#' the user, as you must provide the several small C++ snippets which define the
#' model structure. See examples in ZZZ.
#' @param y Observations as multivariate time series (or matrix) of length \eqn{n}.
#' @param Z_fn,H_fn,T_fn,R_fn An external pointers for the C++ functions which
#' @param Z,H,T,R An external pointers for the C++ functions which
#' define the corresponding model functions.
#' @param Z_gn,T_gn An external pointers for the C++ functions which
#' define the gradients of the corresponding model functions.
#' @param a1 Prior mean for the initial state as a vector of length m.
#' @param P1 Prior covariance matrix for the initial state as m x m matrix.
#' @param theta Parameter vector passed to all model functions.
#' @param prior_pdf An external pointer for the C++ function which
#' @param known_params Vector of known parameters passed to all model functions.
#' @param known_tv_params Matrix of known parameters passed to all model functions.
#' @param n_states Number of states in the model.
#' @param n_etas Dimension of the noise term of the transition equation.
#' @param log_prior_pdf An external pointer for the C++ function which
#' computes the log-prior density given theta.
#' @param time_varying Optional logical vector of length 4, denoting whether the values of
#' Z, H, T, and R vary with respect to time variable (given identical states).
#' If used, can speed up some computations.
#' @param state_varying Optional logical vector of length 2.
#' Do H and T vary (linearly) with respect to states?
#' @param state_names Names for the states.
#' @return Object of class \code{gssm}.
#' @export
Expand Down
47 changes: 31 additions & 16 deletions R/particle_smoother.R
Original file line number Diff line number Diff line change
Expand Up @@ -5,13 +5,22 @@
#'
#' @param object Model.
#' @param nsim Number of samples.
#' @param filter Choice of particle filter algorithm. For Gaussian models,
#' possible choices are \code{"bsf"} (bootstrap particle filter) and
#' \code{"apf"} (auxiliary particle filter). In addition, for non-Gaussian or
#' non-linear models, \code{"psi"} uses psi-particle filter, and
#' for non-linear models options \code{"ekf"} (extended Kalman particle filter)
#' is also available.
#' @param smoothing_method Either \code{"fs"} (filter-smoother), or \code{"fbs"}
#' (forward-backward smoother).
#' @param filter_type For Gaussian models, \code{"bsf"} for bootstrap
#' filter. Also for non-Gaussian or non-linear models, \code{"psi"} uses
#' psi-particle filter, and for non-linear models options \code{"ekf"}
#' (extended Kalman particle filter) and \code{"apf"} (auxiliary particle filter)
#' are available.
#' @param max_iter Maximum number of iterations used in Gaussian approximation. Used psi-PF.
#' @param conv_tol Tolerance parameter used in Gaussian approximation. Used psi-PF.
#' @param iekf_iter If zero (default), first approximation for non-linear
#' Gaussian models is obtained from extended Kalman filter. If
#' \code{iekf_iter > 0}, iterated extended Kalman filter is used with
#' \code{iekf_iter} iterations.
#' @param optimal For Gaussian models, use optimal auxiliary particle filter?
#' Default is \code{TRUE}.
#' @param seed Seed for RNG.
#' @param ... Ignored.
#' @export
Expand All @@ -22,8 +31,9 @@ particle_smoother <- function(object, nsim, ...) {
#' @method particle_smoother gssm
#' @rdname particle_smoother
#' @export
particle_smoother.gssm <- function(object, nsim, smoothing_method = "fs",
filter_type = "bsf", seed = sample(.Machine$integer.max, size = 1),
particle_smoother.gssm <- function(object, nsim,
filter_type = "bsf", smoothing_method = "fs",
seed = sample(.Machine$integer.max, size = 1),
optimal = TRUE, ...) {

filter_type <- match.arg(filter_type, c("bsf", "apf"))
Expand All @@ -32,7 +42,7 @@ particle_smoother.gssm <- function(object, nsim, smoothing_method = "fs",
}
smoothing_method <- match.arg(smoothing_method, c("fs", "fbs"))
if(filter_type == "bsf") {
bsf_smoother(model, nsim, seed, TRUE, 1L)
bsf_smoother(object, nsim, seed, TRUE, 1L)
} else {
out <- aux_smoother(object, nsim, seed, TRUE, 1L, optimal)
}
Expand All @@ -47,8 +57,9 @@ particle_smoother.gssm <- function(object, nsim, smoothing_method = "fs",

#' @method particle_smoother bsm
#' @export
particle_smoother.bsm <- function(object, nsim, smoothing_method = "fs",
filter_type = "bsf", seed = sample(.Machine$integer.max, size = 1),
particle_smoother.bsm <- function(object, nsim,
filter_type = "bsf", smoothing_method = "fs",
seed = sample(.Machine$integer.max, size = 1),
optimal = TRUE, ...) {

filter_type <- match.arg(filter_type, c("bsf", "apf"))
Expand All @@ -72,8 +83,9 @@ particle_smoother.bsm <- function(object, nsim, smoothing_method = "fs",
#' @rdname particle_smoother
#' @method particle_smoother ngssm
#' @export
particle_smoother.ngssm <- function(object, nsim, smoothing_method = "fs",
filter_type = "bsf", seed = sample(.Machine$integer.max, size = 1),
particle_smoother.ngssm <- function(object, nsim,
filter_type = "bsf", smoothing_method = "fs",
seed = sample(.Machine$integer.max, size = 1),
max_iter = 100, conv_tol = 1e-8, ...) {

smoothing_method <- match.arg(smoothing_method, c("fs", "fbs"))
Expand Down Expand Up @@ -124,8 +136,9 @@ particle_smoother.ng_bsm <- function(object, nsim, filter_type = "psi",
}
#' @method particle_smoother svm
#' @export
particle_smoother.svm <- function(object, nsim, smoothing_method = "fs",
filter_type = "psi", seed = sample(.Machine$integer.max, size = 1),
particle_smoother.svm <- function(object, nsim,
filter_type = "psi", smoothing_method = "fs",
seed = sample(.Machine$integer.max, size = 1),
max_iter = 100, conv_tol = 1e-8, ...) {

smoothing_method <- match.arg(smoothing_method, c("fs", "fbs"))
Expand All @@ -150,8 +163,10 @@ particle_smoother.svm <- function(object, nsim, smoothing_method = "fs",

#' @method particle_smoother nlg_ssm
#' @export
particle_smoother.nlg_ssm <- function(object, nsim, smoothing_method = "fs",
filter_type = "psi", seed = sample(.Machine$integer.max, size = 1),
particle_smoother.nlg_ssm <- function(object, nsim,
filter_type = "psi",
smoothing_method = "fs",
seed = sample(.Machine$integer.max, size = 1),
max_iter = 100, conv_tol = 1e-8, iekf_iter = 0, ...) {

smoothing_method <- match.arg(smoothing_method, c("fs", "fbs"))
Expand Down
1 change: 1 addition & 0 deletions R/plot.predict_bssm.R
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
#' Plot predictions based on bssm package
#'
#' @importFrom ggplot2 autoplot ggplot geom_line geom_ribbon
#' @importFrom utils head tail
#' @method autoplot predict_bssm
#' @param object Object of class \code{predict_bssm}.
#' @param plot_mean Draw mean predictions. Default is \code{TRUE}.
Expand Down
1 change: 1 addition & 0 deletions R/predict.R
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@
#' samples of the model parameters to right place.
#' @param return_MCSE For Gaussian models, if \code{TRUE}, the Monte Carlo
#' standard errors of the intervals are also returned.
#' @param seed Seed for RNG.
#' @param ... Ignored.
#' @return List containing the mean predictions, quantiles and Monte Carlo
#' standard errors .
Expand Down
4 changes: 4 additions & 0 deletions R/print_mcmc.R
Original file line number Diff line number Diff line change
Expand Up @@ -82,6 +82,10 @@ print.mcmc_output <- function(x, ...) {
#' returns the expanded sample based on the counts.
#'
#' @param x Output from \code{\link{run_mcmc}}.
#' @param variable Expand parameters \code{"theta"} or states \code{"state"}.
#' @param times Vector of indices. In case of states, what time points to expand? Default is all.
#' @param states Vector of indices. In case of states, what states to expand? Default is all.
#' @param by_states If \code{TRUE} (default), return list by states. Otherwise by time.
#' @param ... Ignored.
#' @export
expand_sample <- function(x, variable = "theta", times, states, by_states = TRUE) {
Expand Down
Loading

0 comments on commit 63af9d3

Please sign in to comment.