Skip to content

Instantly share code, notes, and snippets.

@crsh
Last active November 6, 2015 11:18
Show Gist options
  • Save crsh/debcbc02eca4e0a2bd46 to your computer and use it in GitHub Desktop.
Save crsh/debcbc02eca4e0a2bd46 to your computer and use it in GitHub Desktop.
Simulation of uninformativ datasets for Bayesian t-tests
n_bf <- 50
# Simulate data
library("BayesFactor")
previous_data <- rnorm(2, mean = 0, sd = 3)
uninformative_bf <- function(x, previous_data = NULL, mu = 0) {
bf <- as.vector(ttestBF(c(previous_data, x), mu = mu))
(1 - bf)^2
}
bfs <- as.vector(ttestBF(previous_data, mu = 0))
for(i in 1:(n_bf-1)) {
next_observation <- optim(rnorm(2, mean = 0, sd = 5), uninformative_bf, previous_data = previous_data, method = "BFGS")
previous_data <- c(previous_data, next_observation$par)
bfs <- c(bfs, as.vector(ttestBF(previous_data, mu = 0)))
}
# Plot data
plot(
NA
, NA
, xlim = c(1, n_bf)
, ylim = c(1/10, 10)
, axes = FALSE
, xlab = "n"
, ylab = bquote("BF"[10])
, log = "y"
, main = expression(paste("Bayesian t-test (", mu == 0, ", ", r == sqrt(2)/2, ")"))
)
for(i in c(1/10, 1/3, 1, 3, 10)) lines(x = c(0, n_bf), y = c(i, i), col = "grey", lty = "dashed")
lines(x = c(0, n_bf), y = c(1, 1), lwd = 2)
points(1:n_bf, bfs, pch = 21, col = "steelblue", bg = "lightblue")
axis(1)
axis(2, at = c(1/10, 1/3, 1, 3, 10), labels = c("1/10", "1/3", "1", "3", "10"), las = 1)
@crsh
Copy link
Author

crsh commented Nov 6, 2015

A particularly good example:
bayesian_ttest

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment