## ----include = FALSE---------------------------------------------------------- knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) ## ----setup-------------------------------------------------------------------- library(rmcmc) ## ----------------------------------------------------------------------------- delta_sq <- 0.1 target_distribution <- list( log_density = function(x) -sum(sqrt(delta_sq + x^2)), gradient_log_density = function(x) -x / sqrt(delta_sq + x^2) ) ## ----------------------------------------------------------------------------- sample_bimodal <- function(dimension) { sigma <- 0.1 return( sample(c(-1, 1), dimension, TRUE) * sqrt(1 - sigma^2) + stats::rnorm(dimension) * sigma ) } ## ----------------------------------------------------------------------------- normal_noise_proposal <- barker_proposal() bimodal_noise_proposal <- barker_proposal(sample_auxiliary = sample_bimodal) ## ----------------------------------------------------------------------------- set.seed(7861932L) dimension <- 100 initial_state <- rnorm(dimension) adapters <- list(scale_adapter()) n_warm_up_iteration <- 10000 n_main_iteration <- 10000 ## ----------------------------------------------------------------------------- bimodal_noise_results <- sample_chain( target_distribution = target_distribution, proposal = bimodal_noise_proposal, initial_state = initial_state, n_warm_up_iteration = n_warm_up_iteration, n_main_iteration = n_main_iteration, adapters = adapters ) ## ----------------------------------------------------------------------------- normal_noise_results <- sample_chain( target_distribution = target_distribution, proposal = normal_noise_proposal, initial_state = initial_state, n_warm_up_iteration = n_warm_up_iteration, n_main_iteration = n_main_iteration, adapters = adapters ) ## ----------------------------------------------------------------------------- expected_square_jumping_distance <- function(traces) { n_iteration <- nrow(traces) mean(rowSums(traces[2:n_iteration, ] - traces[1:n_iteration - 1, ])^2) } ## ----------------------------------------------------------------------------- cat( sprintf( "Expected square jumping distance using normal noise is %.2f", expected_square_jumping_distance(bimodal_noise_results$traces[, 1:dimension]) ) ) ## ----------------------------------------------------------------------------- cat( sprintf( "Expected square jumping distance using bimodal noise is %.2f", expected_square_jumping_distance(normal_noise_results$traces[, 1:dimension]) ) ) ## ----------------------------------------------------------------------------- library(posterior) ## ----------------------------------------------------------------------------- cat( sprintf( "Estimated ESS of mean(target_log_density) using bimodal noise is %.0f", ess_mean( extract_variable( as_draws_matrix(bimodal_noise_results$traces), "target_log_density" ) ) ) ) ## ----------------------------------------------------------------------------- cat( sprintf( "Estimated ESS of mean(target_log_density) using normal noise is %.0f", ess_mean( extract_variable( as_draws_matrix(normal_noise_results$traces), "target_log_density" ) ) ) )