Working adstock Finder

working_adstock_finder.Rmd

```{r}
# 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")))
}

my_adstock <- function(x, rate){
  y <- NULL
  my_vector <- vector("numeric")
  for(i in seq(x)) {
    
       # Adstock: y[i] = x[i] + f[1]*y[i-1] + … + f[p]*y[i-p]
       y[1] <- x[i]
       y[i + 1] <- x[i + 1] + rate * y[i + 1 - 1]
       # store as vector
       my_vector[i] <- y[i]
  }
  # return results of vector
  my_vector
}

# 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))

fake_df
```

```{r}
#multivariate adstock function
AdstockRateMV <- function(Impact, Ads, maxiter = 100){
  # parameters
  params = letters[2:(ncol(Ads)+1)]
  # ad variables
  ads = paste0("ad_", params)
  # rate variables
  rates = paste0("rate_", params)
  # create partial formula
  param_fm = paste(
    paste(params, "*adstock(", ads,  ",", rates, ")", sep = ""),
    collapse = " + "
  )
  # create whole formula
  fm = as.formula(paste("Impact ~ a +", param_fm))
  # starting values for nls
  start = c(rep(1, length(params) + 1), rep(.1, length(rates)))
  names(start) = c("a", params, rates)
  # input data
  Ads_df = Ads
  names(Ads_df) = ads
  Data = cbind(Impact, Ads_df)
  # fit model
{
    library(minpack.lm)
    lower = c(rep(-Inf, length(params) + 1), rep(0, length(rates)))
    upper = c(rep(Inf, length(params) + 1), rep(1, length(rates)))
    modFit <- nlsLM(fm, data = Data, start = start,
                   lower = lower, upper = upper,
                   control = nls.lm.control(maxiter = maxiter))
  }
  # model coefficients
  AdstockInt = round(summary(modFit)$coefficients[1, 1])
  AdstockCoef = round(summary(modFit)$coefficients[params, 1], 2)
  AdstockRate = round(summary(modFit)$coefficients[rates, 1], 2)
  # print formula with coefficients
  param_fm_coefs = paste(
    paste(round(AdstockCoef, 2), " * adstock(", names(Ads),  ", ", round(AdstockRate, 2), ")", sep = ""),
    collapse = " + "
  )
  fm_coefs = as.formula(paste("Impact ~ ", AdstockInt, " +", param_fm_coefs))
  # rename rates with original variable names
  names(AdstockRate) = paste0("rate_", names(Ads))
  # calculate percent error
  mape = mean(abs((Impact-predict(modFit))/Impact) * 100)
  # return outputs
  return(list(fm = fm_coefs, base = AdstockInt, rates = AdstockRate, mape = mape))
}

Impact = fake_df$sales
Ads = fake_df %>% select(ad1, ad2, ad3 )
(mod = AdstockRateMV(Impact, Ads))
```