README

Rtwalk is an R implementation of the t-walk, a general-purpose Markov Chain Monte Carlo (MCMC) sampler for continuous distributions, ideal for Bayesian inference problems.

Overview

This package provides an implementation of the t-walk algorithm, as originally proposed by Christen & Fox (2010). The t-walk is a robust, self-adjusting MCMC sampler, which means it does not require the tedious manual tuning of proposal parameters. It is designed to efficiently explore a wide range of target distributions, maintaining good performance even in high-dimensional or multimodal problems.

The implementation supports reproducible sequential and parallel chains, convergence summaries, visualization helpers, constrained parameter spaces, and numerically stable Metropolis–Hastings acceptance probabilities.

Installation

install.packages("Rtwalk")

Example Usage

Here is a simple example of how to use Rtwalk to sample from a bimodal distribution.

library(Rtwalk)

# Define the log-posterior density of the target distribution
# In this case, a mixture of two normal distributions
log_posterior_bimodal <- function(x) {
  log_component_1 <- log(0.5) + dnorm(x, mean = -3, sd = 0.5, log = TRUE)
  log_component_2 <- log(0.5) + dnorm(x, mean = 3, sd = 0.5, log = TRUE)
  max_log <- max(log_component_1, log_component_2)

  max_log + log(
    exp(log_component_1 - max_log) + exp(log_component_2 - max_log)
  )
}

initial_point_1 <- -3
initial_point_2 <- 3

result <- twalk(
  log_posterior = log_posterior_bimodal,
  n_iter = 50000,
  x0 = initial_point_1,
  xp0 = initial_point_2
)

burnin <- floor(nrow(result$samples) * 0.2)
samples <- result$samples[
  seq.int(burnin + 1L, nrow(result$samples)),
  ,
  drop = FALSE
]

par(mfrow = c(1, 2))
hist(samples, breaks = 50, freq = FALSE, 
     main = "Posterior Distribution", xlab = "Parameter Value")
lines(density(samples), col = "blue", lwd = 2)

plot(samples, type = 'l', col = "grey30", 
     main = "Trace Plot", xlab = "Iteration")

The primary time-ordered trajectory is stored in result$samples. The auxiliary t-walk trajectory is available as result$companion_samples. result$all_samples is retained for backward compatibility, but should not be interpreted as one time-ordered chain.

unlist(result[c("acceptance_rate", "move_rate", "no_move_rate")])

Here, acceptance_rate is the Metropolis–Hastings acceptance rate, move_rate counts accepted proposals that changed the state, and no_move_rate counts accepted identity proposals from the original t-walk transition kernel.

Parallel chains

Multiple reproducible chains can be run in parallel. Calling set.seed() before twalk() controls the parallel random-number streams:

set.seed(123)

parallel_result <- twalk(
  log_posterior = log_posterior_bimodal,
  n_iter = 25000,
  x0 = initial_point_1,
  xp0 = initial_point_2,
  n_chains = 2,
  n_cores = 2,
  show_progress = FALSE
)

summary(parallel_result, burnin_frac = 0.2)

Citation

This package is an implementation of the algorithm described in the following paper. If you use Rtwalk in your research, please cite the original work:

License