Error when using adstock function in non-linear models


#1

I am currently trying to reproduce results I found on this post in brms. The goal is to see if I can recover adstock rates, coefficients, and intercept using a non-linear model from simulated data.

The only non-linear models I have tried are the ones in vignette(brms_nonlinear) so I am betting I am making an embarrassingly simple mistake.

The error I am getting is

Error in stanc(model_code = paste(program, collapse = “\n”), model_name = model_cppname, :
failed to parse Stan model ‘file290d3c2c341a’ due to the above error.

I am guessing this is due to the adstock function being part of my model.

My code is as follows:

# load libraries
library(purrr)
library(tidyverse)

# positive_round function 
positive_round <- function(...) round(pmax(..., 0), 0)

# adstock function
adstock<-function(x,rate=0){
  return(as.numeric(stats::filter(x=x,filter=rate,method="recursive")))
}

# base and ad constants 
base     <-  100 
n_weeks  <-  104
avg      <-  20
sdev     <-  10 

# adstock rates
ad1_rate <- .7
ad2_rate <- .4
ad3_rate <- .5

set.seed(42)
# 3 ads with mean
fake_df <- rep(avg, 3) %>%
  # normal distribution with 2 years of observations
  map(~ rnorm(n = n_weeks, mean = .x, sd = sdev)) %>% 
  # make sure there are no 0s and make values integers
  map(~ positive_round(0,.x)) %>% 
  # convert lists to data frame
  map_dfc(~as_tibble(.x)) %>% 
  # sensible names
  setNames(c("ad1", "ad2", "ad3")) %>% 
  # sales = sum of adstock and some noise
  mutate(sales = base + 
           adstock(ad1, ad1_rate) + 
           adstock(ad2, ad2_rate) + 
           adstock(ad3, ad3_rate) +
           rnorm(n(), sd = 5)) %>% 
  # round sales to whole number
  mutate(sales =  round(sales, 0))

# Model building
library(brms)
prior1 <- prior(normal(.5, .5), nlpar = "rate1") + 
          prior(normal(.5, .5), nlpar = "rate2") +
          prior(normal(.5, .5), nlpar = "rate3") +
          prior(normal(1, .5), nlpar = "b1") +
          prior(normal(1, .5), nlpar = "b2") +
          prior(normal(1, .5), nlpar = "b3") +
          prior(normal(100, 50), nlpar = "basesales") 
  
  
fit1 <- brm(bf(sales ~ basesales + 
                 b1 * adstock(ad1, rate1) + 
                 b2 * adstock(ad2, rate2) + 
                 b3 * adstock(ad3, rate3), 
               basesales + b1 + b2 + b3 + rate1 + rate2 + rate3 ~ 1, 
               nl = TRUE),
            data = fake_df, 
            prior = prior1)

summary(fit1)

#2

Can you include what the error message for the “above error” was? That will help us figure out what’s going on.

Thanks!


#3

Sorry

SYNTAX ERROR, MESSAGE(S) FROM PARSER:

No matches for: 

  adstock(real, real)

Function adstock not found.
  error in 'model290d7a10b5d2_file290d3c2c341a' at line 52, column 71
  -------------------------------------------------
    50:   for (n in 1:N) { 
    51:     // compute non-linear predictor 
    52:     mu[n] = mu_basesales[n] + mu_b1[n] * adstock(C_1[n] , mu_rate1[n]) + mu_b2[n] * adstock(C_2[n] , mu_rate2[n]) + mu_b3[n] * adstock(C_3[n] , mu_rate3[n]);
                                                                              ^
    53:   } 
  -------------------------------------------------

Error in stanc(model_code = paste(program, collapse = "\n"), model_name = model_cppname,  : 
  failed to parse Stan model 'file290d3c2c341a' due to the above error.

Looks like the issue is the adstock function. What is a possible work-around for this?


#4

It looks like the Stan code brms writes includes the adstock() function as specified in your formula, that is, bf() is interpreting the formula as the formula to use inside of Stan rather than a formula where the terms are evaluated in R (which I think it does when you don’t use bf() or if nl=FALSE). So yeah, that’s not going to work. adstock() would need to be a user-defined Stan function for this to work and not an R function. So either adstock() needs to be defined in the Stan language or, if possible, used only in R to precompute what you need without having to call adstock() inside of Stan.


#5

Thanks for the help. Unfortunately the rate in the adstock function is something that I need to optimize. So I don’t think that precomputing in R will work. The only option would be to define the function in Stan.


#6

Indeed, bf() interpretes the right-hand side of a formula literally if argument nl is set to TRUE. That’s the main idea of non-linear formulas in brms, basically. So yes, you need to define adstock() in Stan for this model to work.


#8

I have not worked with Stan directly. Since I have to define adstock() in Stan does that mean that I need to write the entire program in Stan?


#9

just the function as a string which you then pass to argument stan_funs of brm().


#10

Given that adstock is a wrapper around another R function (filter in this case) I’d recommend starting by rewriting it using more primitive R functions than filter. In my experience, doing that will often make it easier to translate R into Stan.


#12

After some time in the Stan User’s Manual and a few iterations I was able to produce the following

my_stan_adstock <- "
  vector my_stan_adstock(vector x, real rate, int N){
  vector[N] y;
  vector[N] my_vector;
      for(i in 1:N){ 
       y[0] = x[i];
       y[i] = x[i] + rate * y[i - 1];
       my_vector[i] = y[i];
     }
  return my_vector;
  }
"

This chuck of code seems like it should work, but when I use stan_funs in brms I get the following error:

SYNTAX ERROR, MESSAGE(S) FROM PARSER:

No matches for: 

  my_stan_adstock(real, real, int)

Available argument signatures for my_stan_adstock:

  my_stan_adstock(vector, real, int)

  error in 'model290d28a7de7f_file290d3c8a7f3e' at line 63, column 87
  -------------------------------------------------
    61:   for (n in 1:N) { 
    62:     // compute non-linear predictor 
    63:     mu[n] = mu_basesales[n] + mu_b1[n] * my_stan_adstock(C_1[n] , mu_ad1rate[n] , 104) + mu_b2[n] * my_stan_adstock(C_2[n] , mu_ad2rate[n] , 104) + mu_b3[n] * my_stan_adstock(C_3[n] , mu_ad3rate[n] , 104);
                                                                                              ^
    64:   } 
  -------------------------------------------------

If I understand this correctly, after translation the function is expected to be my_stan_adstock(real, real, int), but I am providing my_stan_adstock(vector, real, int). I would change that vector to a real, but I need a vector to calculate the output.

If anyone is interested in reproducing this error. My code is in a gist here.

My question at this point is how do I take in the whole vector as opposed to the real?


#13

It looks like the brms parser translates your formula into calls of the function within a loop. So this is going to depend on how flexible this feature in brms is. There would need to be a way to inform brms that you need to index some of the arguments differently. @paul.buerkner will know if it’s currently possible.

Either way, it looks like brms is getting you close to the right Stan code so you could also use make_stancode()it to generate the code and make_standata() to get the data. Then tweak the code so that it works with your function and use it with rstan to fit the model.


#14

We have a vignette coming imminently that includes a section about using brms for that purpose.


#15

brms expects functions to work in a pointwise manner (that is for each observation n). That is, your function should take the input necessary for the nth observation and return a real value for this observations. Is it possible to write your function this way?


#16

@alexhallam Looking back at your function are you sure you don’t need i-1 in some of the indexes? If so then the point wise function would need the arguments for the current i and also arguments for the previous output to be fed back in. I’m not sure if it’s possible with brms or not. As a fallback you can always modify the brms code.


#17

Thanks @jonah and @paul.buerkner for your help. I do need the y[i - 1] in there. I updated the code to reflect that. If brms expects the code in a pointwise manner then I may have to do some rethinking to solve this problem. I think that you provided great suggestions:

  1. Using make_stancode() and make_standata()
  2. Possible function rewrite.
  3. Modify some brms code for this specific purpose.

I will work through these suggestions. I am also excited to see the vignette that comes out on this topic.

Thanks again for the help. I am going to keep chipping away at this problem.


#18

You could think of adding another covariate to adstock which contains y[i - 1] and 0 for the first observation.