Hey! Maybe I’m completely missing the point, but I don’t see what’s wrong.

Let’s simulate some data:

```
n_dim <- 50
matCovSquaredFreqs <- diag(rnorm(n_dim, 1, 1)^2)
stan_data <- list(n_dim = n_dim, matCovSquaredFreqs = matCovSquaredFreqs)
```

Running your model gives reasonable results. (By the way, you should never use `multi_normal()`

with a diagonal matrix – that’s very inefficient. Actually, you could get the same results with any RNG, e.g. `rnorm()`

in R – I guess I’m really missing the point here?)

### Running with adapt_delta = 0.8 (default)

```
p0.8 <- sampling(sm, data = stan_data) # sm is the compiled stan model
tibble(Stan = summary(p0.8)$summary[-51,"sd"],
truth = sqrt(diag(matCovSquaredFreqs))) %>%
ggplot(aes(x=truth, y=Stan)) +
geom_abline() +
geom_point() +
ggtitle("True SD's vs. Stan results (delta_adapt = 0.8)")
bayesplot::nuts_params(p0.8) %>%
filter(str_detect(Parameter, "accept|step")) %>%
ggplot(aes(x = Iteration, y = Value, color = factor(Chain))) +
geom_line(show.legend = F) +
facet_grid(Parameter~Chain, scales = "free") +
ggtitle("NUTS parameters with delta_adapt = 0.8 (default)")
```

### Running with adapt_delta = 0.9999

```
p0.9999 <- sampling(sm, data = stan_data, control = list(adapt_delta = 0.9999))
tibble(Stan = summary(p0.9999)$summary[-51,"sd"],
truth = sqrt(diag(matCovSquaredFreqs))) %>%
ggplot(aes(x=truth, y=Stan)) +
geom_abline() +
geom_point() +
ggtitle("True SD's vs. Stan results (adapt_delta = 0.9999)")
bayesplot::nuts_params(p0.9999) %>%
filter(str_detect(Parameter, "accept|step")) %>%
ggplot(aes(x = Iteration, y = Value, color = factor(Chain))) +
geom_line(show.legend = F) +
facet_grid(Parameter~Chain, scales = "free") +
ggtitle("NUTS parameters with adapt_delta = 0.9999")
```

Hope this helps? If this went into a completely wrong direction, please feel free to ignore me. ;)

Cheers,

Max