rmcelreath / rethinking

Statistical Rethinking course and book package
2.14k stars 601 forks source link

Random walk in ulam #423

Open lukaseamus opened 7 months ago

lukaseamus commented 7 months ago

I have been trying to implement a random walk in ulam. I'd be interested in how to do this in general after reading discussions such as this and this on how to build random walks in Stan and its other R interfaces such as brms.

Here is a MRE for a specific case trying to model the decay parameter k of the exponential decay function $f(x)=ae^{-kx}$ as a function of x using random walk:

# generate data
set.seed(10)
x <- 1:250
r <- rnorm(x, mean = 0, sd = sqrt(x^-1.5))
y <- 2 * exp(-0.04 * x) + r

require(tidyverse)
df <- tibble(x = x, y = y)
ggplot(df, aes(x, y)) +
  geom_point() +
  theme_minimal()

# prepare data for ulam
l <- as.list(df)
l$N <- length(l$x)

# run model
require(rethinking)
rw.mod <- ulam(
  alist(
    y ~ dnorm(mean = mu, sd = sigma),
    mu <- a * exp(-k * x),
    k <- 0.04,
    for (i in 2:N) {
      k[i] <- k[i - 1] + r
    },
    r ~ dnorm(mean = 0, sd = tau),
    a ~ dlnorm(meanlog = log(2), sdlog = 0.5),
    c(sigma, tau) ~ dexp(rate = 1)
  ), data = l, chains = 8, cores = parallel::detectCores(), iter = 1e4
)

This yields the following error: Error in get_dist_template(the_dist[1]) : No template for distribution ':'

I also tried specifying the vector 2:N in the for loop as seq(from = 2, to = N, by = 1), in which case I get the same error: Error in get_dist_template(the_dist[1]) : No template for distribution 'seq'

Any solutions or pointers would be greatly appreciated!

rmcelreath commented 7 months ago

Hm, ulam doesn't have any convenient interface for loops or detailed manipulation of vectors. I think I could probably get this to work, but it would be a kludge.

I will save the example to the "requested features" file though.

lukaseamus commented 7 months ago

Thank you for the reply @rmcelreath, and for adding this as a feature request. I assume the same is true for my previous question about reparameterisation of likelihood functions?