Package {qslice}

Title:	Implementations of Various Slice Samplers
Version:	0.4.0
Description:	Implementations of the quantile slice sampler of Heiner et al. (2026+, in preparation) <doi:10.48550/arXiv.2407.12608> as well as other popular slice samplers are provided. Helper functions for specifying pseudo-target distributions are included, both for diagnostics and for tuning the quantile slice sampler. Other implemented methods include the generalized elliptical slice sampler of Nishihara et al. (2014) https://jmlr.org/papers/v15/nishihara14a.html, latent slice sampler of Li and Walker (2023) <doi:10.1016/j.csda.2022.107652>, and stepping-out slice sampler of Neal (2003) <doi:10.1214/aos/1056562461>, as well as independence Metropolis-Hastings and random-walk Metropolis samplers.
License:	MIT + file LICENSE | Apache License 2.0
Depends:	R (≥ 4.1.0)
Encoding:	UTF-8
Suggests:	knitr, rmarkdown, testthat (≥ 3.0.0)
VignetteBuilder:	knitr
Config/roxygen2/version:	8.1.0
NeedsCompilation:	no
Packaged:	2026-10-21 20:50:55 UTC; dahl
Author:	Matthew Heiner [aut], David B. Dahl [aut, cre], Sam Johnson [aut]
Maintainer:	David B. Dahl <dahl@stat.byu.edu>
Repository:	CRAN
Date/Publication:	2026-10-22 00:00:09 UTC

Area Under the Curve (histogram)

Description

Function auc() is deprecated. Use utility_shrinkslice(..., utility_type = "AUC").

Usage

auc(u = NULL, x = NULL, y = NULL, nbins = 30)

Arguments

Value

Scalar value of the AUC (see utility_shrinkslice).

Independence Metropolis-Hastings

Description

Single update using an independence Metropolis-Hastings sampler with a pseudo-target as the proposal distribution.

Usage

imh_pseudo(x, log_target, pseudo, diagnostics = NULL)

Arguments

Value

A list containing the new state, x, whether the proposed value was accepted, logical accpt, and the number of evaluations of the target function, nEvaluations.

Examples

lf <- function(x) dbeta(x[1], 3, 4, log = TRUE) + dbeta(x[2], 5, 3, log = TRUE)
n_iter <- 10 # set to 1e3 for more complete illustration
draws <- matrix(0.2, nrow = n_iter, ncol = 2)
nAccpt <- 0L
pseudo <- list( list(ld = function(x) dbeta(x, 2, 2, log = TRUE),
                     q = function(u) qbeta(u, 2, 2)),
                list(ld = function(x) dbeta(x, 2, 2, log = TRUE),
                     q = function(u) qbeta(u, 2, 2))
)
for (i in seq.int(2, n_iter)) {
 out <- imh_pseudo(draws[i - 1, ], log_target = lf, pseudo = pseudo)
 draws[i,] <- out$x
 nAccpt <- nAccpt + out$accpt
 cat(i, '\r')
}
nAccpt / (nrow(draws) - 1)
plot(draws[,1], draws[,2], xlim = c(0, 1))
hist(draws[,1], freq = FALSE); curve(dbeta(x, 3, 4), col = "blue", add = TRUE)
hist(draws[,2], freq = FALSE); curve(dbeta(x, 5, 3), col = "blue", add = TRUE)

Pseudo-target from Laplace Approximation

Description

Find the location and scale for an approximating pseudo-target via Laplace approximation.

Usage

lapprox(
  log_target,
  init,
  family = "t",
  params = NULL,
  sc_adj = 1,
  lb = -Inf,
  ub = Inf,
  maxit = 100,
  ...
)

Arguments

Value

A list with the same outputs as pseudo_list; also includes opt, which gives output of optim.

Examples

pseu <- lapprox(function(x) dnorm(x, log = TRUE),
  family = "t",
  params = list(loc = NA, sc = NA, df = 5.0),
  init = 0.5, lb = -1.0)
curve(dnorm(x)/(1 - pnorm(-1)), from = -1, to = 6, col = "blue")
xx <- seq(-1, 6, length = 500)
lines(xx, sapply(xx, FUN = pseu$d), lty = 2, col = "red")

Sequence of conditional pseudo-targets from a realization

Description

Given a realization of a random vector, generate the corresponding sequence of conditional pseudo-targets (Heiner et al., 2026+). The pseudo-target is specified as a sequence of growing conditional distributions.

Usage

pseudo_condseq(x, pseudo_init, loc_fn, sc_fn, lb, ub)

Arguments

Details

See the documentation for slice_quantile_mv_seq for examples.

Value

A list containing a sequence of pseudo-targets, each from pseudo_list.

References

Heiner, M. J., Johnson, S. B., Christensen, J. R., and Dahl, D. B. (2026+), "Quantile Slice Sampling," arXiv preprint arXiv:2407.12608 doi:10.48550/arXiv.2407.12608

Examples

# Funnel distribution from Neal (2003).
K <- 10
n_iter <- 50 # MCMC iterations; set to 10e3 for more complete illustration
n <- 100 # number of iid samples from the target; set to 10e3 for more complete illustration
Y <- matrix(NA, nrow = n, ncol = K) # iid samples from the target
Y[, 1] <- rnorm(n, 0.0, 3.0)
for (i in 1:n) {
  Y[i, 2:K] <- rnorm(K - 1, 0.0, exp(0.5 * Y[i, 1]))
}
ltarget <- function(x) {
  dnorm(x[1], 0.0, 3.0, log = TRUE) +
    sum(dnorm(x[2:K], 0.0, exp(0.5 * x[1]), log = TRUE))
}
pseudo_control <- list(
  loc_fn = function(x) {
    0.0
  },
  sc_fn = function(x) {
    if (is.null(x)) {
      out <- 3.0
    } else {
      out <- exp(0.5 * x[1])
    }
    out
  },
  pseudo_init = pseudo_list(
    family = "t",
    params = list(loc = 0.0, sc = 3.0, df = 20),
    lb = -Inf,
    ub = Inf
  ),
  lb = rep(-Inf, K),
  ub = rep(Inf, K)
)
x0 <- runif(K)
draws <- matrix(rep(x0, n_iter + 1), nrow = n_iter + 1, byrow = TRUE)
draws_u <- matrix(rep(x0, n_iter), nrow = n_iter, byrow = TRUE)
n_eval <- 0
for (i in 2:(n_iter + 1)) {
  tmp <- slice_quantile_mv_seq(
    draws[i - 1, ],
    log_target = ltarget,
    pseudo_control = pseudo_control
  )
  draws[i, ] <- tmp$x
  draws_u[i - 1, ] <- tmp$u
  n_eval <- n_eval + tmp$nEvaluations
}
# (es <- coda::effectiveSize(coda::as.mcmc(draws)))
# mean(es)
n_eval / n_iter
sapply(1:K, function(k) {
  utility_shrinkslice(u = draws_u[, k], utility_type = "AUC")
})
hist(draws_u[, 1])
plot(draws[, 1], draws[, 2])
points(Y[, 1], Y[, 2], col = "blue", cex = 0.5)

Inverse transform from sequence of conditional pseudo-targets

Description

Given a vector from the unit hypercube, map to the original (back-transformed) vector using a sequence of conditional pseudo-target inverse CDFs. The pseudo-target is specified as a sequence of growing conditional distributions.

Usage

pseudo_condseq_XfromU(
  u,
  pseudo_init,
  loc_fn,
  sc_fn,
  lb,
  ub,
  diagnostics = NULL
)

Arguments

Details

See the documentation for slice_quantile_mv_seq for examples.

Value

A list containing x obtained from the sequence of inverse CDFs, and pseudo_seq, a list of the corresponding sequential pseudo-targets output from pseudo_list.

Examples

# Funnel distribution from Neal (2003).
K <- 10
n_iter <- 50 # MCMC iterations; set to 10e3 for more complete illustration
n <- 100 # number of iid samples from the target; set to 10e3 for more complete illustration
Y <- matrix(NA, nrow = n, ncol = K) # iid samples from the target
Y[, 1] <- rnorm(n, 0.0, 3.0)
for (i in 1:n) {
  Y[i, 2:K] <- rnorm(K - 1, 0.0, exp(0.5 * Y[i, 1]))
}
ltarget <- function(x) {
  dnorm(x[1], 0.0, 3.0, log = TRUE) +
    sum(dnorm(x[2:K], 0.0, exp(0.5 * x[1]), log = TRUE))
}
pseudo_control <- list(
  loc_fn = function(x) {
    0.0
  },
  sc_fn = function(x) {
    if (is.null(x)) {
      out <- 3.0
    } else {
      out <- exp(0.5 * x[1])
    }
    out
  },
  pseudo_init = pseudo_list(
    family = "t",
    params = list(loc = 0.0, sc = 3.0, df = 20),
    lb = -Inf,
    ub = Inf
  ),
  lb = rep(-Inf, K),
  ub = rep(Inf, K)
)
x0 <- runif(K)
draws <- matrix(rep(x0, n_iter + 1), nrow = n_iter + 1, byrow = TRUE)
draws_u <- matrix(rep(x0, n_iter), nrow = n_iter, byrow = TRUE)
n_eval <- 0
for (i in 2:(n_iter + 1)) {
  tmp <- slice_quantile_mv_seq(
    draws[i - 1, ],
    log_target = ltarget,
    pseudo_control = pseudo_control
  )
  draws[i, ] <- tmp$x
  draws_u[i - 1, ] <- tmp$u
  n_eval <- n_eval + tmp$nEvaluations
}
# (es <- coda::effectiveSize(coda::as.mcmc(draws)))
# mean(es)
n_eval / n_iter
sapply(1:K, function(k) {
  utility_shrinkslice(u = draws_u[, k], utility_type = "AUC")
})
hist(draws_u[, 1])
plot(draws[, 1], draws[, 2])
points(Y[, 1], Y[, 2], col = "blue", cex = 0.5)

Specify a pseudo-target within a given class

Description

Create a list of functions to evaluate a pseudo-target in a given class with supplied parameters (usually location and scale). The distribution is optionally truncated to specified bounds (and renormalized). See Heiner et al. (2026+).

Usage

pseudo_list(family, params, lb = -Inf, ub = Inf, name = NULL)

Arguments

Details

The supported classes of pseudo-targets include: t, cauchy, normal, logistic, and beta.

Value

A list with named components:

d: function to evaluate the density (finite boundary points of the support evaluate to 0)

ld: function to evaluate the log density (finite boundary points of the support evaluate to -Inf)

q: function to evaluate the quantile function

p: function to evaluate the distribution function

txt: text description of the distribution

params: repeats the params argument

lb: lower boundary of support

ub: upper boundary of support

References

Heiner, M. J., Johnson, S. B., Christensen, J. R., and Dahl, D. B. (2026+), "Quantile Slice Sampling," arXiv preprint arXiv:2407.12608 doi:10.48550/arXiv.2407.12608

Examples

pseu <- pseudo_list(family = "t", params = list(loc = 0.0, sc = 1.0, df = 5),
                    lb = 0.0, ub = Inf) # half t
str(pseu)
pseu$d(1.5)
pseu$ld(1.5)
pseu$p(1.5)
# should match (F(x) - F(0)) / (1 - F(0)) below
(pt(1.5, df = 5) - pt(0.0, df = 5)) / pt(0.0, df = 5, lower.tail = FALSE)
pseu$q(0.8060963)
curve(pseu$d(x), from = -1.0, to = 5.0, n = 1000)
pseu <- pseudo_list(family = "cauchy", params = list(loc = 0.0, sc = 1.0),
                    lb = 0.0, ub = Inf) # half Cauchy
str(pseu)
pseu$d(1.5)
pseu$ld(1.5)
pseu$p(1.5)
pseu$q(0.6256659)
pseu <- pseudo_list(family = "normal", params = list(loc = 0.0, sc = 1.0),
                    lb = 0.0, ub = Inf) # half normal
str(pseu)
pseu$d(1.5)
pseu$ld(1.5)
pseu$p(1.5)
pseu$q(0.8663856)
pseu <- pseudo_list(family = "logistic", params = list(loc = 0.0, sc = 1.0),
                    lb = 0.0, ub = Inf) # half logistic
str(pseu)
pseu$d(1.5)
pseu$ld(1.5)
pseu$p(1.5)
pseu$q(0.635149)
pseu <- pseudo_list(family = "beta", params = list(shape1 = 3.0, shape2 = 2.0), lb = 0.2)
str(pseu)
pseu$d(0.5)
pseu$ld(0.5)
pseu$p(0.5)
# should match (F(x) - F(0.2)) / (1 - F(0.2)) below
(pbeta(0.5, 3.0, 2.0) - pbeta(0.2, 3.0, 2.0)) / pbeta(0.2, 3.0, 2.0, lower.tail = FALSE)
pseu$q(0.2932771)
curve(pseu$d(x), from = 0.0, to = 1.0, n = 1000)

Optimal pseudo-target for a given target

Description

Find an optimal pseudo-target in a specified family to approximate the given (unnormalized) target (Heiner et al., 2026+). Optimize over the selected utility function.

Usage

pseudo_opt(
  log_target = NULL,
  samples = NULL,
  type = "samples",
  family = "t",
  df = NULL,
  lb = -Inf,
  ub = Inf,
  utility_type = "AUC",
  nbins = 100,
  tol_opt = 1e-06,
  init_fn = NULL,
  utility_options = list(tol_int = 0.001, n_opt = 10, interval_opt = c(1e-06, 0.999999),
    diagnostics = NULL),
  plot = TRUE,
  verbose = FALSE,
  degf = NULL
)

Arguments

Details

Optionally supply samples from the target distribution.

Value

A list with named components:

pseudo: a list with functions corresponding to the selected pseudo-target; output of pseudo_list.

utility: value of the utility function using the selected pseudo-target; output of utility_pseudo.

utility_type: repeats the input specifying the utility type.

opt: output of optim.

Other outputs repeating inputs.

References

Heiner, M. J., Johnson, S. B., Christensen, J. R., and Dahl, D. B. (2026+), "Quantile Slice Sampling," arXiv preprint arXiv:2407.12608 doi:10.48550/arXiv.2407.12608

Examples

oldpar <- par(mfrow = c(1,2))
(pseu <- pseudo_opt(samples = rnorm(1e3), type = "samples",
               family = "t", utility_type = "AUC",
               nbins = 10, plot = TRUE,
               verbose = FALSE))
(pseu <- pseudo_opt(log_target = function(x) dnorm(x, log = TRUE),
                type = "grid",
                family = "logistic", utility_type = "MSW",
                nbins = 30, plot = TRUE,
                verbose = FALSE))
(pseu <- pseudo_opt(log_target = function(x) dbeta(x, 4, 2, log = TRUE),
                lb = 0, ub = 1,
                type = "function",
                family = "cauchy", utility_type = "AUC",
                init_fn = function() {c(rbeta(1, 3.0, 3.0), rgamma(1, 3.0))},
                nbins = 30, plot = TRUE,
                verbose = FALSE))
par(oldpar)

Univariate Random Walk Metropolis with Normal Proposal

Description

Single update for the univariate random-walk Metropolis sampler with a normal proposal distribution.

Usage

rwm_norm(x, log_target, cand_sd, support = c(-Inf, Inf), diagnostics = NULL)

Arguments

Value

A list with three elements:

x is the new state.

accept Logical indicating whether the proposal was accepted.

nEvaluations is the number of evaluations of the target function used to obtain the new state.

Examples

lf <- function(x) dbeta(x, 3, 4, log = TRUE)
draws <- numeric(10) + 0.5 # set to numeric(1e3) for more complete illustration
nEvaluations <- 0L
nAccept <- 0L
for (i in seq.int(2, length(draws))) {
  out <- rwm_norm(draws[i - 1], log_target = lf, cand_sd = 0.3, support = c(0.0, 1.0))
  draws[i] <- out$x
  nEvaluations <- nEvaluations + out$nEvaluations
  nAccept <- nAccept + out$accept
}
nEvaluations / (length(draws) - 1)
nAccept / (length(draws) - 1)
plot(density(draws), xlim = c(0, 1))
curve(exp(lf(x)), 0, 1, col = "blue", add = TRUE)

Univariate Elliptical Slice Sampler

Description

Algorithm 1 of Nishihara et al. (2014) of the elliptical slice sampler of Murray et al. (2010).

Usage

slice_elliptical(x, log_lik, mu, sigma, diagnostics = NULL)

Arguments

Value

A list with two elements:

x is the new state.

nEvaluations is the number of evaluations of the log-likelihood function used to obtain the new state.

References

Murray, I., Adams, R., and MacKay, D., (2010), "Elliptical Slice Sampling," in Proceedings of the Thirteenth International Conference on Artificial Intelligence and Statistics, JMLR Workshop and Conference Proceedings. https://proceedings.mlr.press/v9/murray10a

Nishihara, R., Murray, I., and Adams, R. P. (2014), "Parallel MCMC with Generalized Elliptical Slice Sampling," Journal of Machine Learning Research, 15, 2087-2112. https://jmlr.org/papers/v15/nishihara14a.html

Examples

log_lik <- function(x) dbeta(x, 8, 2, log = TRUE)
lf <- function(x) log_lik(x) + dnorm(x, mean = mu, sd = sig, log = TRUE)
mu <- 0.2
sig <- 0.3
draws <- numeric(10) + 0.5 # set to numeric(1e3) + 0.5 for more complete illustration
nEvaluations <- 0L
for (i in seq.int(2, length(draws))) {
  out <- slice_elliptical(draws[i - 1], log_lik = log_lik, mu = mu, sigma = sig)
  draws[i] <- out$x
  nEvaluations <- nEvaluations + out$nEvaluations
}
nEvaluations / (length(draws) - 1)
plot(density(draws), xlim = c(-0.2, 1.2))
curve(exp(lf(x)) * 4.5, -0.2, 1.2, col = "blue", add = TRUE) # 4.5 to match curve height

Multivariate Elliptical Slice Sampler

Description

Algorithm 1 of Nishihara et al. (2014) of the elliptical slice sampler of Murray et al. (2010).

Usage

slice_elliptical_mv(x, log_lik, mu, Sig, is_chol = FALSE, diagnostics = NULL)

Arguments

Value

A list containing two elements: x is the new state and nEvaluations is the number of evaluations of the log-likelihood function used to obtain the new state.

References

Murray, I., Adams, R., and MacKay, D., (2010), "Elliptical Slice Sampling," in Proceedings of the Thirteenth International Conference on Artificial Intelligence and Statistics, JMLR Workshop and Conference Proceedings. https://proceedings.mlr.press/v9/murray10a

Nishihara, R., Murray, I., and Adams, R. P. (2014), "Parallel MCMC with Generalized Elliptical Slice Sampling," Journal of Machine Learning Research, 15, 2087-2112. https://jmlr.org/papers/v15/nishihara14a.html

Examples

mu <- c(0.5, 0.5)
Sig <- matrix(c(0.5, 0.25, 0.25, 0.5), nrow = 2)
SigL <- t(chol(Sig))
log_lik <- function(x) dbeta(x[1], 3, 4, log = TRUE) + dbeta(x[2], 5, 3, log = TRUE)
lf <- function(x) log_lik(x) - 0.5 * drop(crossprod(forwardsolve(SigL, x - mu)))
n_iter <- 10 # set to 1e3 for more complete illustration
draws <- matrix(0.3, nrow = n_iter, ncol = 2)
nEvaluations <- 0L
for (i in seq.int(2, n_iter)) {
  out <- slice_elliptical_mv(draws[i - 1,], log_lik = log_lik,
              mu = mu, Sig = SigL, is_chol = TRUE)
  draws[i,] <- out$x
  nEvaluations <- nEvaluations + out$nEvaluations
}
nEvaluations / (n_iter - 1)
plot(draws[,1], draws[,2], xlim = c(0, 1))
hist(draws[,1], freq = FALSE); curve(dbeta(x, 3, 4), col = "blue", add = TRUE)
hist(draws[,2], freq = FALSE); curve(dbeta(x, 5, 3), col = "blue", add = TRUE)

Generalized Elliptical Slice Sampler (univariate)

Description

Single update using the generalized elliptical slice sampler of Nishihara et al. (2014).

Usage

slice_genelliptical(x, log_target, mu, sigma, df, diagnostics = NULL)

Arguments

Value

A list containing two elements:

x is the new state.

nEvaluations is the number of evaluations of the target function used to obtain the new state.

References

Nishihara, R., Murray, I., and Adams, R. P. (2014), "Parallel MCMC with Generalized Elliptical Slice Sampling," Journal of Machine Learning Research, 15, 2087-2112. https://jmlr.org/papers/v15/nishihara14a.html

Examples

lf <- function(x) dbeta(x, 3, 4, log = TRUE)
draws <- numeric(10) + 0.5 # set to numeric(1e3) + 0.5 for more complete illustration
nEvaluations <- 0L
for (i in seq.int(2, length(draws))) {
  out <- slice_genelliptical(draws[i - 1], log_target = lf,
                                      mu = 0.5, sigma = 1, df = 5)
  draws[i] <- out$x
  nEvaluations <- nEvaluations + out$nEvaluations
}
nEvaluations / (length(draws) - 1)
plot(density(draws), xlim = c(0, 1))
curve(exp(lf(x)), 0, 1, col = "blue", add = TRUE)

Generalized Elliptical Slice Sampler (Multivariate)

Description

Generalized Elliptical Slice Sampler, Algorithm 2 of Nishihara et al. (2014)

Usage

slice_genelliptical_mv(
  x,
  log_target,
  mu,
  Sig,
  df,
  is_chol = FALSE,
  diagnostics = NULL
)

Arguments

Value

A list containing two elements: x is the new state and nEvaluations is the number of evaluations of the target function used to obtain the new state.

References

Nishihara, R., Murray, I., and Adams, R. P. (2014), "Parallel MCMC with Generalized Elliptical Slice Sampling," Journal of Machine Learning Research, 15, 2087-2112. https://jmlr.org/papers/v15/nishihara14a.html

Examples

lf <- function(x) dbeta(x[1], 3, 4, log = TRUE) + dbeta(x[2], 5, 3, log = TRUE)
n_iter <- 10 # set to 1e4 for more complete illustration
draws <- matrix(0.3, nrow = n_iter, ncol = 2)
nEvaluations <- 0L
for (i in seq.int(2, n_iter)) {
  out <- slice_genelliptical_mv(draws[i - 1,], log_target = lf,
              mu = c(0.5, 0.5), Sig = matrix(c(0.5, 0.25, 0.25, 0.5), nrow = 2),
              df = 5)
  draws[i,] <- out$x
  nEvaluations <- nEvaluations + out$nEvaluations
}
nEvaluations / (n_iter - 1)
plot(draws[,1], draws[,2], xlim = c(0, 1))
hist(draws[,1], freq = FALSE); curve(dbeta(x, 3, 4), col = "blue", add = TRUE)
hist(draws[,2], freq = FALSE); curve(dbeta(x, 5, 3), col = "blue", add = TRUE)

Multivariate Slice Sampler with Shrinking Hyperrectangle

Description

Multivariate slice sampler in Algorithm 8 of Neal (2003) using the "shrinkage" procedure.

Usage

slice_hyperrect(
  x,
  log_target,
  w = NULL,
  L = NULL,
  R = NULL,
  diagnostics = NULL
)

Arguments

Value

A list containing two elements: "x" is the new state and "nEvaluations" is the number of evaluations of the target function used to obtain the new state.

References

Neal, R. M. (2003), "Slice sampling," The Annals of Statistics, 31, 705-767. doi:10.1214/aos/1056562461

Examples

lf <- function(x) dbeta(x[1], 3, 4, log = TRUE) + dbeta(x[2], 5, 3, log = TRUE)
n_iter <- 10 # set to 1e4 for more complete illustration
draws <- matrix(0.2, nrow = n_iter, ncol = 2)
nEvaluations <- 0L
for (i in seq.int(2, n_iter)) {
 out <- slice_hyperrect(draws[i - 1, ], log_target = lf, w = c(0.5, 0.5))
 draws[i,] <- out$x
 nEvaluations <- nEvaluations + out$nEvaluations
 cat(i, '\r')
}
nEvaluations / (nrow(draws) - 1)
plot(draws[,1], draws[,2], xlim = c(0, 1))
hist(draws[,1], freq = FALSE); curve(dbeta(x, 3, 4), col = "blue", add = TRUE)
hist(draws[,2], freq = FALSE); curve(dbeta(x, 5, 3), col = "blue", add = TRUE)

Latent Slice Sampler

Description

Single update using the latent slice sampler of Li and Walker (2023).

Usage

slice_latent(x, s, log_target, rate, diagnostics = NULL)

Arguments

Value

A list containing three elements:

x is the new state of the target variable.

s is the new state of the latent scale variable.

nEvaluations is the number of evaluations of the target function used to obtain the new state.

References

Li, Y. and Walker, S. G. (2023), "A latent slice sampling algorithm," Computational Statistics and Data Analysis, 179, 107652. doi:10.1016/j.csda.2022.107652

Examples

lf <- function(x) dbeta(x, 3, 4, log = TRUE)
draws <- numeric(10) + 0.5 # set to numeric(1e3) + 0.5 for more complete illustration
nEvaluations <- 0L
s <- 0.5
for (i in seq.int(2, length(draws))) {
  out <- slice_latent(draws[i - 1], s, log_target = lf, rate = 0.3)
  draws[i] <- out$x
  s <- out$s
  nEvaluations <- nEvaluations + out$nEvaluations
}
nEvaluations / (length(draws) - 1)
plot(density(draws), xlim = c(0, 1))
curve(exp(lf(x)), 0, 1, col = "blue", add = TRUE)

Quantile Slice Sampler

Description

Single update using a quantile slice sampler of Heiner et al. (2026+).

Usage

slice_quantile(x, log_target, pseudo, diagnostics = NULL)

Arguments

Details

The sampler works with \psi, the pseudo-target CDF of the state, in double precision. Because floating-point numbers are dense near zero but not near one, states whose pseudo-target upper-tail probability is below about 1e-16 cannot be represented. Use a pseudo-target whose right tail is at least as heavy as that of the target (e.g., Student t).

Value

A list containing three elements:

x is the new state.

u is the value of the CDF of the pseudo-target associated with the returned value (also referred to as psi).

nEvaluations is the number of evaluations of the target function used to obtain the new state.

References

Heiner, M. J., Johnson, S. B., Christensen, J. R., and Dahl, D. B. (2026+), "Quantile Slice Sampling," arXiv preprint arXiv:2407.12608 doi:10.48550/arXiv.2407.12608

Examples

lf <- function(x) dbeta(x, 3, 4, log = TRUE)
pseu <- pseudo_list(family = "beta", params = list(shape1 = 2.0, shape2 = 2.0))
draws <- numeric(10) + 0.5 # set to numeric(1e3) + 0.5 for more complete illustration
nEvaluations <- 0L
for (i in seq.int(2, length(draws))) {
  out <- slice_quantile(draws[i - 1], log_target = lf, pseudo = pseu)
  draws[i] <- out$x
  nEvaluations <- nEvaluations + out$nEvaluations
}
nEvaluations / (length(draws) - 1)
plot(density(draws), xlim = c(0, 1))
curve(exp(lf(x)), 0, 1, col = "blue", add = TRUE)

Multivariate Quantile Slice Sampler

Description

Quantile slice sampler for a random vector (Heiner et al., 2026+). The pseudo-target is specified through independent univariate distributions.

Usage

slice_quantile_mv(x, log_target, pseudo, diagnostics = NULL)

Arguments

Value

A list containing three elements: "x" is the new state, "u" is the vector of values of the pseudo-target CDFs associated with the returned value, and "nEvaluations" is the number of evaluations of the target function used to obtain the new state.

References

Heiner, M. J., Johnson, S. B., Christensen, J. R., and Dahl, D. B. (2026+), "Quantile Slice Sampling," arXiv preprint arXiv:2407.12608 doi:10.48550/arXiv.2407.12608

Examples

lf <- function(x) dbeta(x[1], 3, 4, log = TRUE) + dbeta(x[2], 5, 3, log = TRUE)
ps_shsc <- list(c(2, 2), c(2, 1))
ps <- list(
  list(ld = function(x) dbeta(x, ps_shsc[[1]][1], ps_shsc[[1]][2], log = TRUE),
       p = function(x) pbeta(x, ps_shsc[[1]][1], ps_shsc[[1]][2]),
       q = function(x) qbeta(x, ps_shsc[[1]][1], ps_shsc[[1]][2]) ),
  list(ld = function(x) dbeta(x, ps_shsc[[2]][1], ps_shsc[[2]][2], log = TRUE),
       p = function(x) pbeta(x, ps_shsc[[2]][1], ps_shsc[[2]][2]),
       q = function(x) qbeta(x, ps_shsc[[2]][1], ps_shsc[[2]][2]) )
  )
n_iter <- 10 # set to 1e4 for more complete illustration
draws <- matrix(0.2, nrow = n_iter, ncol = 2)
draws_u <- draws
draws_u[1,] <- sapply(1:length(ps), function(k) ps[[k]]$p(draws[1,k]))
nEvaluations <- 0L
for (i in seq.int(2, n_iter)) {
  out <- slice_quantile_mv(draws[i - 1, ], log_target = lf, pseudo = ps)
  draws[i,] <- out$x
  draws_u[i,] <- out$u
  nEvaluations <- nEvaluations + out$nEvaluations
  cat(i, '\r')
}
nEvaluations / (nrow(draws) - 1)
plot(draws[,1], draws[,2], xlim = c(0, 1))
hist(draws[,1], freq = FALSE); curve(dbeta(x, 3, 4), col = "blue", add = TRUE)
hist(draws[,2], freq = FALSE); curve(dbeta(x, 5, 3), col = "blue", add = TRUE)
plot(draws_u[,1], draws_u[,2], xlim = c(0, 1))
utility_shrinkslice(u = draws_u[,1], nbins = 30, type = "samples",
plot = TRUE, utility_type = "AUC")
utility_shrinkslice(u = draws_u[,2], nbins = 30, type = "samples",
plot = TRUE, utility_type = "AUC")

Multivariate Quantile Slice Sampler from a sequence of conditional pseudo-targets

Description

Quantile slice sampler for a random vector (Heiner et al., 2026+). The pseudo-target is specified as a sequence of growing conditional distributions.

Usage

slice_quantile_mv_seq(x, log_target, pseudo_control, diagnostics = NULL)

Arguments

Value

A list containing three elements: "x" is the new state, "u" is a vector of values of the sequence of conditional CDFs of the pseudo-targets associated with the returned value, and "nEvaluations" is the number of evaluations of the target function used to obtain the new state.

References

Heiner, M. J., Johnson, S. B., Christensen, J. R., and Dahl, D. B. (2026+), "Quantile Slice Sampling," arXiv preprint arXiv:2407.12608 doi:10.48550/arXiv.2407.12608

Examples

# Funnel distribution from Neal (2003).
K <- 10
n_iter <- 50 # MCMC iterations; set to 10e3 for more complete illustration
n <- 100 # number of iid samples from the target; set to 10e3 for more complete illustration
Y <- matrix(NA, nrow = n, ncol = K) # iid samples from the target
Y[, 1] <- rnorm(n, 0.0, 3.0)
for (i in 1:n) {
  Y[i, 2:K] <- rnorm(K - 1, 0.0, exp(0.5 * Y[i, 1]))
}
ltarget <- function(x) {
  dnorm(x[1], 0.0, 3.0, log = TRUE) +
    sum(dnorm(x[2:K], 0.0, exp(0.5 * x[1]), log = TRUE))
}
pseudo_control <- list(
  loc_fn = function(x) {
    0.0
  },
  sc_fn = function(x) {
    if (is.null(x)) {
      out <- 3.0
    } else {
      out <- exp(0.5 * x[1])
    }
    out
  },
  pseudo_init = pseudo_list(
    family = "t",
    params = list(loc = 0.0, sc = 3.0, df = 20),
    lb = -Inf,
    ub = Inf
  ),
  lb = rep(-Inf, K),
  ub = rep(Inf, K)
)
x0 <- runif(K)
draws <- matrix(rep(x0, n_iter + 1), nrow = n_iter + 1, byrow = TRUE)
draws_u <- matrix(rep(x0, n_iter), nrow = n_iter, byrow = TRUE)
n_eval <- 0
for (i in 2:(n_iter + 1)) {
  tmp <- slice_quantile_mv_seq(
    draws[i - 1, ],
    log_target = ltarget,
    pseudo_control = pseudo_control
  )
  draws[i, ] <- tmp$x
  draws_u[i - 1, ] <- tmp$u
  n_eval <- n_eval + tmp$nEvaluations
}
# (es <- coda::effectiveSize(coda::as.mcmc(draws)))
# mean(es)
n_eval / n_iter
sapply(1:K, function(k) {
  utility_shrinkslice(u = draws_u[, k], utility_type = "AUC")
})
hist(draws_u[, 1])
plot(draws[, 1], draws[, 2])
points(Y[, 1], Y[, 2], col = "blue", cex = 0.5)

Slice sampler using the Stepping Out and Shrinkage Procedures

Description

Single update for the univariate slice sampler of Neal (2003) using the "stepping out" procedure, followed by the "shrinkage" procedure.

Usage

slice_stepping_out(x, log_target, w, max = Inf, diagnostics = NULL)

Arguments

Value

A list with two elements:

x is the new state.

nEvaluations is the number of evaluations of the target function used to obtain the new state.

References

Neal, R. M. (2003), "Slice sampling," The Annals of Statistics, 31, 705-767. doi:10.1214/aos/1056562461

Examples

lf <- function(x) dbeta(x, 3, 4, log = TRUE)
draws <- numeric(10) + 0.5 # set to numeric(1e3) for more complete illustration
nEvaluations <- 0L
for (i in seq.int(2, length(draws))) {
  out <- slice_stepping_out(draws[i - 1], log_target = lf, w = 0.7, max = Inf)
  draws[i] <- out$x
  nEvaluations <- nEvaluations + out$nEvaluations
}
nEvaluations / (length(draws) - 1)
plot(density(draws), xlim = c(0, 1))
curve(exp(lf(x)), 0, 1, col = "blue", add = TRUE)

Utility for a given univariate target and pseudo-target

Description

Takes a pseudo-target and target (or samples from the target) and evaluates the utility function for the transformed target, which can be one of Area Under the Curve (AUC) and Mean Slice Width (MSW). See Heiner et al. (2026+).

Usage

utility_pseudo(
  pseudo,
  log_target = NULL,
  samples = NULL,
  type = "samples",
  x = NULL,
  nbins = 30,
  plot = TRUE,
  utility_type = "AUC",
  options = list(tol_int = 0.001, n_opt = 10, interval_opt = c(1e-06, 0.999999),
    diagnostics = NULL)
)

Arguments

Details

Optionally plot the target and pseudo-target densities as well as the transformed target.

Value

Scalar value of the utility function evaluation.

References

Heiner, M. J., Johnson, S. B., Christensen, J. R., and Dahl, D. B. (2026+), "Quantile Slice Sampling," arXiv preprint arXiv:2407.12608 doi:10.48550/arXiv.2407.12608

Examples

pseu <- pseudo_list(family = "t", params = list(loc = 0.25, sc = 1.0, df = 1))
ltargdens <- function(x) dnorm(x, log = TRUE)
oldpar <- par(mfrow = c(1,2))
utility_pseudo(pseudo = pseu, log_target = ltargdens, type = "function", utility_type = "MSW")
utility_pseudo(pseudo = pseu, log_target = ltargdens, type = "function", utility_type = "AUC")
samp <- rnorm(10e3)
utility_pseudo(pseudo = pseu, samples = samp, nbins = 30, type = "samples", utility_type = "MSW")
utility_pseudo(pseudo = pseu, samples = samp, nbins = 30, type = "samples", utility_type = "AUC")
par(oldpar)
pseu <- pseudo_list(family = "beta", params = list(shape1 = 1.5, shape2 = 1.5))
ltargdens <- function(x) dbeta(x, shape1 = 2.0, shape2 = 3.0, log = TRUE)
oldpar <- par(mfrow = c(1,2))
utility_pseudo(pseudo = pseu, log_target = ltargdens, type = "function", utility_type = "MSW")
utility_pseudo(pseudo = pseu, log_target = ltargdens, type = "function", utility_type = "AUC")
samp <- rbeta(10e3, shape1 = 2.0, shape2 = 3.0)
utility_pseudo(pseudo = pseu, samples = samp, nbins = 30, type = "samples", utility_type = "MSW")
utility_pseudo(pseudo = pseu, samples = samp, nbins = 30, type = "samples", utility_type = "AUC")
par(oldpar)

Utility for a Transformed Target

Description

Evaluates the utility function for a transformed target, which can be one of Area Under the Curve (AUC) and Mean Slice Width (MSW).

Usage

utility_shrinkslice(
  h = NULL,
  x = NULL,
  u = NULL,
  type = "samples",
  nbins = 30,
  plot = FALSE,
  utility_type = "AUC",
  options = list(tol_int = 0.001, n_opt = 10, interval_opt = c(1e-06, 0.999999))
)

Arguments

Details

Calculating the utility with a function and type = 'integration' employs numerical optimization (for AUC) and numerical integration (for both AUC and MSW), which may not succeed for highly multimodal or narrow target densities.

Value

Scalar value of the utility function evaluation.

References

Heiner, M. J., Johnson, S. B., Christensen, J. R., and Dahl, D. B. (2026+), "Quantile Slice Sampling," arXiv preprint arXiv:2407.12608 doi:10.48550/arXiv.2407.12608

Examples

a <- 4
b <- 2
u <- rbeta(10e3, a, b)
h <- function(x) dbeta(x, a, b)
x <- sort(runif(20))
type <- "AUC"
type <- "MSW"

utility_shrinkslice(
  h = h,
  x = NULL,
  u = NULL,
  type = "integration",
  plot = TRUE,
  utility_type = type
)

utility_shrinkslice(
  h = h,
  x = NULL,
  u = NULL,
  nbins = 50,
  type = "grid",
  plot = TRUE,
  utility_type = type
)

utility_shrinkslice(
  h = h,
  x = x,
  u = NULL,
  type = "grid",
  plot = TRUE,
  utility_type = type
)

utility_shrinkslice(
  h = NULL,
  x = NULL,
  u = u,
  nbins = 30,
  type = "samples",
  plot = TRUE,
  utility_type = type
)

utility_shrinkslice(
  h = NULL,
  x = x,
  u = u,
  type = "samples",
  plot = TRUE,
  utility_type = type
)