RZigZag

Description

Implements various piecewise deterministic Monte Carlo methods, including the Zig-Zag Sampler (Bierkens, Fearnhead, Roberts, 2019, https://arxiv.org/abs/1607.03188) and the Bouncy Particle Sampler (BPS, Bouchard-Côté et al., 2017, https://arxiv.org/abs/1510.02451). Typical usage consists of first creating a "skeleton" consisting of "events", which can be used directly for plotting trajectories. The skeleton may be post-processed to extract information, such as as moment and covariance estimates, discrete samples at fixed time intervals along the trajectory, effective sample size and asymptotic variance.

Details

This package currently consists of the following functions for generating skeletons: ZigZagLogistic for logistic regression, ZigZagGaussian for multivariate Gaussian, ZigZagIIDGaussian for a IID Gaussian target using Zig-Zag, ZigZagStudentT for spherically symmetric or factorized Student-t distribution, BPSGaussian for multivariate Gaussian using BPS, BPSIIDGaussian for a IID Gaussian target using BPS, BPSStudentT for BPS applied to a spherically symmetric or factorized Student-t distribution. Furthermore the package contains the following functions for post-processing: EstimateESS (to estimate asymptotic variance and effective sample size for individual coordinates), EstimateMoment, EstimateCovarianceMatrix and DiscreteSamples.

Author(s)

Joris Bierkens

With thanks to Matt Moores, https://mattstats.wordpress.com/, for his help in getting from C++ code to a CRAN-ready Rcpp based package.

BPSGaussian

Description

Applies the BPS Sampler to a Gaussian target distribution, as detailed in Bouchard-Côté et al, 2017. Assume potential of the form

U(x) = (x - mu)^T V (x - mu)/2,

i.e. a Gaussian with mean vector mu and covariance matrix inv(V)

Usage

BPSGaussian(V, mu, n_iter = -1L, finalTime = -1, x0 = numeric(0),
  v0 = numeric(0), refresh_rate = 1, unit_velocity = FALSE)

Arguments

Value

Returns a list with the following objects:

Times: Vector of switching times

Positions: Matrix whose columns are locations of switches. The number of columns is identical to the length of skeletonTimes. Be aware that the skeleton points themselves are NOT samples from the target distribution.

Velocities: Matrix whose columns are velocities just after switches. The number of columns is identical to the length of skeletonTimes.

Examples

V <- matrix(c(3,1,1,3),nrow=2)
mu <- c(2,2)
x0 <- c(0,0)
result <- BPSGaussian(V, mu, n_iter = 100, x0 = x0)
plot(result$Positions[1,], result$Positions[2,],type='l',asp=1)

BPSIIDGaussian

Description

Applies the Bouncy Particle Sampler to a IID Gaussian distribution

Usage

BPSIIDGaussian(variance, dim = 1L, n_iter = -1L, finalTime = -1,
  x0 = numeric(0), v0 = numeric(0), refresh_rate = 1,
  unit_velocity = FALSE)

Arguments

Value

Returns a list with the following objects:

Times: Vector of switching times

Positions: Matrix whose columns are locations of switches. The number of columns is identical to the length of skeletonTimes. Be aware that the skeleton points themselves are NOT samples from the target distribution.

Velocities: Matrix whose columns are velocities just after switches. The number of columns is identical to the length of skeletonTimes.

Examples

result <- BPSIIDGaussian(1, 2, 1000)
plot(result$Positions[2,], result$Positions[1,],type='l',asp=1)

BPSStudentT

Description

Applies the Zig-Zag Sampler to a Student T distribution (IID or spherically symmetric)

Usage

BPSStudentT(dof, dim = 1L, n_iter = -1L, finalTime = -1,
  x0 = numeric(0), v0 = numeric(0), sphericallySymmetric = TRUE,
  refresh_rate = 1, unit_velocity = FALSE)

Arguments

Value

Returns a list with the following objects:

Times: Vector of switching times

Positions: Matrix whose columns are locations of switches. The number of columns is identical to the length of skeletonTimes. Be aware that the skeleton points themselves are NOT samples from the target distribution.

Velocities: Matrix whose columns are velocities just after switches. The number of columns is identical to the length of skeletonTimes.

Examples

dim = 2
dof = 4
result <- BPSStudentT(dof, dim, n_iter=1000,sphericallySymmetric = TRUE)
plot(result$Positions[1,], result$Positions[2,],type='l',asp=1)

DiscreteSamples

Description

Extract discrete samples from a skeleton

Usage

DiscreteSamples(skeletonList, n_samples, coordinate = -1L)

Arguments

Value

Returns a list containing the extracted samples and the times (on the continuous time scale) at which the samples are extracted

EstimateCovarianceMatrix

Description

Estimates the covariance matrix of a piecewise deterministic skeleton

Usage

EstimateCovarianceMatrix(skeletonList, coordinate = -1L,
  zeroMeans = FALSE)

Arguments

Value

Returns a list containing the estimated moment

EstimateESS

Description

Estimates the effective sample size (ESS) of a piecewise deterministic skeleton

Usage

EstimateESS(skeletonList, n_batches = 100L, coordinate = -1L,
  zeroMeans = FALSE)

Arguments

Value

Returns a list containing the estimated asymptotic variance, ESS and estimated covariance matrix

EstimateMoment

Description

Estimates the p-th moment of a piecewise deterministic skeleton

Usage

EstimateMoment(skeletonList, p, coordinate = -1L)

Arguments

Value

Returns a list containing the estimated moment

ZigZagGaussian

Description

Applies the Zig-Zag Sampler to a Gaussian target distribution, as detailed in Bierkens, Fearnhead, Roberts, The Zig-Zag Process and Super-Efficient Sampling for Bayesian Analysis of Big Data, 2016. Assume potential of the form

U(x) = (x - mu)^T V (x - mu)/2,

i.e. a Gaussian with mean vector mu and covariance matrix inv(V)

Usage

ZigZagGaussian(V, mu, n_iter = -1L, finalTime = -1, x0 = numeric(0),
  v0 = numeric(0))

Arguments

Value

Returns a list with the following objects:

Times: Vector of switching times

Positions: Matrix whose columns are locations of switches. The number of columns is identical to the length of skeletonTimes. Be aware that the skeleton points themselves are NOT samples from the target distribution.

Velocities: Matrix whose columns are velocities just after switches. The number of columns is identical to the length of skeletonTimes.

Examples

V <- matrix(c(3,1,1,3),nrow=2)
mu <- c(2,2)
result <- ZigZagGaussian(V, mu, 100)
plot(result$Positions[1,], result$Positions[2,],type='l',asp=1)

ZigZagIIDGaussian

Description

Applies the Zig-Zag Sampler to a IID Gaussian distribution

Usage

ZigZagIIDGaussian(variance, dim = 1L, n_iter = -1L, finalTime = -1,
  x0 = numeric(0), v0 = numeric(0))

Arguments

Value

Returns a list with the following objects:

Times: Vector of switching times

Positions: Matrix whose columns are locations of switches. The number of columns is identical to the length of skeletonTimes. Be aware that the skeleton points themselves are NOT samples from the target distribution.

Velocities: Matrix whose columns are velocities just after switches. The number of columns is identical to the length of skeletonTimes.

Examples

result <- ZigZagIIDGaussian(1, 2, 1000)
plot(result$Positions[2,], result$Positions[1,],type='l',asp=1)

ZigZagLogistic

Description

Applies the Zig-Zag Sampler to logistic regression, as detailed in Bierkens, Fearnhead, Roberts, The Zig-Zag Process and Super-Efficient Sampling for Bayesian Analysis of Big Data, 2019.

Usage

ZigZagLogistic(dataX, dataY, n_iter = -1L, finalTime = -1,
  x0 = numeric(0), v0 = numeric(0), cv = FALSE)

Arguments

Value

Returns a list with the following objects:

Times: Vector of switching times

Positions: Matrix whose columns are locations of switches. The number of columns is identical to the length of skeletonTimes. Be aware that the skeleton points themselves are NOT samples from the target distribution.

Velocities: Matrix whose columns are velocities just after switches. The number of columns is identical to the length of skeletonTimes.

Examples

require("RZigZag")

generate.logistic.data <- function(beta, n.obs) {
  dim <- length(beta)
  dataX <- cbind(rep(1.0,n.obs), matrix(rnorm((dim -1) * n.obs), ncol = dim -1));
  vals <- dataX %*% as.vector(beta)
    generateY <- function(p) { rbinom(1, 1, p)}
  dataY <- sapply(1/(1 + exp(-vals)), generateY)
    return(list(dataX = dataX, dataY = dataY))
}

beta <- c(1,2)
data <- generate.logistic.data(beta, 1000)
result <- ZigZagLogistic(data$dataX, data$dataY, 1000)
plot(result$Positions[1,], result$Positions[2,],type='l',asp=1)

ZigZagStudentT

Description

Applies the Zig-Zag Sampler to a Student T distribution (IID or spherically symmetric)

Usage

ZigZagStudentT(dof, dim = 1L, n_iter = -1L, finalTime = -1,
  x0 = numeric(0), v0 = numeric(0), sphericallySymmetric = TRUE)

Arguments

Value

Returns a list with the following objects:

Times: Vector of switching times

Positions: Matrix whose columns are locations of switches. The number of columns is identical to the length of skeletonTimes. Be aware that the skeleton points themselves are NOT samples from the target distribution.

Velocities: Matrix whose columns are velocities just after switches. The number of columns is identical to the length of skeletonTimes.

Examples

dim = 2
dof = 4
result <- ZigZagStudentT(dof, dim, n_iter=1000, sphericallySymmetric = TRUE)
plot(result$Positions[1,], result$Positions[2,],type='l',asp=1)