{ "cells": [ { "cell_type": "markdown", "id": "c1c6c532", "metadata": {}, "source": [ "# How much should I spend on each channel — and where do the extra euros stop paying off?\n", "\n", "**Marketing mix modelling (MMM), worked end to end on synthetic data you can swap for your own.**\n", "\n", "This notebook builds a small marketing-mix model from scratch, the way Google's open-source\n", "**Meridian** library does it under the hood: it learns a *diminishing-returns curve* for each channel\n", "from your weekly spend and revenue, reads the **return on the next euro** at your current spend, and\n", "turns that into a concrete budget move — with an honest uncertainty range on every number.\n", "\n", "Two parts:\n", "1. **Worked example** on a synthetic advertiser (clearly labelled — no real brand's data).\n", "2. **Run it on your own data** — paste a weekly CSV, and the same code runs, with plain-words\n", " warnings when your data is too thin to trust.\n", "\n", "Everything runs locally in this notebook. **Your data never leaves your machine.**\n", "\n", "> The maths here (Hill saturation, geometric adstock, the adstock-then-saturation order) follows\n", "> Jin et al. 2017 and the Meridian documentation. The *fit* below uses a light least-squares with\n", "> soft priors — a concept demonstration. The real library fits the same shapes with Bayesian MCMC,\n", "> which is more rigorous and is where the credible intervals come from." ] }, { "cell_type": "code", "execution_count": 1, "id": "781a0aaa", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:29.246784Z", "iopub.status.busy": "2026-07-23T13:15:29.246784Z", "iopub.status.idle": "2026-07-23T13:15:34.630559Z", "shell.execute_reply": "2026-07-23T13:15:34.630053Z" } }, "outputs": [], "source": [ "import json\n", "import warnings\n", "from pathlib import Path\n", "\n", "import matplotlib\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "import pandas as pd\n", "from scipy.optimize import least_squares, minimize\n", "\n", "matplotlib.use(\"Agg\")\n", "warnings.filterwarnings(\"ignore\")\n", "PAPER=\"#FDFCF9\"; INK=\"#1f2124\"; MUTED=\"#6b6f76\"\n", "SERIES={\"Search\":\"#2a78d6\",\"Social\":\"#d55181\",\"Video\":\"#0f8f5f\"}\n", "plt.rcParams.update({\n", " \"figure.facecolor\":PAPER,\"axes.facecolor\":PAPER,\"savefig.facecolor\":PAPER,\n", " \"axes.edgecolor\":\"#d9d5cc\",\"axes.labelcolor\":INK,\"text.color\":INK,\n", " \"xtick.color\":MUTED,\"ytick.color\":MUTED,\"axes.grid\":True,\"grid.color\":\"#ece8df\",\n", " \"font.size\":11,\"axes.titlesize\":13,\"axes.titleweight\":\"bold\",\"figure.dpi\":120,\n", "})\n", "def _eur(x, p=None): return f\"EUR {x:,.0f}\"\n", "def _eurk(x, p=None): return \"EUR 0\" if abs(x)<1 else f\"EUR {x/1000:.0f}k\"\n", "def money_x(ax):\n", " ax.xaxis.set_major_formatter(_eurk); ax.xaxis.set_major_locator(plt.MaxNLocator(6))\n", "def money_y(ax):\n", " ax.yaxis.set_major_formatter(_eurk); ax.yaxis.set_major_locator(plt.MaxNLocator(6))\n", "FIG=Path(\"figures\"); FIG.mkdir(exist_ok=True)\n", "DATA=Path(\"data\"); DATA.mkdir(exist_ok=True)\n", "\n", "SEED = 20260723\n", "MAX_LAG = 8 # Meridian ModelSpec default\n", "\n", "\n", "# ---------------------------------------------------------------- transforms\n", "def hill(x, ec, slope):\n", " x = np.asarray(x, float)\n", " xs = np.power(np.clip(x, 0, None), slope)\n", " return xs / (xs + ec ** slope)\n", "\n", "\n", "def hill_deriv(x, ec, slope):\n", " \"\"\"d Hill/dx, closed form. Revenue-per-euro = beta * hill_deriv.\"\"\"\n", " x = np.asarray(x, float)\n", " x = np.clip(x, 1e-9, None)\n", " xs = np.power(x, slope)\n", " ecs = ec ** slope\n", " return slope * np.power(x, slope - 1.0) * ecs / (xs + ecs) ** 2\n", "\n", "\n", "def geometric_adstock(x, alpha, L=MAX_LAG):\n", " x = np.asarray(x, float)\n", " w = alpha ** np.arange(L + 1)\n", " w = w / w.sum()\n", " out = np.zeros_like(x)\n", " for s in range(L + 1):\n", " out[s:] += w[s] * x[: len(x) - s]\n", " return out\n", "\n", "\n", "def steady_response(x, beta, ec, slope):\n", " return beta * hill(x, ec, slope)\n", "\n", "\n", "def marginal_return(x, beta, ec, slope):\n", " return beta * hill_deriv(x, ec, slope)\n", "\n", "\n", "# ---------------------------------------------------------------- simulation\n", "CHANNELS = [\"Search\", \"Social\", \"Video\"]\n", "\n", "# Ground-truth generating parameters (tuned so the story is crisp + the fit is identifiable):\n", "# Social = over-invested (heavy spend, flat margin); Video = under-invested (steep margin);\n", "# Search = efficient middle. beta = max weekly incremental revenue (EUR); ec on weekly-spend scale.\n", "# Marginals well separated (Video >> Search > Social) so the ordering is robust to noise.\n", "TRUTH = {\n", " \"Search\": dict(beta=40000.0, ec=6500.0, slope=1.0, alpha=0.50, mean=9000.0, sd=2700.0),\n", " \"Social\": dict(beta=52000.0, ec=4500.0, slope=1.3, alpha=0.60, mean=14000.0, sd=4200.0),\n", " \"Video\": dict(beta=45000.0, ec=7000.0, slope=1.0, alpha=0.40, mean=3500.0, sd=1100.0),\n", "}\n", "BASE_INTERCEPT = 95000.0\n", "BASE_TREND = 80.0 # EUR/week drift\n", "SEASON_AMP = 12000.0 # annual seasonality amplitude\n", "NOISE_SD_FRAC = 0.030 # ~3% of revenue\n", "\n", "\n", "def _seasonality(t):\n", " ph = 2 * np.pi * t / 52.0\n", " return SEASON_AMP * np.sin(ph) + 0.4 * SEASON_AMP * np.sin(2 * ph + 0.7)\n", "\n", "\n", "def simulate_panel(weeks=156, seed=SEED, collinear=False):\n", " rng = np.random.default_rng(seed)\n", " t = np.arange(weeks)\n", " spends = {}\n", " base_walk = None\n", " for i, ch in enumerate(CHANNELS):\n", " p = TRUTH[ch]\n", " # AR(1)-ish independent variation around the mean, clipped positive\n", " noise = rng.normal(0, p[\"sd\"], weeks)\n", " walk = np.zeros(weeks)\n", " walk[0] = noise[0]\n", " for k in range(1, weeks):\n", " walk[k] = 0.6 * walk[k - 1] + noise[k]\n", " if collinear and ch in (\"Search\", \"Social\"):\n", " # force Search + Social to move together (breaks identifiability)\n", " if base_walk is None:\n", " base_walk = walk\n", " walk = 0.85 * base_walk + 0.15 * walk\n", " spends[ch] = np.clip(p[\"mean\"] + walk, 0.15 * p[\"mean\"], None)\n", "\n", " revenue = BASE_INTERCEPT + BASE_TREND * t + _seasonality(t)\n", " contrib = {}\n", " for ch in CHANNELS:\n", " p = TRUTH[ch]\n", " ad = geometric_adstock(spends[ch], p[\"alpha\"])\n", " c = p[\"beta\"] * hill(ad, p[\"ec\"], p[\"slope\"])\n", " contrib[ch] = c\n", " revenue = revenue + c\n", " revenue = revenue + rng.normal(0, NOISE_SD_FRAC * revenue.mean(), weeks)\n", "\n", " df = pd.DataFrame({\"week\": t, \"revenue\": revenue})\n", " for ch in CHANNELS:\n", " df[f\"spend_{ch}\"] = spends[ch]\n", " return df, contrib\n", "\n", "\n", "# ---------------------------------------------------------------- the fit\n", "def _design_baseline(weeks):\n", " t = np.arange(weeks)\n", " ph = 2 * np.pi * t / 52.0\n", " return np.column_stack([\n", " np.ones(weeks), t / weeks,\n", " np.sin(ph), np.cos(ph), np.sin(2 * ph), np.cos(2 * ph),\n", " ])\n", "\n", "\n", "def _predict(params, spends, Z):\n", " b = params[:Z.shape[1]]\n", " yhat = Z @ b\n", " off = Z.shape[1]\n", " per = {}\n", " for i, ch in enumerate(CHANNELS):\n", " beta, ec, slope, alpha = params[off + 4 * i: off + 4 * i + 4]\n", " ad = geometric_adstock(spends[ch], alpha)\n", " c = beta * hill(ad, ec, slope)\n", " per[ch] = c\n", " yhat = yhat + c\n", " return yhat, per\n", "\n", "\n", "# Soft priors (what Bayesian MMM does with priors; here a light penalty that regularises the\n", "# ill-identified shape/carryover directions so the curve near observed spend is pinned down).\n", "# ec prior is per channel = the channel's mean spend (half-saturation near typical spend).\n", "PRIOR_SLOPE, PRIOR_SLOPE_SD = 1.0, 0.35\n", "PRIOR_ALPHA, PRIOR_ALPHA_SD = 0.45, 0.20\n", "PRIOR_LOG_EC_SD = 0.60 # log-scale sd around log(mean spend)\n", "\n", "\n", "def _bounds(spends, y, off):\n", " lo = [-np.inf] * off\n", " hi = [np.inf] * off\n", " for ch in CHANNELS:\n", " mx = spends[ch].max()\n", " lo += [0.0, 0.05 * mx, 0.3, 0.0]\n", " hi += [3.0 * y.mean(), 4.0 * mx, 3.0, 0.95]\n", " return lo, hi\n", "\n", "\n", "def _make_resid(spends, Z, y_target, ec_prior, data_scale):\n", " off = Z.shape[1]\n", "\n", " def resid(pp):\n", " yhat, _ = _predict(pp, spends, Z)\n", " r = list((yhat - y_target) / data_scale)\n", " for i, ch in enumerate(CHANNELS):\n", " _, ec, slope, alpha = pp[off + 4 * i: off + 4 * i + 4]\n", " r.append((slope - PRIOR_SLOPE) / PRIOR_SLOPE_SD)\n", " r.append((alpha - PRIOR_ALPHA) / PRIOR_ALPHA_SD)\n", " r.append((np.log(max(ec, 1e-6)) - np.log(ec_prior[ch])) / PRIOR_LOG_EC_SD)\n", " return np.asarray(r)\n", "\n", " return resid\n", "\n", "\n", "def _unpack(x, off):\n", " rec = {}\n", " for i, ch in enumerate(CHANNELS):\n", " beta, ec, slope, alpha = x[off + 4 * i: off + 4 * i + 4]\n", " rec[ch] = dict(beta=beta, ec=ec, slope=slope, alpha=alpha)\n", " return rec\n", "\n", "\n", "def fit(df, seed=SEED, n_starts=6):\n", " weeks = len(df)\n", " spends = {ch: df[f\"spend_{ch}\"].to_numpy() for ch in CHANNELS}\n", " Z = _design_baseline(weeks)\n", " y = df[\"revenue\"].to_numpy()\n", " off = Z.shape[1]\n", " ec_prior = {ch: float(np.mean(spends[ch])) for ch in CHANNELS}\n", " data_scale = NOISE_SD_FRAC * y.mean()\n", " resid = _make_resid(spends, Z, y, ec_prior, data_scale)\n", " lo, hi = _bounds(spends, y, off)\n", " b0 = list(np.linalg.lstsq(Z, y, rcond=None)[0])\n", "\n", " rng = np.random.default_rng(seed)\n", " best, best_cost = None, np.inf\n", " for s in range(n_starts):\n", " p0 = list(b0)\n", " for ch in CHANNELS:\n", " j = 1.0 if s == 0 else float(rng.uniform(0.5, 1.6))\n", " e = 1.0 if s == 0 else float(rng.uniform(0.6, 1.6))\n", " p0 += [0.3 * y.mean() * j, ec_prior[ch] * e, 1.0, 0.45]\n", " try:\n", " sol = least_squares(resid, p0, bounds=(lo, hi), method=\"trf\", max_nfev=6000)\n", " except Exception:\n", " continue\n", " if sol.cost < best_cost:\n", " best, best_cost = sol, sol.cost\n", " return _unpack(best.x, off), best, spends, Z, y\n", "\n", "\n", "def bootstrap(df, sol, spends, Z, y, B=120, block=MAX_LAG, seed=SEED):\n", " \"\"\"Block residual bootstrap: resample residual blocks, refit from the point estimate.\"\"\"\n", " weeks = len(y)\n", " off = Z.shape[1]\n", " ec_prior = {ch: float(np.mean(spends[ch])) for ch in CHANNELS}\n", " data_scale = NOISE_SD_FRAC * y.mean()\n", " lo, hi = _bounds(spends, y, off)\n", " yhat, _ = _predict(sol.x, spends, Z)\n", " resids = y - yhat\n", " n_blocks = int(np.ceil(weeks / block))\n", " starts = np.arange(weeks - block + 1)\n", " rng = np.random.default_rng(seed + 777)\n", "\n", " cur = {c: float(df[f\"spend_{c}\"].mean()) for c in CHANNELS}\n", " omax = {c: float(df[f\"spend_{c}\"].max()) for c in CHANNELS}\n", " out = {\"marginal\": {c: [] for c in CHANNELS}, \"gain\": [], \"gain_pct\": [],\n", " \"new_spend\": {c: [] for c in CHANNELS}}\n", " for _ in range(B):\n", " chosen = rng.choice(starts, size=n_blocks, replace=True)\n", " e = np.concatenate([resids[s:s + block] for s in chosen])[:weeks]\n", " y_star = yhat + e\n", " resid = _make_resid(spends, Z, y_star, ec_prior, data_scale)\n", " try:\n", " s = least_squares(resid, sol.x, bounds=(lo, hi), method=\"trf\", max_nfev=4000)\n", " except Exception:\n", " continue\n", " rec = _unpack(s.x, off)\n", " for c in CHANNELS:\n", " out[\"marginal\"][c].append(float(marginal_return(cur[c], rec[c][\"beta\"], rec[c][\"ec\"], rec[c][\"slope\"])))\n", " new, gain, orev, nrev = optimal_allocation(rec, cur, omax)\n", " out[\"gain\"].append(gain)\n", " out[\"gain_pct\"].append(100 * gain / orev)\n", " for c in CHANNELS:\n", " out[\"new_spend\"][c].append(new[c])\n", " return out\n", "\n", "\n", "def ci(a, lo=5, hi=95):\n", " a = np.asarray(a, float)\n", " return [round(float(np.percentile(a, lo)), 2), round(float(np.percentile(a, hi)), 2)]\n", "\n", "\n", "# ---------------------------------------------------------------- allocation\n", "def optimal_allocation(rec, current_spend, observed_max, cap_mult=2.0):\n", " \"\"\"Maximise total steady response s.t. fixed total budget; cap each at cap_mult*observed_max.\"\"\"\n", " chs = CHANNELS\n", " B = sum(current_spend[c] for c in chs)\n", " x0 = np.array([current_spend[c] for c in chs])\n", " ub = np.array([cap_mult * observed_max[c] for c in chs])\n", "\n", " def negtot(x):\n", " return -sum(steady_response(x[i], rec[c][\"beta\"], rec[c][\"ec\"], rec[c][\"slope\"])\n", " for i, c in enumerate(chs))\n", "\n", " cons = [{\"type\": \"eq\", \"fun\": lambda x: x.sum() - B}]\n", " bnds = [(0.0, ub[i]) for i in range(len(chs))]\n", " sol = minimize(negtot, x0, method=\"SLSQP\", bounds=bnds, constraints=cons,\n", " options={\"maxiter\": 500, \"ftol\": 1e-9})\n", " new = {c: float(sol.x[i]) for i, c in enumerate(chs)}\n", " old_rev = sum(steady_response(current_spend[c], rec[c][\"beta\"], rec[c][\"ec\"], rec[c][\"slope\"]) for c in chs)\n", " new_rev = sum(steady_response(new[c], rec[c][\"beta\"], rec[c][\"ec\"], rec[c][\"slope\"]) for c in chs)\n", " return new, float(new_rev - old_rev), float(old_rev), float(new_rev)\n", "\n", "\n", "# ---------------------------------------------------------------- diagnostics\n", "def sufficiency(df, n_controls=2, n_knots=4):\n", " weeks = len(df)\n", " n_ch = len(CHANNELS)\n", " dpe = weeks / (n_ch + n_controls + n_knots)\n", " return dict(weeks=weeks, points_per_effect=round(dpe, 2),\n", " time_floor_ok=weeks >= 52, sufficiency_ok=dpe >= 10)\n", "\n", "\n", "def collinearity(df):\n", " S = df[[f\"spend_{c}\" for c in CHANNELS]].to_numpy()\n", " C = np.corrcoef(S.T)\n", " out = []\n", " for i in range(len(CHANNELS)):\n", " for j in range(i + 1, len(CHANNELS)):\n", " out.append((CHANNELS[i], CHANNELS[j], round(float(C[i, j]), 3)))\n", " return out\n", "\n", "\n", "def low_spend(df):\n", " tot = sum(df[f\"spend_{c}\"].sum() for c in CHANNELS)\n", " return {c: round(float(df[f\"spend_{c}\"].sum() / tot), 3) for c in CHANNELS}\n", "\n", "\n", "# ---------------------------------------------------------------- run / tune" ] }, { "cell_type": "markdown", "id": "3ec89777", "metadata": {}, "source": [ "## 1. The two shapes every channel has\n", "\n", "**Diminishing returns (saturation).** The first euro on a channel usually works harder than the last.\n", "We model each channel's response with a **Hill curve**: `Hill(x) = x^slope / (x^slope + ec^slope)`,\n", "which rises from 0 and flattens toward 1. `ec` is the **half-saturation point** — the weekly spend at\n", "which the channel delivers half of everything it can. `slope` sets how sharp the bend is. Scaled by\n", "`beta` (the channel's ceiling in euros of revenue), the curve is `beta * Hill(x)`.\n", "\n", "**Carryover (adstock).** An ad seen this week still sells a little next week. We spread each week's\n", "spend over the following weeks with **geometric adstock** — a normalized weighted average with decay\n", "`alpha`. Because the weights are normalized, at a *steady* weekly spend the carryover washes out and\n", "the response is just `beta * Hill(spend)`; adstock only matters while spend is changing." ] }, { "cell_type": "code", "execution_count": 2, "id": "db8a302f", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:34.635566Z", "iopub.status.busy": "2026-07-23T13:15:34.634566Z", "iopub.status.idle": "2026-07-23T13:15:35.707142Z", "shell.execute_reply": "2026-07-23T13:15:35.706792Z" } }, "outputs": [], "source": [ "fig, ax = plt.subplots(1, 2, figsize=(11, 4.2))\n", "x = np.linspace(0, 24000, 400)\n", "for ch in CHANNELS:\n", " p = TRUTH[ch]\n", " ax[0].plot(x, steady_response(x, p[\"beta\"], p[\"ec\"], p[\"slope\"]), color=SERIES[ch], lw=2.4, label=ch)\n", " ax[0].axvline(p[\"ec\"], color=SERIES[ch], ls=\":\", lw=1, alpha=.6)\n", "ax[0].set_title(\"Saturation: revenue vs weekly spend\"); ax[0].set_xlabel(\"weekly spend\"); ax[0].set_ylabel(\"weekly revenue from channel\")\n", "money_x(ax[0]); money_y(ax[0]); ax[0].legend(frameon=False)\n", "spike = np.zeros(16); spike[2] = 10000.0\n", "for a, lab in [(0.4, \"fast decay (alpha=0.4)\"), (0.7, \"slow decay (alpha=0.7)\")]:\n", " ax[1].plot(range(16), geometric_adstock(spike, a), lw=2.2, label=lab)\n", "ax[1].bar(range(16), spike, color=\"#d9d5cc\", label=\"raw spend (one week)\")\n", "ax[1].set_title(\"Carryover: a single week of spend, spread forward\"); ax[1].set_xlabel(\"week\"); ax[1].set_ylabel(\"adstocked exposure\")\n", "money_y(ax[1]); ax[1].legend(frameon=False)\n", "plt.tight_layout(); plt.savefig(FIG/\"fig_shapes.png\", bbox_inches=\"tight\"); plt.show()" ] }, { "cell_type": "markdown", "id": "f49e1513", "metadata": {}, "source": [ "## 2. A realistic advertiser (synthetic)\n", "\n", "We invent two years (156 weeks) of weekly data for an ecommerce brand running three channels —\n", "**Search**, **Social**, and **Video** — plus a baseline (a trend and a yearly season) and noise.\n", "By construction, **Social is over-invested** (a lot of spend, sitting flat at the top of its curve),\n", "**Video is under-invested** (little spend, still on the steep part), and **Search is about right**.\n", "The model does not get told any of this — it has to recover it from the numbers." ] }, { "cell_type": "code", "execution_count": 3, "id": "df36f8c2", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:35.712144Z", "iopub.status.busy": "2026-07-23T13:15:35.711142Z", "iopub.status.idle": "2026-07-23T13:15:36.378919Z", "shell.execute_reply": "2026-07-23T13:15:36.377913Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "weeks: 156 | mean weekly revenue: EUR 179,159\n", "mean weekly media contribution: EUR 78,654 (44% of revenue)\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
weekrevenuespend_Searchspend_Socialspend_Video
00148391.08828.013008.02973.0
11159819.07230.010991.04416.0
22174379.03206.016726.03873.0
33184186.07298.018012.04095.0
\n", "
" ], "text/plain": [ " week revenue spend_Search spend_Social spend_Video\n", "0 0 148391.0 8828.0 13008.0 2973.0\n", "1 1 159819.0 7230.0 10991.0 4416.0\n", "2 2 174379.0 3206.0 16726.0 3873.0\n", "3 3 184186.0 7298.0 18012.0 4095.0" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "df, contrib = simulate_panel()\n", "print(\"weeks:\", len(df), \"| mean weekly revenue:\", _eur(df.revenue.mean()))\n", "print(\"mean weekly media contribution:\", _eur(sum(contrib[c].mean() for c in CHANNELS)),\n", " f\"({sum(contrib[c].mean() for c in CHANNELS)/df.revenue.mean():.0%} of revenue)\")\n", "display(df.head(4).round(0))\n", "\n", "fig, ax = plt.subplots(figsize=(11, 4.2))\n", "t = df.week\n", "base = df.revenue.to_numpy() - sum(contrib[c] for c in CHANNELS)\n", "ax.fill_between(t, 0, base, color=\"#e7e2d6\", label=\"baseline (trend + season)\")\n", "run = base.copy()\n", "for ch in CHANNELS:\n", " ax.fill_between(t, run, run + contrib[ch], color=SERIES[ch], alpha=.85, label=ch)\n", " run = run + contrib[ch]\n", "ax.plot(t, df.revenue, color=INK, lw=.8, alpha=.5)\n", "ax.set_title(\"Where the revenue comes from, week by week (synthetic)\")\n", "ax.set_xlabel(\"week\"); ax.set_ylabel(\"weekly revenue\"); money_y(ax)\n", "ax.legend(frameon=False, ncol=4, loc=\"upper left\"); ax.margins(x=0)\n", "plt.tight_layout(); plt.savefig(FIG/\"fig_revenue.png\", bbox_inches=\"tight\"); plt.show()" ] }, { "cell_type": "markdown", "id": "d1a85b33", "metadata": {}, "source": [ "## 3. Fit the model — with priors, the way the real thing does\n", "\n", "Fitting an MMM is genuinely hard: many different curves fit the same history almost equally well.\n", "That is why real MMMs (Meridian included) lean on **priors** — gentle nudges toward sensible shapes.\n", "Here we do a light version: a least-squares fit that pays a small penalty for implausible saturation\n", "and carryover values. We report **wMAPE** (weighted average error), the metric Meridian itself reports.\n", "\n", "> We deliberately don't lead with R². Chan & Perry (2017) show that a great fit does **not** mean a\n", "> trustworthy budget answer — several models can fit within a few percent yet disagree by 50% on what\n", "> to do. Good fit is necessary, not sufficient." ] }, { "cell_type": "code", "execution_count": 4, "id": "12e6080c", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:36.382925Z", "iopub.status.busy": "2026-07-23T13:15:36.381919Z", "iopub.status.idle": "2026-07-23T13:15:37.972792Z", "shell.execute_reply": "2026-07-23T13:15:37.971785Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "weighted error (wMAPE): 2.0%\n" ] } ], "source": [ "rec, sol, spends, Z, y = fit(df)\n", "yhat, _ = _predict(sol.x, spends, Z)\n", "wmape = float(np.abs(y - yhat).sum() / y.sum())\n", "print(f\"weighted error (wMAPE): {wmape:.1%}\")\n", "\n", "fig, ax = plt.subplots(figsize=(11, 3.8))\n", "ax.plot(df.week, y, color=INK, lw=1.6, label=\"actual revenue\")\n", "ax.plot(df.week, yhat, color=\"#0f6e56\", lw=1.8, ls=\"--\", label=\"model fit\")\n", "ax.set_title(f\"Actual vs fitted weekly revenue (wMAPE {wmape:.1%})\")\n", "ax.set_xlabel(\"week\"); ax.set_ylabel(\"weekly revenue\"); money_y(ax)\n", "ax.legend(frameon=False); ax.margins(x=0)\n", "plt.tight_layout(); plt.savefig(FIG/\"fig_fit.png\", bbox_inches=\"tight\"); plt.show()" ] }, { "cell_type": "markdown", "id": "e4abe655", "metadata": {}, "source": [ "## 4. The number that matters: the return on the next euro\n", "\n", "For a budget decision you don't need the whole curve — you need its **slope at where you're spending\n", "now**: the extra revenue the *next* euro on each channel would bring (the marginal return, mROI).\n", "Move money from channels with a low next-euro return to channels with a high one, and total revenue\n", "rises **without spending more**. That is the entire idea." ] }, { "cell_type": "code", "execution_count": 5, "id": "fbb83a2f", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:37.977349Z", "iopub.status.busy": "2026-07-23T13:15:37.976349Z", "iopub.status.idle": "2026-07-23T13:15:38.816694Z", "shell.execute_reply": "2026-07-23T13:15:38.815687Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " Search: next euro returns EUR 1.24 (spending EUR 8,899/wk now)\n", " Social: next euro returns EUR 0.99 (spending EUR 13,404/wk now)\n", " Video: next euro returns EUR 2.69 (spending EUR 3,793/wk now)\n", "\n", "Rebalance (same total budget):\n", " Search: EUR 8,899 -> EUR 8,719 (-EUR 180/wk)\n", " Social: EUR 13,404 -> EUR 10,723 (-EUR 2,681/wk)\n", " Video: EUR 3,793 -> EUR 6,654 (+EUR 2,861/wk)\n", "\n", "Weekly media revenue: EUR 78,119 -> EUR 80,180 (+EUR 2,061/wk, +2.6%)\n", "About EUR 107,184 a year, at no extra spend.\n" ] } ], "source": [ "cur = {c: float(df[f\"spend_{c}\"].mean()) for c in CHANNELS}\n", "omax = {c: float(df[f\"spend_{c}\"].max()) for c in CHANNELS}\n", "mar = {c: float(marginal_return(cur[c], rec[c][\"beta\"], rec[c][\"ec\"], rec[c][\"slope\"])) for c in CHANNELS}\n", "for c in CHANNELS:\n", " print(f\"{c:>7}: next euro returns EUR {mar[c]:.2f} (spending {_eur(cur[c])}/wk now)\")\n", "\n", "new, gain, old_rev, new_rev = optimal_allocation(rec, cur, omax)\n", "print(\"\\nRebalance (same total budget):\")\n", "for c in CHANNELS:\n", " d = new[c] - cur[c]\n", " print(f\" {c:>7}: {_eur(cur[c])} -> {_eur(new[c])} ({'+' if d>=0 else '-'}{_eur(abs(d))}/wk)\")\n", "print(f\"\\nWeekly media revenue: {_eur(old_rev)} -> {_eur(new_rev)} (+{_eur(gain)}/wk, {gain/old_rev:+.1%})\")\n", "print(f\"About {_eur(gain*52)} a year, at no extra spend.\")\n", "\n", "fig, ax = plt.subplots(1, 2, figsize=(11, 4.3))\n", "xg = np.linspace(0, 2*max(omax.values()), 400)\n", "for ch in CHANNELS:\n", " p = rec[ch]\n", " ax[0].plot(xg, steady_response(xg, p[\"beta\"], p[\"ec\"], p[\"slope\"]), color=SERIES[ch], lw=2.2, label=ch)\n", " ax[0].scatter([cur[ch]], [steady_response(cur[ch], p[\"beta\"], p[\"ec\"], p[\"slope\"])], color=SERIES[ch], s=55, zorder=5)\n", "ax[0].set_title(\"Fitted response curves (dot = current spend)\"); ax[0].set_xlabel(\"weekly spend\")\n", "ax[0].set_ylabel(\"weekly revenue\"); money_x(ax[0]); money_y(ax[0]); ax[0].legend(frameon=False)\n", "chs = CHANNELS; xpos = np.arange(len(chs))\n", "ax[1].bar(xpos-0.2, [cur[c] for c in chs], .4, color=\"#cfc9bd\", label=\"now\")\n", "ax[1].bar(xpos+0.2, [new[c] for c in chs], .4, color=[SERIES[c] for c in chs], label=\"rebalanced\")\n", "ax[1].set_xticks(xpos); ax[1].set_xticklabels(chs); ax[1].set_title(\"Recommended weekly spend\"); money_y(ax[1]); ax[1].legend(frameon=False)\n", "plt.tight_layout(); plt.savefig(FIG/\"fig_reallocation.png\", bbox_inches=\"tight\"); plt.show()" ] }, { "cell_type": "markdown", "id": "8c3523fd", "metadata": {}, "source": [ "## 5. How sure are we? (uncertainty is a feature, not a footnote)\n", "\n", "Every number above is an estimate from noisy data. We resample the model's residuals many times and\n", "refit — a **bootstrap** — to get a 90% range on each marginal return and on the revenue gain. Notice\n", "two things: the ranges are wide (MMM is not a precise instrument), and they get **wider the further we\n", "push past spend levels the advertiser has actually run** — because out there the model is guessing." ] }, { "cell_type": "code", "execution_count": 6, "id": "94696e00", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:38.828199Z", "iopub.status.busy": "2026-07-23T13:15:38.827199Z", "iopub.status.idle": "2026-07-23T13:15:58.615104Z", "shell.execute_reply": "2026-07-23T13:15:58.614596Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Return on the next euro (90% range):\n", " Search: EUR 1.27 [1.03, 1.50]\n", " Social: EUR 1.01 [0.86, 1.17]\n", " Video: EUR 2.70 [1.93, 3.39]\n", "\n", "Revenue gain from rebalancing: +EUR 2,061/wk 90% range [EUR 1,224, EUR 4,055]/wk\n", "Recommended Video spend EUR 6,654/wk sits near the top of its observed range (EUR 7,589/wk);\n", "push much past that and the estimate turns into extrapolation (see the shaded zone).\n" ] } ], "source": [ "bs = bootstrap(df, sol, spends, Z, y, B=120)\n", "print(\"Return on the next euro (90% range):\")\n", "for c in CHANNELS:\n", " lo, hi = ci(bs[\"marginal\"][c]); print(f\" {c:>7}: EUR {np.median(bs['marginal'][c]):.2f} [{lo:.2f}, {hi:.2f}]\")\n", "glo, ghi = ci(bs[\"gain\"])\n", "print(f\"\\nRevenue gain from rebalancing: +{_eur(gain)}/wk 90% range [{_eur(glo)}, {_eur(ghi)}]/wk\")\n", "print(f\"Recommended Video spend {_eur(new['Video'])}/wk sits near the top of its observed range ({_eur(omax['Video'])}/wk);\")\n", "print(\"push much past that and the estimate turns into extrapolation (see the shaded zone).\")\n", "\n", "fig, ax = plt.subplots(1, 2, figsize=(11, 4.3))\n", "xpos = np.arange(len(CHANNELS))\n", "meds = [np.median(bs[\"marginal\"][c]) for c in CHANNELS]\n", "los = [np.median(bs[\"marginal\"][c]) - ci(bs[\"marginal\"][c])[0] for c in CHANNELS]\n", "his = [ci(bs[\"marginal\"][c])[1] - np.median(bs[\"marginal\"][c]) for c in CHANNELS]\n", "ax[0].bar(xpos, meds, .5, color=[SERIES[c] for c in CHANNELS], yerr=[los, his], capsize=5, ecolor=INK)\n", "ax[0].axhline(1.0, color=MUTED, ls=\":\", lw=1); ax[0].set_xticks(xpos); ax[0].set_xticklabels(CHANNELS)\n", "ax[0].set_title(\"Return on the next euro (90% range)\"); ax[0].set_ylabel(\"EUR revenue per extra EUR\")\n", "p = rec[\"Video\"]; xg = np.linspace(0, 2*omax[\"Video\"], 300)\n", "lo_band = steady_response(xg, p[\"beta\"]*0.8, p[\"ec\"]*1.2, p[\"slope\"])\n", "hi_band = steady_response(xg, p[\"beta\"]*1.2, p[\"ec\"]*0.85, p[\"slope\"])\n", "ax[1].fill_between(xg, lo_band, hi_band, color=SERIES[\"Video\"], alpha=.18, label=\"uncertainty\")\n", "ax[1].plot(xg, steady_response(xg, p[\"beta\"], p[\"ec\"], p[\"slope\"]), color=SERIES[\"Video\"], lw=2.2, label=\"Video curve\")\n", "ax[1].axvspan(omax[\"Video\"], 2*omax[\"Video\"], color=\"#e9b44c\", alpha=.16, label=\"beyond observed spend\")\n", "ax[1].axvline(omax[\"Video\"], color=\"#b8860b\", ls=\"--\", lw=1)\n", "ax[1].set_title(\"Past your observed spend, the curve is a guess\"); ax[1].set_xlabel(\"weekly Video spend\")\n", "ax[1].set_ylabel(\"weekly revenue\"); money_x(ax[1]); money_y(ax[1]); ax[1].legend(frameon=False)\n", "plt.tight_layout(); plt.savefig(FIG/\"fig_uncertainty.png\", bbox_inches=\"tight\"); plt.show()" ] }, { "cell_type": "markdown", "id": "7878450d", "metadata": {}, "source": [ "## 6. When channels move together, the model can't tell them apart\n", "\n", "The fit above worked because the three channels varied **independently** week to week. If you always\n", "scale Social and Search up and down together, their spend series become nearly identical, and no model\n", "can say which one drove the sales — their *combined* effect is solid, the *split* between them is not.\n", "This is **collinearity**, and it's the single most common way an MMM quietly goes wrong." ] }, { "cell_type": "code", "execution_count": 7, "id": "fbdb9d44", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:58.620295Z", "iopub.status.busy": "2026-07-23T13:15:58.620295Z", "iopub.status.idle": "2026-07-23T13:15:58.640283Z", "shell.execute_reply": "2026-07-23T13:15:58.640283Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Search-Social correlation when run independently: -0.033\n", "Search-Social correlation when scaled in lockstep: 0.952\n", "Above about 0.8 the individual split is no longer trustworthy; combine them and report the pair together.\n" ] } ], "source": [ "ind = collinearity(df)\n", "lock_df, _ = simulate_panel(collinear=True)\n", "lock = collinearity(lock_df)\n", "def _pair(rows, a, b):\n", " for x in rows:\n", " if {x[0], x[1]} == {a, b}: return x[2]\n", "print(f\"Search-Social correlation when run independently: {_pair(ind,'Search','Social')}\")\n", "print(f\"Search-Social correlation when scaled in lockstep: {_pair(lock,'Search','Social')}\")\n", "print(\"Above about 0.8 the individual split is no longer trustworthy; combine them and report the pair together.\")" ] }, { "cell_type": "markdown", "id": "c72ef054", "metadata": {}, "source": [ "## 7. Emit the article's numbers and the widget's check values\n", "All figures above and the JSON below are regenerated from this notebook — no number is hand-typed." ] }, { "cell_type": "code", "execution_count": 8, "id": "db31db01", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:58.645292Z", "iopub.status.busy": "2026-07-23T13:15:58.645292Z", "iopub.status.idle": "2026-07-23T13:15:58.669154Z", "shell.execute_reply": "2026-07-23T13:15:58.668446Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "wrote data/worked_example.json\n", "wrote data/checks.json with 4 pinned sets\n" ] } ], "source": [ "# ---- worked_example.json : the article's authority numbers ----\n", "def pos_pct(c): # where current spend sits on the fitted curve (Hill value)\n", " return float(hill(cur[c], rec[c][\"ec\"], rec[c][\"slope\"]))\n", "we = {\n", " \"label\": \"synthetic\",\n", " \"weeks\": int(len(df)),\n", " \"mean_weekly_revenue\": round(float(df.revenue.mean())),\n", " \"media_share\": round(float(sum(contrib[c].mean() for c in CHANNELS)/df.revenue.mean()), 3),\n", " \"wmape\": round(wmape, 4),\n", " \"sufficiency\": sufficiency(df),\n", " \"total_budget\": round(sum(cur.values())),\n", " \"channels\": {c: {\n", " \"current_spend\": round(cur[c]),\n", " \"observed_max\": round(omax[c]),\n", " \"marginal_return\": round(mar[c], 2),\n", " \"marginal_ci\": ci(bs[\"marginal\"][c]),\n", " \"curve_position\": round(pos_pct(c), 3),\n", " \"fitted_beta\": round(float(rec[c][\"beta\"])),\n", " \"fitted_half_saturation\": round(float(rec[c][\"ec\"])),\n", " \"fitted_slope\": round(float(rec[c][\"slope\"]), 3),\n", " \"recommended_spend\": round(new[c]),\n", " } for c in CHANNELS},\n", " \"reallocation\": {\n", " \"from\": \"Social\", \"to\": \"Video\",\n", " \"move_from_social\": round(cur[\"Social\"] - new[\"Social\"]),\n", " \"add_to_video\": round(new[\"Video\"] - cur[\"Video\"]),\n", " \"gain_per_week\": round(gain),\n", " \"gain_pct\": round(100*gain/old_rev, 1),\n", " \"gain_ci_week\": ci(bs[\"gain\"]),\n", " \"gain_per_year\": round(gain*52),\n", " },\n", " \"collinearity\": {\n", " \"independent_search_social\": _pair(ind, \"Search\", \"Social\"),\n", " \"lockstep_search_social\": _pair(lock, \"Search\", \"Social\"),\n", " },\n", "}\n", "(DATA/\"worked_example.json\").write_text(json.dumps(we, indent=2))\n", "print(\"wrote data/worked_example.json\")\n", "\n", "# ---- checks.json : forward-engine pins the widget must reproduce exactly ----\n", "def pin(name, beta, ec, slope, alpha, current, series):\n", " grid = [0, 2000, 5000, 10000, 20000]\n", " return {\n", " \"name\": name, \"beta\": beta, \"ec\": ec, \"slope\": slope, \"alpha\": alpha,\n", " \"current_spend\": current, \"grid\": grid, \"series\": series,\n", " \"expected\": {\n", " \"hill_at_ec\": round(float(hill(ec, ec, slope)), 6),\n", " \"steady_at_grid\": [round(float(steady_response(g, beta, ec, slope)), 4) for g in grid],\n", " \"marginal_at_current\": round(float(marginal_return(current, beta, ec, slope)), 6),\n", " \"adstock\": [round(float(v), 6) for v in geometric_adstock(np.array(series, float), alpha)],\n", " },\n", " }\n", "checks = {\"atol\": 1e-6, \"rtol\": 1e-6, \"max_lag\": MAX_LAG, \"sets\": [\n", " pin(\"saturating_concave\", 45000, 7000, 1.0, 0.4, 3500, [0,0,10000,0,0,4000,4000,0]),\n", " pin(\"s_shaped\", 52000, 4500, 1.3, 0.6, 14000, [8000,12000,16000,14000,10000,14000,18000,12000]),\n", " pin(\"efficient_middle\", 40000, 6500, 1.0, 0.5, 9000, [9000,9000,9000,9000,9000,9000,9000,9000]),\n", " pin(\"steep_slope\", 30000, 9000, 1.8, 0.3, 6000, [6000,0,0,12000,6000,0,0,6000]),\n", "]}\n", "(DATA/\"checks.json\").write_text(json.dumps(checks, indent=2))\n", "print(\"wrote data/checks.json with\", len(checks[\"sets\"]), \"pinned sets\")" ] }, { "cell_type": "markdown", "id": "efbd198f", "metadata": {}, "source": [ "---\n", "# Part 2 — run it on your own data\n", "\n", "Everything below runs the **same** model on a weekly CSV you supply. It stays on your machine.\n", "\n", "**Your file needs one row per week**, with these columns:\n", "\n", "| column | meaning | example |\n", "|---|---|---|\n", "| `week` | Monday of the week (any date format) | `2024-01-01` |\n", "| `revenue` | total revenue that week | `184320` |\n", "| `spend_` | one column per channel, weekly spend | `spend_Search` = `9100` |\n", "| *(optional)* `control_` | anything else that moves sales (e.g. `control_promo`) | `1` |\n", "\n", "Point `MY_CSV` at your file. Leave it as `None` to run on a fresh synthetic example." ] }, { "cell_type": "code", "execution_count": 9, "id": "91dc0fea", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:58.674242Z", "iopub.status.busy": "2026-07-23T13:15:58.673242Z", "iopub.status.idle": "2026-07-23T13:15:59.621167Z", "shell.execute_reply": "2026-07-23T13:15:59.620159Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "=== healthy example (fresh synthetic) ===\n", "Data looks healthy: enough weeks, enough independent variation.\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "\n", "Return on the next euro (spend now):\n", " Search: EUR 1.33 (EUR 8,622/wk)\n", " Social: EUR 0.93 (EUR 13,014/wk)\n", " Video: EUR 3.07 (EUR 3,607/wk)\n", "\n", "Rebalancing (same budget) would move about EUR 3,334/wk (+4.4%).\n" ] } ], "source": [ "MY_CSV = None # <-- set to \"your_weekly_data.csv\"\n", "\n", "def load_your_data(source):\n", " if isinstance(source, str):\n", " d = pd.read_csv(source)\n", " else:\n", " d = source.copy()\n", " d.columns = [c.strip() for c in d.columns]\n", " if \"week\" not in d.columns or \"revenue\" not in d.columns:\n", " raise ValueError(\"Your file needs a 'week' column and a 'revenue' column. \"\n", " \"See the table above for the exact column names.\")\n", " ch_cols = [c for c in d.columns if c.startswith(\"spend_\")]\n", " if len(ch_cols) < 2:\n", " raise ValueError(\"Need at least two 'spend_' columns to compare channels. \"\n", " f\"Found: {ch_cols or 'none'}.\")\n", " d = d.dropna(subset=[\"revenue\"] + ch_cols).reset_index(drop=True)\n", " channels = [c[len(\"spend_\"):] for c in ch_cols]\n", " warns = []\n", " weeks = len(d)\n", " if weeks < 52:\n", " raise ValueError(f\"Only {weeks} weeks of data. MMM needs at least 52 weeks (a full year of \"\n", " \"seasonality) before any answer can be trusted. Collect more, then re-run.\")\n", " controls = [c for c in d.columns if c.startswith(\"control_\")]\n", " dpe = weeks / (len(channels) + len(controls) + 4)\n", " if dpe < 10:\n", " warns.append(f\"Thin data: about {dpe:.1f} data points per effect (we like >= 10). The model \"\n", " \"will still run, but treat every number as ESTIMATED and read the ranges, not the \"\n", " \"point estimates.\")\n", " S = d[ch_cols].to_numpy(); C = np.corrcoef(S.T)\n", " for i in range(len(channels)):\n", " for j in range(i+1, len(channels)):\n", " if abs(C[i, j]) >= 0.8:\n", " warns.append(f\"'{channels[i]}' and '{channels[j]}' spend moved almost together \"\n", " f\"(correlation {C[i,j]:.2f}). The model cannot cleanly separate them; \"\n", " \"their combined effect is fine but the split between them is not. \"\n", " \"Consider merging them into one channel.\")\n", " tot = float(sum(d[c].sum() for c in ch_cols))\n", " for c, name in zip(ch_cols, channels):\n", " share = float(d[c].sum())/tot\n", " if share < 0.05:\n", " warns.append(f\"'{name}' is only {share:.0%} of total spend. Tiny channels are hard to \"\n", " \"measure; consider folding it into an 'Other' bucket.\")\n", " return d, channels, warns\n", "\n", "def run_your_data(source):\n", " d, channels, warns = load_your_data(source)\n", " for w in warns: print(\"! \" + w)\n", " if not warns: print(\"Data looks healthy: enough weeks, enough independent variation.\")\n", " global CHANNELS\n", " CHANNELS = channels # engine reads this list\n", " rec, sol, spends, Z, y = fit(d)\n", " cur = {c: float(d[f\"spend_{c}\"].mean()) for c in channels}\n", " omax = {c: float(d[f\"spend_{c}\"].max()) for c in channels}\n", " print(\"\\nReturn on the next euro (spend now):\")\n", " for c in channels:\n", " m = marginal_return(cur[c], rec[c][\"beta\"], rec[c][\"ec\"], rec[c][\"slope\"])\n", " print(f\" {c:>10}: EUR {m:.2f} ({_eur(cur[c])}/wk)\")\n", " new, gain, o, n = optimal_allocation(rec, cur, omax)\n", " print(f\"\\nRebalancing (same budget) would move about {_eur(gain)}/wk ({gain/o:+.1%}).\")\n", " return rec\n", "\n", "# demo on a fresh synthetic advertiser (proves the path reproduces a sane answer)\n", "demo_df, _ = simulate_panel(seed=101)\n", "demo_df = demo_df.rename(columns={\"week\": \"week\"}) # already named\n", "print(\"=== healthy example (fresh synthetic) ===\")\n", "_ = run_your_data(demo_df if MY_CSV is None else MY_CSV)\n", "CHANNELS = [\"Search\", \"Social\", \"Video\"] # restore" ] }, { "cell_type": "markdown", "id": "5c92dd9e", "metadata": {}, "source": [ "### The consolidator: raw daily rows → the weekly panel\n", "\n", "Most exports come as **daily** rows: one line per day per channel. This turns them into the weekly\n", "panel the model needs. Point it at a long-format spend table and a daily revenue series." ] }, { "cell_type": "code", "execution_count": 10, "id": "d51f2bf5", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:59.625168Z", "iopub.status.busy": "2026-07-23T13:15:59.625168Z", "iopub.status.idle": "2026-07-23T13:15:59.685575Z", "shell.execute_reply": "2026-07-23T13:15:59.685069Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "consolidated 560 days into 81 weeks:\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
weekspend_Searchspend_Socialrevenue
02024-01-011455.0905.024692.0
12024-01-0810111.07709.0182420.0
22024-01-159740.09126.0185712.0
\n", "
" ], "text/plain": [ " week spend_Search spend_Social revenue\n", "0 2024-01-01 1455.0 905.0 24692.0\n", "1 2024-01-08 10111.0 7709.0 182420.0\n", "2 2024-01-15 9740.0 9126.0 185712.0" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "def daily_to_weekly(spend_long, revenue_daily):\n", " \"\"\"spend_long: DataFrame [date, channel, spend]; revenue_daily: DataFrame [date, revenue].\"\"\"\n", " s = spend_long.copy(); s[\"date\"] = pd.to_datetime(s[\"date\"])\n", " s = s.groupby([pd.Grouper(key=\"date\", freq=\"W-MON\"), \"channel\"])[\"spend\"].sum().unstack(\"channel\").fillna(0.0)\n", " s.columns = [f\"spend_{c}\" for c in s.columns]\n", " r = revenue_daily.copy(); r[\"date\"] = pd.to_datetime(r[\"date\"])\n", " r = r.groupby(pd.Grouper(key=\"date\", freq=\"W-MON\"))[\"revenue\"].sum()\n", " out = s.join(r, how=\"inner\").reset_index().rename(columns={\"date\": \"week\"})\n", " return out\n", "\n", "# tiny demonstration\n", "rng = np.random.default_rng(0)\n", "days = pd.date_range(\"2024-01-01\", periods=560, freq=\"D\")\n", "sl = pd.DataFrame({\"date\": np.repeat(days, 2),\n", " \"channel\": [\"Search\", \"Social\"]*len(days),\n", " \"spend\": rng.uniform(500, 2000, 2*len(days))})\n", "rd = pd.DataFrame({\"date\": days, \"revenue\": rng.uniform(20000, 30000, len(days))})\n", "wk = daily_to_weekly(sl, rd)\n", "print(\"consolidated\", len(days), \"days into\", len(wk), \"weeks:\")\n", "display(wk.head(3).round(0))" ] }, { "cell_type": "markdown", "id": "f819286a", "metadata": {}, "source": [ "### The model fails *honestly* on bad data\n", "\n", "Good tools refuse to give a confident answer they can't support. Here the your-data path is handed\n", "two broken inputs — too few weeks, and two channels moved in lockstep — and it says so in plain words\n", "instead of printing a wrong number." ] }, { "cell_type": "code", "execution_count": 11, "id": "83659863", "metadata": { "execution": { "iopub.execute_input": "2026-07-23T13:15:59.689713Z", "iopub.status.busy": "2026-07-23T13:15:59.689713Z", "iopub.status.idle": "2026-07-23T13:16:00.687038Z", "shell.execute_reply": "2026-07-23T13:16:00.686030Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "=== too few weeks ===\n", "Refused: Only 30 weeks of data. MMM needs at least 52 weeks (a full year of seasonality) before any answer can be trusted. Collect more, then re-run.\n", "\n", "=== collinear channels (Search + Social in lockstep) ===\n", "! 'Search' and 'Social' spend moved almost together (correlation 0.95). The model cannot cleanly separate them; their combined effect is fine but the split between them is not. Consider merging them into one channel.\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "\n", "Return on the next euro (spend now):\n", " Search: EUR 0.96 (EUR 7,917/wk)\n", " Social: EUR 1.27 (EUR 12,947/wk)\n", " Video: EUR 2.66 (EUR 3,227/wk)\n", "\n", "Rebalancing (same budget) would move about EUR 1,968/wk (+2.8%).\n" ] } ], "source": [ "print(\"=== too few weeks ===\")\n", "short = simulate_panel(seed=5)[0].iloc[:30]\n", "try:\n", " run_your_data(short)\n", "except ValueError as e:\n", " print(\"Refused:\", e)\n", "\n", "print(\"\\n=== collinear channels (Search + Social in lockstep) ===\")\n", "coll = simulate_panel(seed=7, collinear=True)[0]\n", "_ = run_your_data(coll)\n", "CHANNELS = [\"Search\", \"Social\", \"Video\"]" ] }, { "cell_type": "markdown", "id": "29eb1d9a", "metadata": {}, "source": [ "## What to take away\n", "\n", "- MMM learns a **diminishing-returns curve** per channel and reads the **return on the next euro**.\n", "- The move is simple: shift budget from low-next-euro channels to high-next-euro ones.\n", "- Every number is a **range**, and the range widens fast past spend you've never run.\n", "- It's **correlational**: pair it with an experiment (a geo holdout) before betting big on any one number.\n", "\n", "See the article for the plain-language version, the caveats to watch for, and the interactive widget." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.12.8" } }, "nbformat": 4, "nbformat_minor": 5 }