Skip to main content
Back to blog
Marketing Mix Modeling
Hierarchical Models
Priors
Bayesian Statistics
PyMC-Marketing

When Your MMM Won't Mix: Sweep the Priors Before You Refit Everything

A SKU × customer × brand PyMC-Marketing case study. The first fit technically ran: it also produced chains that disagreed about price elasticity, hierarchical scales that barely moved, and elasticities with the wrong economic sign. Instead of guessing which prior to tighten and paying for each guess in NUTS hours, this workflow repairs the sampling geometry, screens 18 candidate priors with importance reweighting, and spends the one expensive refit where it counts.

Niall OultonAugust 20, 202613 min read

The model technically ran. That was the problem.

We had three brands, four variants per brand, five retail customers, four media channels, three commercial controls and 156 weeks of data. The model had enough structure to answer the questions the business actually cared about: which media channels were incremental, how promotion interacted with demand, and whether a 500 ml zero-sugar SKU at an online retailer was more price sensitive than the same brand in a national grocer.

It also had enough structure to become a sampling nightmare.

The first fit produced chains that disagreed about price elasticity, hierarchical scale parameters that barely moved, and a handful of SKU/customer combinations with the wrong economic sign. Increasing the number of draws would have made a larger file, not a better model.

This article shows a different workflow using pymc-marketing-sweep: rescue the geometry without changing the prior distribution, use PSIS importance reweighting to screen candidate priors quickly, reject candidates that move outside posterior support, then spend the expensive MCMC refit on the smallest safe change.

Reproducibility note: the business data in this article are fully synthetic. The complete case study (2.7 MB zip) includes deterministic posterior fixtures so every chart renders immediately without a multi-hour MCMC run, plus the full PyMC-Marketing scripts that recreate the workflow using real saved .nc models. The quoted diagnostics below are from the deterministic fixture until you rerun run_demo.sh in a PyMC-Marketing environment.

Weekly units and commercial drivers for example SKU/customer cells in the synthetic soft-drinks panel
The synthetic panel: 60 SKU/customer cells, 156 weeks, correlated price, promotion and media.

The model we actually wanted

The synthetic business is a soft-drinks portfolio with:

  • 3 brands: Aster, Brio and Cinder
  • 4 variants per brand: Original 330 ml, Zero 330 ml, Original 500 ml and Zero 500 ml
  • 5 customers: NorthMart, ValueHub, CityExpress, FreshBasket and eGrocer
  • 156 weekly observations from 2023 through 2025
  • 60 SKU/customer cells and 9,360 panel rows
  • 4 media channels: TV, paid social, paid search and retailer media
  • 3 commercial controls: price discount signal, promotion depth and weighted distribution

The data generator deliberately makes the world uncomfortable in realistic ways. List prices inflate over time. Promotions vary by customer. Zero variants ramp into distribution. Retail media rises during promotion weeks. Paid search reacts to campaigns and promotional pressure. National media pulses across brands rather than arriving as independent white noise.

The target is log_units. Price enters as:

price_discount_signal = -log(actual_price / list_price)

That small choice is useful. If the coefficient on price_discount_signal is 1.7, the conventional log-log price elasticity is approximately −1.7. We can therefore inspect a SKU/customer price elasticity directly from the model’s control coefficient while still using the multidimensional PyMC-Marketing MMM.

The true elasticities in the synthetic data range from roughly −1.24 to −2.48. That is wide enough to make partial pooling useful, but not so wide that every retailer should behave like a different universe.

Why price becomes the canary in the coal mine

Price is rarely varied independently in retail data.

In this simulation, discount depth and the engineered price signal correlate at roughly 0.99 by construction. Retail media also rises around promotional events. Distribution shifts at the same time that new packs launch. The model has to decide whether demand moved because a product became cheaper, because the retailer supported it, because distribution expanded, because national advertising had carry-over, or because several of those things happened together.

Correlation matrix of commercial drivers in the synthetic panel, with discount depth and the price signal correlated at roughly 0.99
Confounding by construction: promotions, retail media and distribution move together.

This is exactly where a hierarchical Bayesian MMM earns its keep. We want the model to borrow information across customers and SKUs, but not force everything to be identical.

The baseline control prior was intentionally generous:

python
from pymc_extras.prior import Prior

gamma_control = Prior(
    "Normal",
    mu=0.25,
    sigma=Prior(
        "HalfNormal",
        sigma=1.0,
        dims=("control", "brand", "variant"),
    ),
    dims=("control", "brand", "customer", "variant"),
    centered=True,
)

Media response was also allowed to vary heavily across channel, brand, customer and variant:

python
from pymc_marketing.special_priors import LogNormalPrior

saturation_beta = LogNormalPrior(
    mean=Prior(
        "Gamma",
        mu=0.22,
        sigma=0.18,
        dims=("channel", "brand"),
    ),
    std=Prior(
        "HalfNormal",
        sigma=0.80,
        dims=("channel", "brand"),
    ),
    dims=("channel", "brand", "customer", "variant"),
    centered=True,
)

Nothing here is absurd in isolation. The problem is the combination: many weakly identified local effects, broad hierarchical scales, correlated commercial drivers and a centred hierarchical parameterisation.

The first fit is not merely noisy. It is untrustworthy

The deterministic failure fixture has 187 divergent transitions, a maximum rank R-hat of about 2.50, and a minimum bulk ESS below 5. Two SKU/customer price elasticities even flip to the wrong sign on posterior mean.

The trace for one Aster Zero 500 ml elasticity tells the story more clearly than any summary table:

Trace plot of one Aster Zero 500 ml price elasticity where the four chains spend long periods in different regions
Four chains, several stories: the centred fit’s trace for one price elasticity.

Those chains are not four noisy estimates of the same posterior. They are spending long periods in different regions.

This distinction matters because a prior sweep based on importance reweighting should not be run blindly on a posterior that has not explored its target distribution. Reweighting bad samples more cleverly does not turn them into good samples.

That is why the workflow has an explicit rescue step.

Step 1: repair the parameterisation, not the prior

For the reference fit we change the centred hierarchical priors to non-centred form while keeping the intended prior distributions the same.

python
# Same Normal hierarchy, different parameterisation
gamma_control = Prior(
    "Normal",
    mu=0.25,
    sigma=Prior(
        "HalfNormal",
        sigma=1.0,
        dims=("control", "brand", "variant"),
    ),
    dims=("control", "brand", "customer", "variant"),
    centered=False,
)

PyMC-Marketing’s own multidimensional documentation calls out non-centred parameterisations as a practical way to address computational difficulties in hierarchical MMMs. It also recommends building complexity gradually and checking model performance as dimensions and partial pooling are added.

The non-centred reference fit is not the final model. It is a trustworthy representation of the same broad-prior target that we can use for the sweep.

In the fixture, maximum R-hat drops to about 1.006 and minimum bulk ESS rises above 1,360.

Bar chart of maximum rank R-hat falling from about 2.5 in the centred fit to about 1.006 in the non-centred reference
Maximum R-hat, centred versus non-centred.
Bar chart of minimum bulk ESS rising from below 5 in the centred fit to above 1,360 in the non-centred reference
Minimum bulk ESS, centred versus non-centred.

This is the point at which sweeping becomes statistically meaningful.

Step 2: stop guessing which prior should be tighter

A common workflow after a bad MMM fit is surprisingly manual:

  • Pick a prior that looks suspicious.
  • Tighten it.
  • Refit for an hour or six.
  • Discover that either nothing improved or the posterior changed too much.
  • Repeat.

With several hierarchical scales, this becomes a combinatorial search paid for in NUTS hours.

pymc-marketing-sweep takes a saved PyMC-Marketing model and asks a cheaper question first: if I changed these priors, does the existing well-sampled posterior contain enough support to approximate the new posterior safely?

For a candidate prior p₁(θ) and the fitted prior p₀(θ), each posterior draw receives an importance ratio proportional to:

p₁(θ) / p₀(θ)

The package then uses Pareto-smoothed importance sampling to stabilise the weights and reports:

  • Pareto-k as an overlap diagnostic
  • importance reweighting ESS
  • ESS as a fraction of available posterior draws
  • posterior distortion on variables you care about
  • a distance from the fitted prior
  • the least invasive candidate that satisfies all configured guardrails

It does not call reweighting ESS “predicted NUTS ESS”. Those are different quantities. The only proof that the sampler is fixed is the final validation refit.

The sweep we ran

We targeted three places where the original model had too much room to move:

  • The hierarchical scale controlling customer variation in commercial controls, including price.
  • The cross-customer/variant scale on media saturation beta.
  • The cell-level intercept prior.

The YAML is intentionally readable by a human or an agent:

yaml
sweep:
  mode: product
  max_candidates: 100

  priors:
    price_control_pooling:
      prior_path: gamma_control.sigma
      posterior_variable: gamma_control_sigma
      parameters:
        sigma:
          values: [0.50, 0.35, 0.25]

    media_response_pooling:
      prior_path: saturation_beta.std
      posterior_variable: saturation_beta_std
      parameters:
        sigma:
          values: [0.50, 0.35, 0.25]

    baseline_regularisation:
      prior_path: intercept
      posterior_variable: intercept_contribution
      parameters:
        sigma:
          values: [1.5, 1.0]

constraints:
  max_pareto_k: auto
  min_reweight_ess: 250
  min_reweight_ess_fraction: 0.08
  max_relative_rmse: 0.05

That is 18 candidate prior configurations from one reference fit.

The command is simply:

bash
python scripts/run_prior_sweep.py

or, when the package is installed:

bash
pmm-sweep run configs/sku_customer_brand_sweep.yaml

The best part of the sweep is what it refuses to recommend

Aggressive regularisation is seductive. A very tight prior can make a posterior look stable because it removes most of the space the sampler was struggling to explore.

That does not mean it is supported by the fitted likelihood.

Several candidates in this example push Pareto-k above 0.8 and collapse importance ESS to around 40 effective draws. The sweep rejects them.

Scatter of sweep candidates by Pareto-k and reweighting ESS, with aggressive candidates failing the guardrails
The sweep frontier: aggressive candidates fail the Pareto-k and ESS guardrails.

The recommended fixture candidate is deliberately boring:

python
recommended_model_config_patch = {
    "gamma_control": Prior(
        "Normal",
        mu=0.25,
        sigma=Prior(
            "HalfNormal",
            sigma=0.50,
            dims=("control", "brand", "variant"),
        ),
        dims=("control", "brand", "customer", "variant"),
        centered=False,
    ),

    "intercept": Prior(
        "Normal",
        mu=7.25,
        sigma=1.50,
        dims=("brand", "customer", "variant"),
    ),

    "saturation_beta": LogNormalPrior(
        mean=Prior(
            "Gamma",
            mu=0.22,
            sigma=0.18,
            dims=("channel", "brand"),
        ),
        std=Prior(
            "HalfNormal",
            sigma=0.50,
            dims=("channel", "brand"),
        ),
        dims=("channel", "brand", "customer", "variant"),
        centered=False,
    ),
}

Relative to the reference posterior, that candidate has:

DiagnosticFixture result
Pareto-k0.035
Reweighting ESS3,859
Reweighting ESS fraction68.9%
Maximum protected-variable relative RMSE0.23%
Candidate statusRecommended

The point is not that 0.50 is a magical value. It is that this configuration is a supported change. We can make it without asking the data to teleport into a posterior region the reference fit never visited.

What happens to price elasticity before the final refit?

The bad centred fit has price-elasticity RMSE of about 0.73 against the known synthetic truth. The non-centred reference reduces that to about 0.18 simply because it actually samples the intended posterior.

The sweep approximation is slightly better again at roughly 0.17.

That may sound modest, but it is exactly what you want from a screening step. Importance reweighting should not manufacture a completely different model. It should tell you which prior change is plausible enough to deserve the real refit.

Estimated versus true price elasticities from the centred fit, with large errors and wrong-sign cells
Price elasticities from the centred fit: RMSE around 0.73, including wrong-sign cells.

Step 3: run one real validation fit

We now rebuild the PyMC-Marketing model with the generated recommended_model_config.py patch and fit it once.

bash
python scripts/fit_fixed_model.py

That is the expensive step, but we pay for it once rather than across the entire grid.

The validation fixture has:

  • 0 divergences
  • maximum R-hat around 1.004
  • minimum bulk ESS above 2,200
  • price-elasticity RMSE around 0.036
  • 0 wrong-sign price elasticities
Estimated versus true price elasticities after the validation refit, tightly aligned along the diagonal
Price elasticities after the validation refit: RMSE around 0.036, no wrong signs.

The effect is easier to see when we follow one brand across every customer and variant.

Aster brand price elasticities across every customer and variant, before and after the repair workflow
One brand, every customer and variant: differences survive, unsupported wandering does not.

The final model still allows differences. eGrocer can be more price sensitive than NorthMart. Zero 500 ml can be more elastic than Original 330 ml. What disappears is unsupported retailer/SKU wandering caused by an over-flexible hierarchy and poor geometry.

Price-elasticity RMSE falling from 0.73 in the bad fit to 0.18 in the reference, 0.17 in the sweep approximation and 0.036 after validation
Elasticity RMSE along the workflow: bad fit, reference, sweep approximation, validation refit.

This is not “tight priors fix convergence”

That would be the wrong lesson.

The actual workflow is:

bad mixing
    ↓
identify whether geometry or model structure is broken
    ↓
get a trustworthy reference posterior
    ↓
sweep plausible prior changes
    ↓
reject changes with poor posterior overlap
    ↓
protect business quantities from excessive movement
    ↓
choose the smallest supported change
    ↓
run one real validation fit

There are problems a prior sweep should not be expected to repair:

  • omitted-variable bias
  • incorrect causal structure
  • a completely misspecified likelihood
  • non-identifiability caused by perfectly collinear inputs
  • a reference posterior that never explored its own target
  • a proposed prior so different that PSIS says the old posterior has no useful overlap

The package treats the last case as a failure, not an invitation to extrapolate harder.

It also gives an agent a much safer job to do

Once this is in a repository, you can ask Codex, Claude or Cursor something like:

Read README.md and the sweep package instructions first.

The SKU/customer/brand MMM has unstable price elasticities.
Do not edit the statistical engine.

Inspect the saved reference model and confirm which posterior variable
corresponds to gamma_control.sigma.

Test HalfNormal scale values 0.50, 0.35 and 0.25 for the control hierarchy.
Also test saturation_beta.std at 0.50, 0.35 and 0.25.
Test intercept sigma at 1.5 and 1.0.

Keep max Pareto-k on auto.
Require at least 8% reweighting ESS.
Do not allow gamma_control or channel contributions to move by more than
5% relative RMSE.

Run the sweep.
Explain which candidates were rejected and why.
Return the exact PyMC-Marketing prior patch for the recommended candidate.
Do not claim the sampler is fixed until you run the validation refit and
check divergences, R-hat and ESS.

That is a much better use of an agent than “keep changing priors until the model looks okay”.

Approximate future predictions are possible too

The sweep package can importance-resample the joint posterior for a candidate prior and temporarily pass those draws back into the loaded PyMC-Marketing model. It then delegates future prediction to PyMC-Marketing’s own predict_posterior method.

The config in this repo includes:

yaml
prediction:
  data: ../data/future_media_price_plan.csv
  draws: 1000

So the same sweep can produce an approximate prediction under the candidate prior before the expensive validation refit.

That is useful for questions such as:

  • Does this regularisation materially change the next-quarter sales forecast?
  • Do contribution shares move enough to change a budget decision?
  • Does the price response become economically sensible without rewriting the media story?

Again, the prediction is an importance-reweighted approximation. Once a candidate is promoted, the production answer should come from the validation fit.

How to run it yourself

Everything in this article ships in the case study zip: the data generator, the model scripts, the figures, and the complete pymc-marketing-sweep package. Create an environment and install the sweep package:

bash
python -m venv .venv
source .venv/bin/activate
pip install -r requirements.txt
pip install -e ./pymc-marketing-sweep

Generate the synthetic data:

bash
python scripts/generate_data.py

Fit the deliberately difficult centred model, then the non-centred broad-prior reference:

bash
python scripts/fit_bad_model.py
python scripts/diagnose_fit.py --model outputs/runtime/bad_model.nc --label bad

python scripts/fit_reference_model.py
python scripts/diagnose_fit.py --model outputs/runtime/reference_model.nc --label reference

Run the sweep:

bash
python scripts/run_prior_sweep.py

Inspect:

outputs/runtime/sweep/
├── sweep_results.csv
├── recommended_scenario.yaml
├── recommended_model_config.py
├── recommended_importance_weights.npz
├── recommended_weighted_posterior_summary.nc
└── recommended_future_prediction_summary.csv

Then perform the only refit that matters:

bash
python scripts/fit_fixed_model.py
python scripts/diagnose_fit.py --model outputs/runtime/fixed_model.nc --label fixed
python scripts/extract_price_elasticities.py

Or run everything with ./run_demo.sh.

Why I like this workflow

A complex MMM has two competing requirements.

It needs enough flexibility to represent the business. A national brand can behave differently by retailer. A 500 ml pack can have a different price response from a 330 ml pack. Retail media can behave differently from TV. Pooling everything would be convenient and wrong.

But every extra dimension creates another place where weak identification can turn into posterior geometry.

The answer is not to avoid hierarchical MMMs. It is to make the iteration loop cheaper and more disciplined.

A prior sweep turns “try a few priors and see” into a constrained search:

  • Is the candidate supported by the existing posterior?
  • How much effective information survives the reweighting?
  • Which business quantities move?
  • Is a more aggressive candidate actually worth the distortion?
  • What exact PyMC-Marketing configuration should be validated next?

That is a much more useful conversation than debating whether sigma=0.5 simply “feels tighter”.

And when the final validation fit comes back clean, with stable media contributions and price elasticities that line up across SKU, customer and brand in economically sensible ways, you know why it improved.

Not because you asked NUTS to work harder.

Because you gave the model a better-shaped problem.

Technical notes

At the time this case study was assembled, the latest tagged PyMC-Marketing release was 0.19.4, while the live documentation was already documenting the modern multidimensional API and the migration path towards the 1.0 line. The bundled sweep package contains compatibility handling for the current multidimensional class and the newer root import path.

PyMC-Marketing’s multidimensional MMM supports additional panel dimensions, hierarchical priors through Prior(..., dims=...), original-scale contribution deterministics, and saved model workflows. The official multidimensional example also explicitly discusses centred versus non-centred hierarchical parameterisations and the computational trade-offs of partial pooling.

SIMBA builds on PyMC-Marketing for transparent Bayesian MMM, including the convergence and diagnostics workflow described here. If your hierarchical MMM will not mix and you want a second pair of eyes, book a call.

Published on August 20, 2026 by Niall Oulton

All posts