View source: R/hyperrectangle_slice.R
slice_quantile_mv_seq | R Documentation |
Quantile slice sampler for a random vector (Heiner et al., 2024+). The pseudo-target is specified as a sequence of growing conditional distributions.
slice_quantile_mv_seq(x, log_target, pseudo_control)
x |
The current state (as a numeric vector). |
log_target |
A function taking numeric vector that evaluates the log-target density, returning a numeric scalar. |
pseudo_control |
A list with
|
A list containing three elements: "x" is the new state, "u" is a vector of values of the sequence of conditional CDFs of the psuedo-targets associated with the returned value, and "nEvaluations is the number of evaluations of the target function used to obtain the new state.
Heiner, M. J., Johnson, S. B., Christensen, J. R., and Dahl, D. B. (2024+), "Quantile Slice Sampling," arXiv preprint arXiv:###
# 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, degf = 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) auc(u = draws_u[,k]))
hist(draws_u[,1])
plot(draws[,1], draws[,2])
points(Y[,1], Y[,2], col = "blue", cex = 0.5)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.