Your First BLM Model

Goal

Build, fit, and interpret a single-market Bayesian linear model (BLM) using the DSAMbayes R API.

This page is a hands-on tutorial. For the broader methodological questions behind prior-setting and diagnostics, use the workflow pages:

Prerequisites

  • DSAMbayes installed locally (see Install and Setup).
  • Familiarity with R and lm()-style formulas.

Dataset

This walkthrough uses the tracked synthetic dataset at data/timeseries/demo_data_synthetic.csv. It contains weekly observations for a single series with columns for:

  • Response: revenue, weekly synthetic KPI.
  • Media: channel0_spend, channel1_spend, channel2_spend, channel3_spend.
  • Controls: events, newsletters.
  • Date: time, weekly date index.
library(DSAMbayes)
df <- read.csv("data/timeseries/demo_data_synthetic.csv")
str(df)

Step 1: Construct the model

blm() creates an unfitted model object. No fitting happens yet.

model <- blm(
  revenue ~
    events + newsletters +
    channel0_spend + channel1_spend + channel2_spend + channel3_spend,
  data = df
) %>%
  set_date(time)

Inspect defaults:

peek_prior(model)      # data-dependent normal priors: normal(ybar, sy) intercept, normal(0, sy / sx) per slope
peek_boundary(model)   # (-Inf, Inf), unconstrained

Step 2: Set boundaries

Media channels should have non-negative effects. Use inequality notation:

model <- model %>%
  set_boundary(
    channel0_spend > 0, channel1_spend > 0,
    channel2_spend > 0, channel3_spend > 0
  )

Step 3: (Optional) Override priors

Default priors are weakly informative. Override only with domain knowledge:

model <- model %>%
  set_prior(events ~ normal(0, 1))

See Stage 2: Model and Priors and Minimal-Prior Policy for guidance.

Step 4: Fit with MCMC

fitted_model <- model %>%
  fit(cores = 2, iter = 2000, warmup = 1000, seed = 123)

First-time Stan compilation takes 1–3 minutes. Subsequent runs use a cached binary. With the default 4 chains split across 2 cores, sampling on synthetic data typically completes in under 2 minutes.

Step 5: Sampler diagnostics

chain_diagnostics(fitted_model)
Metric Good Concern
Max Rhat <= 1.01 > 1.01 means chains have not converged
Min ESS (bulk) > 400 < 200 means too few effective samples
Divergences 0 Any non-zero count warrants investigation

These are the computational checks. They do not by themselves prove the model is adequate for interpretation.

Step 6: Extract the posterior

post <- get_posterior(fitted_model)

post is a tibble with one row per draw containing coef (named coefficient list), yhat (fitted values), noise_sd, r2, rmse, and smape.

Summarise coefficients:

library(dplyr); library(tidyr)

coef_summary <- post %>%
  select(coef) %>%
  unnest_wider(coef) %>%
  pivot_longer(everything(), names_to = "term") %>%
  group_by(term) %>%
  summarise(
    mean = mean(value), median = median(value), sd = sd(value),
    ci_low = quantile(value, 0.025), ci_high = quantile(value, 0.975),
    .groups = "drop"
  )
print(coef_summary, n = 30)

What to look for:

  • Media coefficients should be positive (boundaries enforce this).
  • Wide credible intervals mean the prior dominates, the data cannot identify the effect precisely.

Step 7: Assess model fit

fit_tbl <- fitted(fitted_model)
cat("Median R²:", median(r2(fitted_model)), "\n")
cat("Median RMSE:", median(rmse(fitted_model)), "\n")

For a well-specified MMM on weekly data, in-sample R² above 0.85 is typical.

Step 8: Response decomposition

tbl_map <- data.frame(
  actual.vars = c(
    "intercept", "events", "newsletters",
    "channel0_spend", "channel1_spend", "channel2_spend", "channel3_spend"
  ),
  group = c(rep("BASE", 3), rep("media", 4)),
  ref_points = c("none", rep("mean", 6))
)
decomp_result <- decomp(fitted_model, tbl_map = tbl_map)
head(decomp_result$DecompedData)

$DecompedData shows each term’s contribution (coefficient × design-matrix column) to the predicted KPI at each time point. The result is a native dsambayes_decomposition S3 object.

Step 9: MAP for rapid iteration

During development, use MAP for fast point estimates:

map_model <- model %>% fit_map(n_runs = 10)
get_posterior(map_model)

Use MCMC for final reporting; MAP for formula iteration. A MAP result is one point estimate, not posterior draws: it has no valid credible intervals or MCMC diagnostics. Review the restart diagnostics if objectives differ materially; use fit() when uncertainty, convergence, or decision risk matters.

Common pitfalls

Pitfall Symptom Fix
Forgetting to set boundaries Media coefficients go negative Add set_boundary(m_x > 0)
Too few iterations High Rhat, low ESS Increase iter and warmup
Missing controls High residual autocorrelation Add trend, seasonality, or holiday terms
Scaling confusion Coefficients look wrong model.scale: true is default; get_posterior() back-transforms automatically

Next steps