harmonicHMC: Sample from a truncated Gaussian distribution with the...

In hdtg: Generate Samples from Multivariate Truncated Normal Distributions

Sample from a truncated Gaussian distribution with the harmonic HMC

Description

Generate MCMC samples from a d-dimensional truncated Gaussian distribution with constraints Fx+g >= 0 using the Harmonic Hamiltonian Monte Carlo sampler (Harmonic-HMC).

Usage

harmonicHMC(
  nSample,
  burnin = 0,
  mean,
  choleskyFactor,
  constrainDirec,
  constrainBound,
  init,
  time = c(pi/8, pi/2),
  precFlg,
  seed = NULL,
  extraOutputs = c()
)

Arguments

nSample	number of samples after burn-in.
burnin	number of burn-in samples (default = 0).
mean	a d-dimensional mean vector.
choleskyFactor	upper triangular matrix R from Cholesky decomposition of precision or covariance matrix into R^TR.
constrainDirec	the k-by-d F matrix (k is the number of linear constraints).
constrainBound	the k-dimensional g vector.
init	a d-dimensional vector of the initial value. init must satisfy all constraints.
time	HMC integration time for each iteration. Can either be a scalar value for a fixed time across all samples, or a length 2 vector of a lower and upper bound for uniform distribution from which the time is drawn from for each iteration.
precFlg	logical. whether choleskyFactor is from precision (TRUE) or covariance matrix (FALSE).
seed	random seed (default = 1).
extraOutputs	vector of strings. "numBounces" and/or "bounceDistances" can be requested, with the latter containing the distances in-between bounces for each sample and hence incurring significant computational and memory costs.

Value

samples: nSample-by-d matrix of samples or, if extraOutputs is non-empty, a list of samples and the extra outputs.

References

\insertRef

pakman2014exacthdtg

Examples

set.seed(1)
d <- 10
A <- matrix(runif(d^2)*2 - 1, ncol=d)
Sigma <- t(A) %*% A
R <- cholesky(Sigma)
mu <- rep(0, d)
constrainDirec <- diag(d)
constrainBound <- rep(0,d)
initial <- rep(1, d)
results <- harmonicHMC(1000, 1000, mu, R, constrainDirec, constrainBound, initial, precFlg = FALSE)

hdtg documentation built on June 8, 2025, 1:07 p.m.