# Robyn quickstart on the Vexpower practice dataset
# ------------------------------------------------
# Runs Meta's Robyn (R package, v3.12.1 on CRAN as of October 2026) on
# mmm-practice-synthetic.csv: inputs, hyperparameter ranges, the Nevergrad
# search, the Pareto front, and a budget allocation for one chosen model.
#
# Written against the official demo (demo/demo.R in the Robyn repo) and the
# function signatures in the package source. It has not been run end to end by
# Vexpower because R wasn't available on the machine this was written on, so
# treat it as a starting point and check ?robyn_inputs etc. if anything errors.
#
# Setup, once:
#   install.packages("Robyn")
#   install.packages("reticulate")
#   reticulate::virtualenv_create("r-reticulate")
#   reticulate::use_virtualenv("r-reticulate", required = TRUE)
#   reticulate::py_install("nevergrad", pip = TRUE)
# See demo/install_nevergrad.R in the Robyn repo if Nevergrad won't load.
#
# Data: synthetic Vexpower teaching data, 104 weeks with columns week, revenue,
# search_spend, social_spend, video_spend, promotion. Be aware that the file was
# generated from an exactly linear formula with no noise:
#   revenue = 20000 + 1.4*search + 0.8*social + 0.6*video
#             + 3000*sin(pi*t/26) + 1800*promotion
# There is no adstock and no saturation in the truth, so a good run should find
# thetas near zero and push the Hill parameters to whichever end of their
# ranges makes the curves closest to straight. Fitted decay rates here are not
# findings; the file is for learning the workflow, not for testing Robyn.

library(Robyn)
data("dt_prophet_holidays")

Sys.setenv(R_FUTURE_FORK_ENABLE = "true")
options(future.fork.enable = TRUE)

csv_path <- "mmm-practice-synthetic.csv"
out_dir <- file.path(getwd(), "robyn-output")
dir.create(out_dir, showWarnings = FALSE)

dt <- read.csv(csv_path)
dt$week <- as.Date(dt$week)

# 1. Inputs ---------------------------------------------------------------
InputCollect <- robyn_inputs(
  dt_input = dt,
  dt_holidays = dt_prophet_holidays,
  date_var = "week",
  dep_var = "revenue",
  dep_var_type = "revenue",
  prophet_vars = c("trend", "season"), # no "holiday", so no country needed
  context_vars = c("promotion"),
  paid_media_spends = c("search_spend", "social_spend", "video_spend"),
  paid_media_vars = c("search_spend", "social_spend", "video_spend"),
  window_start = "2024-01-01",
  window_end = max(as.character(dt$week)),
  adstock = "geometric"
)

# 2. Hyperparameter ranges -------------------------------------------------
# Names must match hyper_names(); run this to check them.
print(hyper_names(adstock = InputCollect$adstock, all_media = InputCollect$all_media))

hyperparameters <- list(
  search_spend_alphas = c(0.5, 3),
  search_spend_gammas = c(0.3, 1),
  search_spend_thetas = c(0, 0.3), # digital: little carryover
  social_spend_alphas = c(0.5, 3),
  social_spend_gammas = c(0.3, 1),
  social_spend_thetas = c(0, 0.3),
  video_spend_alphas = c(0.5, 3),
  video_spend_gammas = c(0.3, 1),
  video_spend_thetas = c(0.1, 0.6), # video: allow more carryover
  train_size = c(0.5, 0.8)
)
InputCollect <- robyn_inputs(InputCollect = InputCollect, hyperparameters = hyperparameters)
print(InputCollect)

# Optional: calibrate with an experiment (point estimate of incremental revenue
# for the spend in the test window). Uncomment and fill in with a real result.
# calibration_input <- data.frame(
#   channel = "social_spend",
#   liftStartDate = as.Date("2025-03-03"),
#   liftEndDate = as.Date("2025-03-30"),
#   liftAbs = 9000,
#   spend = 7000,
#   confidence = 0.9,
#   metric = "revenue",
#   calibration_scope = "immediate"
# )
# InputCollect <- robyn_inputs(InputCollect = InputCollect, calibration_input = calibration_input)

# 3. Run the search -------------------------------------------------------
# 2,000 iterations x 5 trials is Robyn's recommendation for the demo data.
OutputModels <- robyn_run(
  InputCollect = InputCollect,
  cores = NULL, # all cores but one
  iterations = 2000,
  trials = 5,
  ts_validation = TRUE, # train / validation / test split on time
  add_penalty_factor = FALSE
)
print(OutputModels)
OutputModels$convergence$moo_distrb_plot
OutputModels$convergence$moo_cloud_plot

# 4. Pareto front, clusters and one-pagers ----------------------------------
OutputCollect <- robyn_outputs(
  InputCollect, OutputModels,
  pareto_fronts = "auto",
  csv_out = "pareto",
  clusters = TRUE,
  export = TRUE,
  plot_folder = out_dir,
  plot_pareto = TRUE
)
print(OutputCollect)

# The best model in each cluster. Open their one-pagers in out_dir, compare
# decomposition, ROAS and response curves, and pick the one that matches what
# you know about the business. Don't just take the first row.
print(OutputCollect$clusters$models)
select_model <- OutputCollect$clusters$models$solID[1] # replace with your pick

ExportedModel <- robyn_write(InputCollect, OutputCollect, select_model, export = TRUE)
print(ExportedModel)

# 5. Budget allocation ------------------------------------------------------
# Keep each channel between 70 and 150 percent of its current spend.
AllocatorCollect <- robyn_allocator(
  InputCollect = InputCollect,
  OutputCollect = OutputCollect,
  select_model = select_model,
  date_range = "last_13", # the last 13 weeks
  channel_constr_low = 0.7,
  channel_constr_up = 1.5,
  scenario = "max_response",
  export = TRUE
)
print(AllocatorCollect)
plot(AllocatorCollect)

# 6. Marginal return of the next 1,000 on search (metric_value is the total
# spend over date_range; the practice data spends about 37,000 on search in its
# last 13 weeks)
r1 <- robyn_response(
  InputCollect = InputCollect, OutputCollect = OutputCollect,
  select_model = select_model, metric_name = "search_spend",
  metric_value = 37000, date_range = "last_13"
)
r2 <- robyn_response(
  InputCollect = InputCollect, OutputCollect = OutputCollect,
  select_model = select_model, metric_name = "search_spend",
  metric_value = 38000, date_range = "last_13"
)
(r2$sim_mean_response - r1$sim_mean_response) / (r2$sim_mean_spend - r1$sim_mean_spend)