---
title: Blackbox MCMC Testing
format:
    html:
        copy-code: true
        embed-resources: true
fig-responsive: true
fig-width: 8
project:
    execute-dir: project
---

```{r}
library(tidyverse)
source('blackbox.r')
```

Making sure initially that we can plot the mystery function:
```{r}
xs = seq(-30,30, 0.01)
data = tibble(x = xs, y = mystery(xs))
ggplot(data=data, aes(x=x, y=y)) + geom_point()
```

Setting up our sampler:

```{r}
mcmc_sampler <- function(func, start, step_size, num_steps) {
    samples = numeric(num_steps)
    samples[1] = start
    accepted = 0
    for (i in 2:num_steps) {
        current = samples[i-1]
        proposed = rnorm(1, current, step_size)
        ratio = func(proposed) / func(current)
        if (ratio > runif(1)) {
            samples[i] = proposed
            accepted = accepted + 1
        } else {
            samples[i] = current
        }
    }
    return(list(samples=samples, accept_rate=accepted / length(samples)))
}
```

Running the sampler
```{r}
out = mcmc_sampler(mystery, -10, 4, 50000)
samples = out$samples
accept_rate = out$accept_rate
accept_rate
```

Looking at the trace plot:
```{r}
data = tibble(iter=1:length(samples), x=samples)
ggplot(data=data, aes(x=iter, y=x)) + geom_line()
```

Looking at the lag plot:
```{r}
burned = data %>% filter(iter > 1000)
ggplot(data=burned, aes(x=x, y=lag(x))) + geom_point(alpha=0.1)
```

And looking at the overall histogram:
```{r}
ggplot(data=burned, aes(x=x)) + geom_histogram(bins=100)
```