diff --git a/.github/workflows/mmm-interop.yml b/.github/workflows/mmm-interop.yml
new file mode 100644
index 00000000..aa20a30f
--- /dev/null
+++ b/.github/workflows/mmm-interop.yml
@@ -0,0 +1,129 @@
+name: MMM Interop
+
+# Schema-drift canary for the diff_diff.mmm exporters against the REAL MMM
+# frameworks (pymc-marketing, google-meridian). Lives in its own workflow file
+# because `paths:` is workflow-level: inside notebooks.yml these heavy jobs
+# would fire on every notebook-prose PR and pay full TensorFlow/PyMC installs.
+#
+# Version policy: CI installs FLOATING ranges (like the balance job's
+# "balance>=0.21") so the weekly cron genuinely resolves new upstream releases
+# and acts as an upstream-drift canary; the tutorial notebooks and their local
+# execution venvs stay pinned to exact versions for reproducibility.
+#
+# Fragile edge: google-meridian pins an exact tfp-nightly build. If that wheel
+# ever disappears from PyPI, first bump/resolve google-meridian to the newest
+# release; if none resolves, flip this workflow to cron-only in a hotfix PR —
+# never leave it red on labeled PRs.
+
+on:
+ push:
+ branches: [main]
+ paths:
+ # The exporters' behavior is owned by mmm.py plus the pinned
+ # AggregationResult container shape (aggregation.py); estimator-internal
+ # changes cannot alter the exporter contract without touching these, and
+ # the normal suite owns the container contract itself - so the heavy
+ # framework installs do not fan out to unrelated library changes.
+ - 'diff_diff/mmm.py'
+ - 'diff_diff/aggregation.py'
+ - 'pyproject.toml'
+ - 'tests/test_mmm_interop_pymc.py'
+ - 'tests/test_mmm_interop_meridian.py'
+ - '.github/workflows/mmm-interop.yml'
+ pull_request:
+ branches: [main]
+ types: [opened, synchronize, reopened, labeled, unlabeled]
+ paths:
+ # The exporters' behavior is owned by mmm.py plus the pinned
+ # AggregationResult container shape (aggregation.py); estimator-internal
+ # changes cannot alter the exporter contract without touching these, and
+ # the normal suite owns the container contract itself - so the heavy
+ # framework installs do not fan out to unrelated library changes.
+ - 'diff_diff/mmm.py'
+ - 'diff_diff/aggregation.py'
+ - 'pyproject.toml'
+ - 'tests/test_mmm_interop_pymc.py'
+ - 'tests/test_mmm_interop_meridian.py'
+ - '.github/workflows/mmm-interop.yml'
+ schedule:
+ # Weekly Sunday 7am UTC (offset from notebooks.yml) — upstream-drift canary:
+ # the floating installs resolve the latest framework releases.
+ - cron: '0 7 * * 0'
+
+permissions:
+ contents: read
+
+jobs:
+ mmm-pymc-interop:
+ name: PyMC-Marketing interop smoke
+ # Same ready-for-ci label gate as the notebooks.yml jobs (keep in sync).
+ if: >-
+ github.event_name != 'pull_request'
+ || (contains(github.event.pull_request.labels.*.name, 'ready-for-ci')
+ && (github.event.action != 'labeled' && github.event.action != 'unlabeled'
+ || github.event.label.name == 'ready-for-ci'))
+ runs-on: ubuntu-latest
+
+ steps:
+ - uses: actions/checkout@3d3c42e5aac5ba805825da76410c181273ba90b1 # v7
+
+ - name: Set up Python
+ uses: actions/setup-python@5fda3b95a4ea91299a34e894583c3862153e4b97 # v7.0.0
+ with:
+ python-version: '3.12'
+
+ - name: Install dependencies
+ run: |
+ pip install numpy pandas scipy pytest "pymc-marketing>=1.0,<2"
+ # Add repo root to Python path so tests can import diff_diff
+ # (pip install -e . requires the Rust/maturin toolchain; .pth avoids that)
+ python -c "import site; print(site.getsitepackages()[0])" | xargs -I{} sh -c 'echo "$PWD" > {}/diff_diff_dev.pth'
+
+ - name: Import canary
+ # Its own step so an install/transitive-import failure FAILS the job
+ # instead of being swallowed by pytest.importorskip into an all-skipped
+ # green run.
+ run: python -c "import pymc_marketing; print(pymc_marketing.__version__)"
+
+ - name: Run interop smoke tests
+ env:
+ DIFF_DIFF_BACKEND: python
+ run: |
+ pytest tests/test_mmm_interop_pymc.py -v --tb=short
+
+ mmm-meridian-interop:
+ name: Meridian interop smoke
+ # Same ready-for-ci label gate as the notebooks.yml jobs (keep in sync).
+ if: >-
+ github.event_name != 'pull_request'
+ || (contains(github.event.pull_request.labels.*.name, 'ready-for-ci')
+ && (github.event.action != 'labeled' && github.event.action != 'unlabeled'
+ || github.event.label.name == 'ready-for-ci'))
+ runs-on: ubuntu-latest
+
+ steps:
+ - uses: actions/checkout@3d3c42e5aac5ba805825da76410c181273ba90b1 # v7
+
+ - name: Set up Python
+ uses: actions/setup-python@5fda3b95a4ea91299a34e894583c3862153e4b97 # v7.0.0
+ with:
+ python-version: '3.12'
+
+ - name: Install dependencies
+ run: |
+ pip install numpy pandas scipy pytest "google-meridian>=1.8,<2"
+ # Add repo root to Python path so tests can import diff_diff
+ # (pip install -e . requires the Rust/maturin toolchain; .pth avoids that)
+ python -c "import site; print(site.getsitepackages()[0])" | xargs -I{} sh -c 'echo "$PWD" > {}/diff_diff_dev.pth'
+
+ - name: Import canary
+ # Its own step so an install/transitive-import failure FAILS the job
+ # instead of being swallowed by pytest.importorskip into an all-skipped
+ # green run.
+ run: python -c "import meridian; print(meridian.__version__)"
+
+ - name: Run interop smoke tests
+ env:
+ DIFF_DIFF_BACKEND: python
+ run: |
+ pytest tests/test_mmm_interop_meridian.py -v --tb=short
diff --git a/.github/workflows/notebooks.yml b/.github/workflows/notebooks.yml
index 274805a7..d2e268db 100644
--- a/.github/workflows/notebooks.yml
+++ b/.github/workflows/notebooks.yml
@@ -10,6 +10,9 @@ on:
- '.github/workflows/notebooks.yml'
# the interop drift guard runs only in this workflow (balance installed)
- 'tests/test_t26_composition_drift_calibration_drift.py'
+ # the MMM tutorial drift guards run as a step in execute-notebooks
+ - 'tests/test_t29_mmm_calibration_pymc_drift.py'
+ - 'tests/test_t30_mmm_calibration_meridian_drift.py'
pull_request:
branches: [main]
types: [opened, synchronize, reopened, labeled, unlabeled]
@@ -20,6 +23,9 @@ on:
- '.github/workflows/notebooks.yml'
# the interop drift guard runs only in this workflow (balance installed)
- 'tests/test_t26_composition_drift_calibration_drift.py'
+ # the MMM tutorial drift guards run as a step in execute-notebooks
+ - 'tests/test_t29_mmm_calibration_pymc_drift.py'
+ - 'tests/test_t30_mmm_calibration_meridian_drift.py'
schedule:
# Weekly Sunday 6am UTC — smoke test that notebooks still execute cleanly
- cron: '0 6 * * 0'
@@ -63,6 +69,8 @@ jobs:
--ignore=docs/tutorials/06_power_analysis.ipynb \
--ignore=docs/tutorials/10_trop.ipynb \
--ignore=docs/tutorials/26_composition_drift_calibration.ipynb \
+ --ignore=docs/tutorials/29_mmm_calibration_pymc.ipynb \
+ --ignore=docs/tutorials/30_mmm_calibration_meridian.ipynb \
-v \
--tb=short
# Excluded notebooks (too slow for pure-Python CI without Rust backend):
@@ -72,6 +80,21 @@ jobs:
# 26_composition_drift_calibration — requires the balance package;
# runs in the isolated interop-notebooks job below so this job's
# minimal env keeps enforcing that tutorials add no dependencies
+ # 29/30_mmm_calibration_* — require pymc-marketing / google-meridian,
+ # which cannot even share ONE environment (arviz conflict) and run
+ # multi-minute MCMC; executed locally with committed outputs. Their
+ # drift guards run in the next step, and the mmm-interop.yml workflow
+ # smoke-tests the exporters against the real frameworks (no sampling).
+
+ - name: Run MMM tutorial drift guards
+ # The MMM notebooks are excluded from execution above, so their committed
+ # surface is protected here: DGP/sampler/integration needles + diff-diff
+ # recomputation of the DiD-side numbers (no MMM frameworks needed).
+ env:
+ DIFF_DIFF_BACKEND: python
+ run: |
+ pytest tests/test_t29_mmm_calibration_pymc_drift.py \
+ tests/test_t30_mmm_calibration_meridian_drift.py -v --tb=short
- name: Upload failed notebook outputs
if: failure()
diff --git a/.github/workflows/rust-test.yml b/.github/workflows/rust-test.yml
index 4883441b..5ca666a2 100644
--- a/.github/workflows/rust-test.yml
+++ b/.github/workflows/rust-test.yml
@@ -19,6 +19,7 @@ on:
# this suite (enforced by TestRustTestWorkflowPathFilter).
- '.github/workflows/docs-tests.yml'
- '.github/workflows/notebooks.yml'
+ - '.github/workflows/mmm-interop.yml'
- '.github/workflows/ci-gate.yml'
- '.github/workflows/release-build-check.yml'
# lint.yml carries the ruff/black pin-sync contract locked by
@@ -51,6 +52,7 @@ on:
# this suite (enforced by TestRustTestWorkflowPathFilter).
- '.github/workflows/docs-tests.yml'
- '.github/workflows/notebooks.yml'
+ - '.github/workflows/mmm-interop.yml'
- '.github/workflows/ci-gate.yml'
- '.github/workflows/release-build-check.yml'
# lint.yml carries the ruff/black pin-sync contract locked by
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 91bf6f6a..bad43cff 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -32,6 +32,29 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
REGISTRY's LWDiD section gains the maintainer implementation notes (RI
convention, dual overall-ATT conventions, Sec 4.3 rejection, inference
default, API conformance).
+- **Executed MMM calibration tutorials + CI interop canary.** Two new tutorials
+ run the diff-diff -> MMM hand-off end-to-end against the real frameworks:
+ `29_mmm_calibration_pymc.ipynb` (a staggered geo spend-boost experiment
+ exported with `to_pymc_marketing_lift_test` and fed to pymc-marketing 1.0's
+ `MMM.add_lift_test_measurements`) and `30_mmm_calibration_meridian.ipynb` (a
+ geo-holdout launch exported with `to_meridian_roi_prior` +
+ `meridian_calibration_mask`, with the generated `to_code()` snippet executed
+ verbatim in google-meridian 1.8). Both fit the MMM without and with the
+ calibration and demonstrate the ROI posterior moving to the simulated truth
+ and tightening, with in-notebook sampler-health and acceptance asserts. The
+ heavy frameworks stay out of diff-diff's dependencies: the notebooks execute
+ locally (committed outputs; they cannot even share one environment - the two
+ frameworks pin conflicting `arviz` ranges) and are guarded by
+ diff-diff-only drift tests (`tests/test_t29_*`/`test_t30_*`), while a new
+ `mmm-interop.yml` workflow smoke-tests the exporters against floating
+ framework releases on exporter-touching PRs (paths-scoped to `mmm.py` +
+ `aggregation.py`, the two files that own the exporter contract) and a
+ weekly cron
+ (`tests/test_mmm_interop_pymc.py`, `tests/test_mmm_interop_meridian.py` -
+ schema, dims, value-retention, and Meridian defaults-drift pins; the
+ 1.7.0-pinned `to_code()` templates and `_MERIDIAN_PARAM_DEFAULTS` are
+ execution-validated on Meridian 1.8.0). The practitioner guide's Step 8
+ gains the route-qualified MMM hand-off pointer.
- **`results.aggregate('total')` - the estimator-owned total incremental
outcome** on CallawaySantAnna, EfficientDiD, ImputationDiD, and TwoStageDiD,
promoted into the library-wide aggregation vocabulary. The single
diff --git a/TODO.md b/TODO.md
index 2ef1d085..9a0fa06c 100644
--- a/TODO.md
+++ b/TODO.md
@@ -76,11 +76,11 @@ generic sparse-FE, QR+SVD rank-detection redundancy, `check_finite` bypass — m
| Issue | Location | Origin | Effort | Priority |
|-------|----------|--------|--------|----------|
+| Optional scheduled end-to-end execution gate for the MMM tutorials (29/30): a cron-only workflow (or extension of `mmm-interop.yml`) that executes both notebooks in isolated exact-pin environments, so a stale/invalid committed posterior cannot stay green indefinitely - today the hybrid posture (deliberate: notebooks execute locally with committed outputs; CI smoke-tests the exporters without sampling; drift tests pin source + committed-output needles) leaves the MCMC claims un-re-executed in CI | `.github/workflows/mmm-interop.yml`, `docs/tutorials/29_mmm_calibration_pymc.ipynb`, `docs/tutorials/30_mmm_calibration_meridian.ipynb` | mmm-interop | Mid | Low |
| Committed `fixest::feols` event-study golden for TWFE `event_study=True` (within + pooled specs, unbalanced + covariate panels, matched CR1 cluster convention, per-period effects + vcov block) - the in-suite gates are shared-core cross-checks (TWFE-within == MPD-absorb, pooled == MPD bit-exact), so a defect common to the shared core would pass; the live-R harness (`benchmarks/R/benchmark_multiperiod.R`, `feols(y ~ treated * time_f \| unit)`) validated the within design in `docs/benchmarks.rst` but is not a committed regression test - follow the `fixest_did_twfe_golden.json` committed-golden pattern (pytest.skip when absent) | `tests/test_fixest_did_twfe_parity.py`, `benchmarks/R/` | 3(a) R2 | Mid | Medium |
| Type-blind `n_bootstrap` acceptance in already-validated estimators - HAD bool (`isinstance(..., int)` passes `True`, runs as 1 replicate), dCDH bool+float (its bare `< 0` check passes both `True` and `2.5`), TROP float (`2.5` passes the `>= 2` floor), SyntheticDiD float under all three variance methods + bool/negative under jackknife (its floor check is skipped there) - align these local checks with the `utils.validate_n_bootstrap` type guard (M-081 kept them out of the sweep: it scoped to previously-UNvalidated estimators only) | `diff_diff/had.py`, `diff_diff/chaisemartin_dhaultfoeuille.py`, `diff_diff/trop.py`, `diff_diff/synthetic_did.py` | 2(d) PR-B | Quick | Low |
| Evaluate adding the `BaseEstimator` param surface (get_params/set_params) to the exported classes that never had it - `PowerAnalysis`, `LinearRegression`, `BusinessReport`, `DiagnosticReport`, `TWFEWeightsResult` (a NEW public surface, deliberately out of the 2(c)-i pure-refactor scope; `LinearRegression` is the one `fit`-bearing class excluded from the contract suite's roster-completeness test). | `diff_diff/linalg.py`, `diff_diff/power.py` | mixin PR | Mid | Low |
| Tighten the mypy suppressions that back the enforced-zero posture: burn down `prep_dgp`'s per-module `[index]` override (needs a None-vs-array restructure that preserves the seeded RNG stream), and evaluate re-enabling the globally disabled codes (`arg-type`, `return-value`, `var-annotated`, `assignment`) one at a time — `assignment` alone hid several real annotation drifts found during the 2026-07 triage. | `pyproject.toml` `[tool.mypy]`, `diff_diff/prep_dgp.py` | lint-CI | Mid | Low |
-| MMM interop PR-B: calibration tutorial notebook (fit DiD/CS -> scope -> `to_pymc_marketing_lift_test` / `to_meridian_roi_prior`) + a `llms-practitioner.txt` Step 8 pointer to the exporters as the MMM hand-off. | `docs/tutorials/`, `diff_diff/guides/llms-practitioner.txt` | mmm-interop | Mid | Low |
| Tracking-file contract guard test: reject NEW active deferred-work pointers at `TODO.md` (deferred rows live in `DEFERRED.md`; allowlist for historical/past-tense prose and actionable-row pointers) and assert rows cross-linking a `docs/v4-deprecations.yaml` `M-xxx` id don't restate ledger status. Origin: tracking-split local review R2. | `tests/`, `TODO.md`, `DEFERRED.md` | tracking-split | Quick | Low |
| Real-data CI canary for dataset-backed replication tests: `test_methodology_lwdid.py`'s Prop 99 / Walmart goldens skip (visibly) when loaders fall back to synthetic; add a lane or canary asserting `df.attrs["source"] == "lwdid_ssc_ancillary"` in CI so network regressions cannot silently de-gate the replication tests. Follow-on from the loader-fallback repair (#723), which made provenance explicit but deliberately did not add a network-dependent CI lane. | `tests/test_methodology_lwdid.py`, `.github/workflows/` | LWDiD validation suite | Quick | Low |
| df-provenance completion follow-up (PR C descoping): a Wooldridge public per-row df channel (a per-key dict that can also carry the hc2_bm BM contrast dofs the scalar `_df_analytic_fallback` cannot represent — the unified surface's Wooldridge df column stays NaN until then); ImputationDiD lead-horizon `event_study_df` (needs the `_compute_lead_coefficients`/`_aggregate_event_study` return-path change plus bootstrap-clearing semantics for post rows); and the `inference_df` (DiD/MPD/SA/StackedDiD) vs `df_inference` (CS/dCDH) adapter-naming unification. NOTE (M-024): StackedDiD's `df_survey=None` on analytical fits is DELIBERATE - the container resolver's `df_inference` fallback not matching `inference_df` is load-bearing (normal-theory FLCI; see the REGISTRY StackedDiD M-024 Note and the pin in `tests/test_event_study_consumers.py`) - the unification must preserve that outcome or renegotiate it explicitly. | `diff_diff/wooldridge_results.py`, `diff_diff/imputation.py`, `diff_diff/results_base.py` | #variance-PR-C | Mid | Low |
diff --git a/diff_diff/guides/llms-practitioner.txt b/diff_diff/guides/llms-practitioner.txt
index ebcd8e33..a35b69bf 100644
--- a/diff_diff/guides/llms-practitioner.txt
+++ b/diff_diff/guides/llms-practitioner.txt
@@ -591,6 +591,25 @@ guidance = practitioner_next_steps(results)
# Returns context-aware suggestions for what to do next
```
+### Hand off to an MMM (marketing mix model)
+
+When the DiD measured a geo experiment on a marketing channel, export the
+result as MMM calibration input via `diff_diff.mmm` - the route depends on
+the grain the target MMM row represents:
+
+- The MMM row IS the campaign total (e.g. Meridian ROI priors):
+ `results.aggregate('total')` -> `to_meridian_roi_prior(aggregation_result=...,
+ spend=...)` - the total container arrives pre-scaled; `scale` is not accepted.
+- The MMM row's population/period grain differs from the total's span (e.g. a
+ PyMC-Marketing weekly lift row for a national model):
+ `results.aggregate('simple')` or `aggregate('group')` plus an explicit
+ `scale=` -> `to_pymc_marketing_lift_test(aggregation_result=..., scale=...,
+ x=..., delta_x=...)` - every lift-row field must describe the SAME
+ observation (same channel, population, and period span).
+
+Executed end-to-end examples: tutorials 29 (PyMC-Marketing lift test) and
+30 (Meridian ROI prior + calibration mask).
+
---
## Common Pitfalls
diff --git a/diff_diff/mmm.py b/diff_diff/mmm.py
index d4739c9b..2e7f5ea0 100644
--- a/diff_diff/mmm.py
+++ b/diff_diff/mmm.py
@@ -135,12 +135,16 @@
# Meridian prior parameters this exporter can target, with each one's Meridian
# default LogNormal(mu, sigma) per channel (verified against
-# meridian/model/prior_distribution.py at 1.7.0): roi_m is the return on a
+# meridian/model/prior_distribution.py at 1.7.0; execution-validated unchanged on
+# 1.8.0, 2026-08, and continuously checked by tests/test_mmm_interop_meridian.py's
+# defaults-drift canary): roi_m is the return on a
# channel's full spend (zero-spend counterfactual), mroi_m the marginal return.
# Channels without an experiment keep the default in the vector snippet.
_MERIDIAN_PARAM_DEFAULTS = {"roi_m": (0.2, 0.9), "mroi_m": (0.0, 0.5)}
-# The .to_code() templates are pinned against google-meridian 1.7.0 (2026-06); they
+# The .to_code() templates are pinned against google-meridian 1.7.0 (2026-06) and
+# execution-validated on 1.8.0 (2026-08: the generated snippets exec verbatim into a
+# ModelSpec that Meridian 1.8.0 accepts); they
# are convenience snippets, not a programmatic contract. Meridian's roi_m/mroi_m have
# batch shape n_media_channels; a scalar LogNormal broadcasts to EVERY channel, so the
# scalar template is gated behind an explicit single_channel opt-in.
@@ -817,7 +821,8 @@ def to_code(
roi_calibration_period: Optional[Union[str, np.ndarray]] = None,
full_model_window: bool = False,
) -> str:
- """Ready-to-paste Meridian snippet (channel- and time-scoped; 1.7.0 pinned).
+ """Ready-to-paste Meridian snippet (channel- and time-scoped; 1.7.0 pinned,
+ execution-validated on 1.8.0).
Meridian's ``roi_m``/``mroi_m`` prior has batch shape ``n_media_channels`` and
a scalar distribution broadcasts to EVERY media channel - a TV experiment's
diff --git a/docs/api/mmm.rst b/docs/api/mmm.rst
index 3d3063a5..f5bdbcb6 100644
--- a/docs/api/mmm.rst
+++ b/docs/api/mmm.rst
@@ -8,6 +8,14 @@ Assemble Marketing Mix Model (MMM) calibration inputs from experiment results.
Overview
--------
+.. seealso::
+
+ Tutorials 29 (:doc:`../tutorials/29_mmm_calibration_pymc`) and 30
+ (:doc:`../tutorials/30_mmm_calibration_meridian`) run both hand-offs
+ end-to-end against the real frameworks - a lift-test calibration in
+ PyMC-Marketing and an ROI-prior calibration in Meridian - including the
+ with/without-calibration posterior comparison.
+
MMM practitioners calibrate their models against experimental evidence. The two
dominant Python MMM frameworks consume that evidence in different shapes:
diff --git a/docs/dev-status.md b/docs/dev-status.md
index 337b8860..3a8bc55d 100644
--- a/docs/dev-status.md
+++ b/docs/dev-status.md
@@ -152,3 +152,24 @@ Linux.
covariates); logistic-regression overflow in edge cases (separate from the BLAS bug).
- **Long-term:** revert to the `@` operator when numpy ≥ 2.3 becomes the minimum supported
version.
+
+## MMM tutorial environments (notebooks 29 / 30)
+
+The two executed MMM calibration tutorials need frameworks that CANNOT share one
+environment: pymc-marketing 1.0 requires `arviz>=1.2,<2.0` while google-meridian
+requires `arviz<0.20`. Each notebook therefore has its own Python 3.12 venv and
+Jupyter kernel (TensorFlow does not support the 3.14 dev environment):
+
+- `~/.venvs/diffdiff-mmm-pymc` -> kernel `diffdiff-mmm-pymc`
+ (`pymc-marketing==1.0.0` + matplotlib/ipykernel/jupyter/nbconvert/pytest)
+- `~/.venvs/diffdiff-mmm-meridian` -> kernel `diffdiff-mmm-meridian`
+ (`google-meridian==1.8.0` + the same tooling)
+
+Both use the CI jobs' `diff_diff_dev.pth` shim (the repo root written into
+site-packages) instead of `pip install -e .` (the maturin build backend would
+demand a Rust build). Committed notebooks stay kernelspec-free; every local
+execution names the kernel explicitly (`--nbmake-kernel=...` /
+`nbconvert --ExecutePreprocessor.kernel_name=...`) - never the default
+`python3` kernel. Fragile edge: google-meridian pins an exact `tfp-nightly`
+build; the response protocol lives in `.github/workflows/mmm-interop.yml`'s
+header comment.
diff --git a/docs/doc-deps.yaml b/docs/doc-deps.yaml
index 7a054780..683af42e 100644
--- a/docs/doc-deps.yaml
+++ b/docs/doc-deps.yaml
@@ -1073,7 +1073,7 @@ sources:
- path: docs/methodology/REGISTRY.md
section: "MMM Calibration Export (interop)"
type: methodology
- note: "Lift-test scaling semantics, spend-weighted pooling, lognormal closed form (pinned to Meridian 1.7.0), sign/positivity policies."
+ note: "Lift-test scaling semantics, spend-weighted pooling, lognormal closed form (pinned to Meridian 1.7.0; execution-validated on 1.8.0), sign/positivity policies."
- path: docs/api/mmm.rst
type: api_reference
- path: README.md
@@ -1087,6 +1087,15 @@ sources:
- path: diff_diff/guides/llms-full.txt
section: "MMM Calibration Export"
type: user_guide
+ - path: diff_diff/guides/llms-practitioner.txt
+ section: "Step 8 MMM hand-off"
+ type: user_guide
+ - path: docs/tutorials/29_mmm_calibration_pymc.ipynb
+ type: tutorial
+ note: "pymc-marketing interop: executed lift-test calibration; drift guard tests/test_t29_mmm_calibration_pymc_drift.py"
+ - path: docs/tutorials/30_mmm_calibration_meridian.ipynb
+ type: tutorial
+ note: "meridian interop: executed ROI-prior calibration + to_code snippet; drift guard tests/test_t30_mmm_calibration_meridian_drift.py"
diff_diff/power.py:
drift_risk: low
@@ -1231,6 +1240,15 @@ sources:
type: user_guide
- path: docs/r_comparison.rst
type: user_guide
+ - path: diff_diff/guides/llms-practitioner.txt
+ section: "Step 8 MMM hand-off (aggregate routes)"
+ type: user_guide
+ - path: docs/tutorials/29_mmm_calibration_pymc.ipynb
+ type: tutorial
+ note: "aggregate('simple')+scale lift-row route taught end-to-end"
+ - path: docs/tutorials/30_mmm_calibration_meridian.ipynb
+ type: tutorial
+ note: "aggregate('total') container route taught end-to-end"
diff_diff/results_base.py:
drift_risk: low
diff --git a/docs/index.rst b/docs/index.rst
index 303ac99e..7b61e46c 100644
--- a/docs/index.rst
+++ b/docs/index.rst
@@ -78,7 +78,7 @@ Explore the Documentation
:link: tutorials/index
:link-type: doc
- 28 hands-on notebooks, from basic 2x2 DiD to survey-weighted and
+ 31 hands-on notebooks, from basic 2x2 DiD to survey-weighted and
spillover-aware designs.
.. grid-item-card:: User Guide
diff --git a/docs/methodology/REGISTRY.md b/docs/methodology/REGISTRY.md
index 917dbf84..31d84fb7 100644
--- a/docs/methodology/REGISTRY.md
+++ b/docs/methodology/REGISTRY.md
@@ -6537,7 +6537,7 @@ estimator-focused:
## MMM Calibration Export (interop)
-**Primary sources:** [PyMC-Marketing lift-test calibration](https://www.pymc-marketing.io/en/stable/notebooks/mmm/mmm_lift_test.html) (`MMM.add_lift_test_measurements` schema); [Google Meridian, "Set custom prior distributions using past experiments"](https://developers.google.com/meridian/docs/advanced-modeling/set-custom-priors-past-experiments) and the closed-form lognormal conversion in `meridian/model/prior_distribution.py` (`lognormal_dist_from_mean_std`, Meridian 1.7.0); the roi_m/mroi_m estimand definitions in [Meridian's ROI parameterization docs](https://developers.google.com/meridian/docs/advanced-modeling/roi-mroi-contribution-parameterizations); the `roi_calibration_period` mask shape/semantics contract in [Meridian's "Set the ROI calibration period"](https://developers.google.com/meridian/docs/user-guide/configure-model) (configure-model guide).
+**Primary sources:** [PyMC-Marketing lift-test calibration](https://www.pymc-marketing.io/en/stable/notebooks/mmm/mmm_lift_test.html) (`MMM.add_lift_test_measurements` schema); [Google Meridian, "Set custom prior distributions using past experiments"](https://developers.google.com/meridian/docs/advanced-modeling/set-custom-priors-past-experiments) and the closed-form lognormal conversion in `meridian/model/prior_distribution.py` (`lognormal_dist_from_mean_std`, Meridian 1.7.0; execution-validated on 1.8.0, 2026-08); the roi_m/mroi_m estimand definitions in [Meridian's ROI parameterization docs](https://developers.google.com/meridian/docs/advanced-modeling/roi-mroi-contribution-parameterizations); the `roi_calibration_period` mask shape/semantics contract in [Meridian's "Set the ROI calibration period"](https://developers.google.com/meridian/docs/user-guide/configure-model) (configure-model guide).
**Module:** `diff_diff/mmm.py`
@@ -6559,7 +6559,7 @@ estimator-focused:
*Meridian ROI prior (`to_meridian_roi_prior`):*
- Per experiment `roi = incremental_outcome / spend`, `roi_sd = incremental_outcome_se / spend * se_widening`. Multiple experiments pooled spend-weighted - the spend-weighted average ROI the [Meridian FAQ](https://developers.google.com/meridian/docs/faqs) suggests (citing sec 3.4 of Google's MMM calibration whitepaper): `roi_mean = sum(w_i * roi_i)`, `w_i = spend_i / sum(spend)`.
-- **Estimand (`parameter={roi_m,mroi_m}`, default roi_m):** `roi_m` is the return on the channel's full spend (a zero-spend/full-holdout estimand; `spend` = total channel spend); `mroi_m` is the marginal return of a spend change (`spend` = the spend change). Under saturation these differ, so the caller declares which their `incremental_outcome` measures; the returned prior and `.to_code()` target that parameter with its own Meridian default for non-experiment channels (roi_m `LogNormal(0.2, 0.9)`, mroi_m `LogNormal(0.0, 0.5)`, both verified at 1.7.0). Sign, aggregation, and population are the caller's responsibility - a go-dark contract's numerator is `Y_exposed - Y_zero`, which the caller supplies with the correct sign.
+- **Estimand (`parameter={roi_m,mroi_m}`, default roi_m):** `roi_m` is the return on the channel's full spend (a zero-spend/full-holdout estimand; `spend` = total channel spend); `mroi_m` is the marginal return of a spend change (`spend` = the spend change). Under saturation these differ, so the caller declares which their `incremental_outcome` measures; the returned prior and `.to_code()` target that parameter with its own Meridian default for non-experiment channels (roi_m `LogNormal(0.2, 0.9)`, mroi_m `LogNormal(0.0, 0.5)`, both verified at 1.7.0, execution-validated unchanged on 1.8.0, and continuously checked by the `tests/test_mmm_interop_meridian.py` defaults-drift canary). Sign, aggregation, and population are the caller's responsibility - a go-dark contract's numerator is `Y_exposed - Y_zero`, which the caller supplies with the correct sign.
- Lognormal conversion matches Meridian's `lognormal_dist_from_mean_std` exactly: `sigma = sqrt(log1p((s/m)^2))`, `mu = ln(m) - log1p((s/m)^2)/2` (`log1p` keeps a ~1e-8 relative SE from rounding `sigma` to 0).
- A non-positive pooled ROI mean raises (lognormal positivity), pointing to pooling, wider priors, or Meridian's contribution/coefficient parameterizations. `spend`, `incremental_outcome_se`, and `se_widening` must be finite and positive; `incremental_outcome` finite.
- **Note:** Pooled `roi_sd = sqrt(sum((w_i * sd_i)^2))` treats experiments as independent (no covariance term). Ignoring covariance can misstate the pooled sd - it is anti-conservative when the net weighted covariance is positive, which experiments sharing control units or overlapping windows typically induce - and the direction cannot be determined from marginal SEs alone; the docstring instructs users to widen via `se_widening` (a conservative heuristic, not an exact correction).
diff --git a/docs/references.rst b/docs/references.rst
index a3611ff4..f92547cb 100644
--- a/docs/references.rst
+++ b/docs/references.rst
@@ -371,7 +371,7 @@ MMM Calibration Interop
- **Google.** "Set Custom Prior Distributions Using Past Experiments." *Google Meridian documentation*. https://developers.google.com/meridian/docs/advanced-modeling/set-custom-priors-past-experiments
- Documents the experiment-to-ROI-prior calibration workflow and caveats; the lognormal ``(mu, sigma)`` closed form emitted by ``to_meridian_roi_prior`` matches Meridian's ``prior_distribution.lognormal_dist_from_mean_std`` (verified against Meridian 1.7.0 source).
+ Documents the experiment-to-ROI-prior calibration workflow and caveats; the lognormal ``(mu, sigma)`` closed form emitted by ``to_meridian_roi_prior`` matches Meridian's ``prior_distribution.lognormal_dist_from_mean_std`` (verified against Meridian 1.7.0 source; execution-validated on 1.8.0, 2026-08).
- **Google.** "Set the ROI Calibration Period." *Google Meridian documentation (configure-model guide)*. https://developers.google.com/meridian/docs/user-guide/configure-model
diff --git a/docs/tutorials/28_rdd_scholarship_illusion.ipynb b/docs/tutorials/28_rdd_scholarship_illusion.ipynb
index f3dd143a..73a047db 100644
--- a/docs/tutorials/28_rdd_scholarship_illusion.ipynb
+++ b/docs/tutorials/28_rdd_scholarship_illusion.ipynb
@@ -471,7 +471,7 @@
"An RD estimate is only as credible as its design checks. Three standard ones, all\n",
"runnable with the estimator itself.\n",
"\n",
- "### 4a. Covariate balance\n",
+ "### Diagnostic A: covariate balance\n",
"\n",
"Continuity has a testable implication: **predetermined covariates cannot jump at the\n",
"cutoff** - the scholarship offer cannot change your parents' income or your high-school\n",
@@ -527,7 +527,7 @@
"absence supports it, but no balance test can prove continuity of the (unobservable)\n",
"potential-outcome functions.\n",
"\n",
- "### 4b. Placebo cutoffs\n",
+ "### Diagnostic B: placebo cutoffs\n",
"\n",
"At any threshold where nothing is assigned, the estimator should find nothing.\n",
"Practical guides examine one or more artificial cutoffs across the support;\n",
@@ -585,7 +585,7 @@
"cutoffs cannot rule out discontinuities everywhere else - but finding none where none\n",
"should exist is exactly the supporting evidence this check exists to collect.\n",
"\n",
- "### 4c. Bandwidth sensitivity\n",
+ "### Diagnostic C: bandwidth sensitivity\n",
"\n",
"The bandwidth is the estimator's main tuning parameter, so referees will ask: does the\n",
"answer depend on it? Overriding the selector with fixed `h` values maps out the\n",
diff --git a/docs/tutorials/29_mmm_calibration_pymc.ipynb b/docs/tutorials/29_mmm_calibration_pymc.ipynb
new file mode 100644
index 00000000..bfa2c12c
--- /dev/null
+++ b/docs/tutorials/29_mmm_calibration_pymc.ipynb
@@ -0,0 +1,878 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "id": "964fea06",
+ "metadata": {},
+ "source": [
+ "# MMM Calibration with PyMC-Marketing Lift Tests\n",
+ "\n",
+ "Marketing mix models fit on observational spend history inherit its biases: when\n",
+ "search budgets chase demand, the model credits search with sales that demand would\n",
+ "have produced anyway. A geo experiment breaks that circularity - and diff-diff's\n",
+ "MMM exporters carry the experiment's result into\n",
+ "[PyMC-Marketing](https://www.pymc-marketing.io/) as a *lift test* with one function\n",
+ "call.\n",
+ "\n",
+ "This tutorial runs the full pipeline **end-to-end with the real frameworks**:\n",
+ "\n",
+ "1. simulate a geo-week world where search spend chased demand historically, then a\n",
+ " staggered geo **spend-boost experiment** was run under a budget freeze;\n",
+ "2. measure the boost's incremental sales with `CallawaySantAnna` -\n",
+ " `results.aggregate('total')` for the campaign headline, and\n",
+ " `results.aggregate('simple')` for the per-geo-week ATT the lift row needs;\n",
+ "3. export the measurement with `to_pymc_marketing_lift_test` via the container\n",
+ " route - the only number you supply is the count of treated geos (the\n",
+ " grain-reconciliation `scale`);\n",
+ "4. fit the national PyMC-Marketing MMM **without and with** the lift test and watch\n",
+ " the search-ROI posterior move to the truth and tighten.\n",
+ "\n",
+ "**Requirements.** This notebook needs `pymc-marketing==1.0.0` (Python 3.12), which\n",
+ "is *not* a diff-diff dependency: `pip install \"pymc-marketing==1.0.0\" diff-diff`.\n",
+ "It cannot share an environment with the Meridian tutorial (notebook 30) - the two\n",
+ "frameworks pin incompatible `arviz` ranges - so each tutorial gets its own\n",
+ "environment/kernel. Runtime is a few minutes (two small MCMC fits).\n",
+ "\n",
+ "**Reproducibility.** MCMC draws depend on framework and transitive-dependency\n",
+ "versions and on the platform. The committed outputs were produced with\n",
+ "pymc-marketing 1.0.0; re-executing under other versions may shift posterior\n",
+ "numbers. The `assert` cells check the qualitative claims (sampler health, and\n",
+ "that calibration moves the posterior toward the simulated truth and tightens it)\n",
+ "with slack margins - if one trips on your machine, it means exactly that claim\n",
+ "did not reproduce in your environment."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "1b76eaeb",
+ "metadata": {},
+ "source": [
+ "## A world where spend chased demand\n",
+ "\n",
+ "Fifteen geos, 104 weeks, two channels. The **search** channel is the experiment\n",
+ "channel and is deliberately simple: zero carryover (an immediate-response channel)\n",
+ "and a *linear* true response - every search dollar returns `BETA_SEARCH = 2.0`\n",
+ "dollars of sales, so the true incremental ROI is 2.0 by construction. The **tv**\n",
+ "channel keeps adstock and saturation for realism.\n",
+ "\n",
+ "The bias mechanism: for the first 52 weeks, search budgets *chased seasonal\n",
+ "demand* (`spend = base * (1 + 0.5 * season)` while the same season moves sales\n",
+ "directly). From week 52 budgets were **frozen** - standard geo-test hygiene, and\n",
+ "what makes the lift test's baseline spend `x` a single well-defined number - and\n",
+ "from week 60 a staggered boost of \\$50/week rolled out across nine treated geos\n",
+ "(cohorts at weeks 60, 68, 76), leaving six geos as controls."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "id": "dba18972",
+ "metadata": {
+ "execution": {
+ "iopub.execute_input": "2026-08-16T19:57:05.565006Z",
+ "iopub.status.busy": "2026-08-16T19:57:05.564773Z",
+ "iopub.status.idle": "2026-08-16T19:57:06.439770Z",
+ "shell.execute_reply": "2026-08-16T19:57:06.439363Z"
+ }
+ },
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "\n",
+ "\n",
+ "
\n",
+ " \n",
+ " \n",
+ " geo \n",
+ " week \n",
+ " sales \n",
+ " search_spend \n",
+ " tv_spend \n",
+ " first_treat \n",
+ " \n",
+ " \n",
+ " \n",
+ " \n",
+ " 0 \n",
+ " 0 \n",
+ " 0 \n",
+ " 3543.352739 \n",
+ " 100.000000 \n",
+ " 203.578696 \n",
+ " 60 \n",
+ " \n",
+ " \n",
+ " 1 \n",
+ " 0 \n",
+ " 1 \n",
+ " 3687.103795 \n",
+ " 106.026834 \n",
+ " 212.798263 \n",
+ " 60 \n",
+ " \n",
+ " \n",
+ " 2 \n",
+ " 0 \n",
+ " 2 \n",
+ " 3753.001771 \n",
+ " 111.965783 \n",
+ " 209.345368 \n",
+ " 60 \n",
+ " \n",
+ " \n",
+ " 3 \n",
+ " 0 \n",
+ " 3 \n",
+ " 3840.767622 \n",
+ " 117.730244 \n",
+ " 207.410011 \n",
+ " 60 \n",
+ " \n",
+ " \n",
+ " 4 \n",
+ " 0 \n",
+ " 4 \n",
+ " 3905.212324 \n",
+ " 123.236159 \n",
+ " 207.098347 \n",
+ " 60 \n",
+ " \n",
+ " \n",
+ "
\n",
+ "
\n",
+ "\n",
+ "
\n",
+ " \n",
+ " \n",
+ " channel \n",
+ " x \n",
+ " delta_x \n",
+ " delta_y \n",
+ " sigma \n",
+ " \n",
+ " \n",
+ " \n",
+ " \n",
+ " 0 \n",
+ " search \n",
+ " 1500.0 \n",
+ " 450.0 \n",
+ " 932.760714 \n",
+ " 82.569488 \n",
+ " \n",
+ " \n",
+ "
\n",
+ "
"
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "fig, ax = plt.subplots(figsize=(8, 4.5))\n",
+ "bins = np.linspace(0.5, 6.5, 80)\n",
+ "ax.hist(np.asarray(roi_plain).ravel(), bins=bins, density=True, alpha=0.55,\n",
+ " label=\"without lift test\", color=\"#c44e52\")\n",
+ "ax.hist(np.asarray(roi_cal).ravel(), bins=bins, density=True, alpha=0.65,\n",
+ " label=\"with lift test\", color=\"#4c72b0\")\n",
+ "ax.axvline(BETA_SEARCH, color=\"black\", lw=2, label=\"simulated truth (2.0)\")\n",
+ "did_roi = df_lift.delta_y[0] / df_lift.delta_x[0]\n",
+ "ax.axvline(did_roi, color=\"#55a868\", lw=2, ls=\"--\", label=\"DiD lift measurement\")\n",
+ "ax.set_xlabel(\"campaign incremental ROI (sales per boost dollar)\")\n",
+ "ax.set_ylabel(\"posterior density\")\n",
+ "ax.set_title(\"What the lift test buys: search-ROI posterior, before vs after calibration\")\n",
+ "ax.legend()\n",
+ "plt.tight_layout()\n",
+ "plt.show()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "4326b1bb",
+ "metadata": {},
+ "source": [
+ "The uncalibrated posterior is wide and biased high - the model still credits\n",
+ "search with demand it chased historically. The calibrated posterior concentrates\n",
+ "near the truth: the lift likelihood anchors the response curve at the one spend\n",
+ "change we actually experimented on."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "a83be6fe",
+ "metadata": {},
+ "source": [
+ "## Acceptance checks\n",
+ "\n",
+ "Hard checks on the claims this notebook makes: sampler health first, then the\n",
+ "calibration claims with slack margins against draw noise. The saturation\n",
+ "parameters sit on a ridge (`beta` and `lam` trade off along the fitted curve), so\n",
+ "the demanding sampler settings above (`target_accept=0.999`, deeper trees, long\n",
+ "tuning) are what it takes to traverse that geometry cleanly - with them, BOTH fits\n",
+ "sample with zero divergences, zero maximum-tree-depth events, and whole-model\n",
+ "maximum R-hat under PyMC's own 1.01 warning bar, so the wide uncalibrated\n",
+ "posterior is a faithful picture of the model's uncertainty rather than a sampling\n",
+ "artifact. ArviZ's chain-aware `az.rhat` and bulk/tail `az.ess` are additionally\n",
+ "checked on the derived campaign-ROI posterior - the quantity the notebook\n",
+ "reports."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "id": "64294128",
+ "metadata": {
+ "execution": {
+ "iopub.execute_input": "2026-08-16T20:03:05.883166Z",
+ "iopub.status.busy": "2026-08-16T20:03:05.883107Z",
+ "iopub.status.idle": "2026-08-16T20:03:06.202559Z",
+ "shell.execute_reply": "2026-08-16T20:03:06.202135Z"
+ }
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "plain: divergences=0, treedepth events=0, whole-model max Rhat=1.0069, Rhat(campaign ROI)=1.0016, ESS bulk=458 tail=366\n",
+ "calibrated: divergences=0, treedepth events=0, whole-model max Rhat=1.0062, Rhat(campaign ROI)=1.0014, ESS bulk=2104 tail=1197\n",
+ "error: 1.52 -> 0.21; 90% width: 3.59 -> 0.59\n",
+ "all acceptance checks passed\n"
+ ]
+ }
+ ],
+ "source": [
+ "import arviz as az\n",
+ "\n",
+ "ess_floors = {\"plain\": 300, \"calibrated\": 800}\n",
+ "for name, m, roi in [(\"plain\", mmm_plain, roi_plain), (\"calibrated\", mmm_cal, roi_cal)]:\n",
+ " div = int(m.idata.sample_stats[\"diverging\"].sum())\n",
+ " mtd = int(m.idata.sample_stats[\"reached_max_treedepth\"].sum())\n",
+ " with warnings.catch_warnings():\n",
+ " # constant/deterministic posterior variables have undefined R-hat;\n",
+ " # arviz emits RuntimeWarnings computing them - excluded below.\n",
+ " warnings.simplefilter(\"ignore\", RuntimeWarning)\n",
+ " rh_all = az.rhat(m.idata.posterior)\n",
+ " finite_rhats = np.concatenate(\n",
+ " [np.asarray(rh_all[v]).ravel() for v in rh_all.data_vars]\n",
+ " )\n",
+ " max_rhat_model = float(finite_rhats[np.isfinite(finite_rhats)].max())\n",
+ " rhat_roi = float(az.rhat(roi).values)\n",
+ " ess_bulk = float(az.ess(roi).values)\n",
+ " ess_tail = float(az.ess(roi, method=\"tail\").values)\n",
+ " print(f\"{name}: divergences={div}, treedepth events={mtd}, \"\n",
+ " f\"whole-model max Rhat={max_rhat_model:.4f}, \"\n",
+ " f\"Rhat(campaign ROI)={rhat_roi:.4f}, ESS bulk={ess_bulk:.0f} tail={ess_tail:.0f}\")\n",
+ " assert div == 0, f\"{name} fit: {div} divergences - sampling is unhealthy here\"\n",
+ " assert mtd == 0, f\"{name} fit: {mtd} maximum-tree-depth events - raise max_treedepth\"\n",
+ " assert max_rhat_model < 1.01, (\n",
+ " f\"{name} fit: some model parameter has not converged (max Rhat {max_rhat_model:.4f})\"\n",
+ " )\n",
+ " assert rhat_roi < 1.01, f\"{name} fit: campaign-ROI chains disagree (Rhat {rhat_roi:.4f})\"\n",
+ " assert ess_bulk >= ess_floors[name], (\n",
+ " f\"{name} fit: campaign-ROI bulk ESS too low ({ess_bulk:.0f})\"\n",
+ " )\n",
+ " assert ess_tail >= 300, f\"{name} fit: campaign-ROI tail ESS too low ({ess_tail:.0f})\"\n",
+ "\n",
+ "flat_plain = np.asarray(roi_plain).ravel()\n",
+ "flat_cal = np.asarray(roi_cal).ravel()\n",
+ "err_plain = abs(float(flat_plain.mean()) - BETA_SEARCH)\n",
+ "err_cal = abs(float(flat_cal.mean()) - BETA_SEARCH)\n",
+ "width_plain = float(np.quantile(flat_plain, 0.95) - np.quantile(flat_plain, 0.05))\n",
+ "width_cal = float(np.quantile(flat_cal, 0.95) - np.quantile(flat_cal, 0.05))\n",
+ "print(f\"error: {err_plain:.2f} -> {err_cal:.2f}; 90% width: {width_plain:.2f} -> {width_cal:.2f}\")\n",
+ "\n",
+ "assert err_cal < err_plain, (\n",
+ " \"calibrated ROI mean is not closer to the simulated truth - the lift test did \"\n",
+ " \"not correct the bias in this environment\"\n",
+ ")\n",
+ "assert width_cal < 0.8 * width_plain, (\n",
+ " \"calibrated 90% interval is not clearly narrower - the lift test did not add \"\n",
+ " \"information in this environment\"\n",
+ ")\n",
+ "print(\"all acceptance checks passed\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "f7ee1dbf",
+ "metadata": {},
+ "source": [
+ "## Where to go next\n",
+ "\n",
+ "- **Meridian**: notebook 30 runs the same hand-off into Google Meridian's ROI\n",
+ " priors, using a geo-holdout *launch* experiment (the other classic design) and\n",
+ " `to_meridian_roi_prior` + `meridian_calibration_mask`.\n",
+ "- **API reference**: the `diff_diff.mmm` module documents both exporters'\n",
+ " contracts - including the two input routes (a pre-scaled `aggregate('total')`\n",
+ " container, or `aggregate('simple')`/`aggregate('group')` plus an explicit\n",
+ " `scale=`) and the additive-outcome caveat.\n",
+ "- **Estimator choices**: the DiD side works with any diff-diff estimator that\n",
+ " supports `aggregate('total')` (CallawaySantAnna, EfficientDiD, ImputationDiD,\n",
+ " TwoStageDiD)."
+ ]
+ }
+ ],
+ "metadata": {
+ "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.13"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/docs/tutorials/30_mmm_calibration_meridian.ipynb b/docs/tutorials/30_mmm_calibration_meridian.ipynb
new file mode 100644
index 00000000..927af01d
--- /dev/null
+++ b/docs/tutorials/30_mmm_calibration_meridian.ipynb
@@ -0,0 +1,970 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "id": "d86bf740",
+ "metadata": {},
+ "source": [
+ "# MMM Calibration with Meridian ROI Priors\n",
+ "\n",
+ "Google [Meridian](https://developers.google.com/meridian) calibrates a marketing\n",
+ "mix model by putting an experiment-informed prior on a channel's ROI. diff-diff\n",
+ "closes the loop from a geo experiment to that prior:\n",
+ "`results.aggregate('total')` gives the campaign's total incremental outcome in\n",
+ "the estimator's own units, `to_meridian_roi_prior` turns it into Meridian's\n",
+ "lognormal `roi_m` prior, and `meridian_calibration_mask` scopes the prior to the\n",
+ "experiment window. The prior object even generates a ready-to-paste Meridian\n",
+ "snippet - which this notebook executes verbatim.\n",
+ "\n",
+ "The experiment here is the second classic geo design (notebook 29 showed the\n",
+ "spend-boost lift test): a **geo-holdout launch**. The channel runs *only* in\n",
+ "launched geos; never-launched geos are the holdout. That matters for the\n",
+ "estimand: Meridian's `roi_m` is the return on a channel's *full spend against a\n",
+ "zero-spend counterfactual*, and a launch experiment measures exactly that - the\n",
+ "DiD total *is* the incremental outcome versus zero spend.\n",
+ "\n",
+ "**Requirements.** This notebook needs `google-meridian==1.8.0` (Python 3.12,\n",
+ "TensorFlow-based), which is *not* a diff-diff dependency:\n",
+ "`pip install \"google-meridian==1.8.0\" diff-diff`. It cannot share an environment\n",
+ "with the PyMC-Marketing tutorial (notebook 29) - the two frameworks pin\n",
+ "incompatible `arviz` ranges - so each tutorial gets its own environment/kernel.\n",
+ "Runtime is a few minutes (two small MCMC fits).\n",
+ "\n",
+ "**Reproducibility.** MCMC draws depend on framework and transitive-dependency\n",
+ "versions and on the platform. The committed outputs were produced with\n",
+ "google-meridian 1.8.0; re-executing under other versions may shift posterior\n",
+ "numbers. The `assert` cells check the qualitative claims with slack margins - if\n",
+ "one trips on your machine, it means exactly that claim did not reproduce in your\n",
+ "environment."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "7e0b5e03",
+ "metadata": {},
+ "source": [
+ "## A launch the national model cannot see\n",
+ "\n",
+ "Twelve geos, 104 weeks. The **search** channel launches in eight geos (staggered\n",
+ "cohorts at weeks 30, 44, 58) at a constant \\$150/week and never runs in the four\n",
+ "holdout geos; its true response is linear with `BETA_SEARCH = 2.5` - the true ROI.\n",
+ "The **tv** channel runs everywhere with adstock and saturation.\n",
+ "\n",
+ "Two realistic wrinkles create the calibration problem:\n",
+ "\n",
+ "- the launches ride a **demand ramp** (the brand expanded during its growth\n",
+ " season), and\n",
+ "- weekly **national shocks** (macro, weather, promos) hit every geo at once.\n",
+ "\n",
+ "The geo-level DiD is immune to both - common time effects difference out against\n",
+ "the holdout geos. But the MMM below is *national* (as most production MMMs are),\n",
+ "and aggregation destroys the holdout contrast: nationally, search spend simply\n",
+ "ramps up while sales ramp up, so the model overcredits search - and the common\n",
+ "shocks it must absorb make its evidence noisy too. The experiment's information\n",
+ "survives only if we carry it in as a prior."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "id": "c0d3841c",
+ "metadata": {
+ "execution": {
+ "iopub.execute_input": "2026-08-16T19:52:58.412094Z",
+ "iopub.status.busy": "2026-08-16T19:52:58.412020Z",
+ "iopub.status.idle": "2026-08-16T19:52:59.370306Z",
+ "shell.execute_reply": "2026-08-16T19:52:59.369837Z"
+ }
+ },
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "\n",
+ "\n",
+ "
\n",
+ " \n",
+ " \n",
+ " geo \n",
+ " week \n",
+ " sales \n",
+ " search_spend \n",
+ " tv_spend \n",
+ " first_treat \n",
+ " \n",
+ " \n",
+ " \n",
+ " \n",
+ " 0 \n",
+ " 0 \n",
+ " 0 \n",
+ " 5288.883191 \n",
+ " 0.0 \n",
+ " 274.984892 \n",
+ " 30 \n",
+ " \n",
+ " \n",
+ " 1 \n",
+ " 0 \n",
+ " 1 \n",
+ " 5442.518730 \n",
+ " 0.0 \n",
+ " 265.642663 \n",
+ " 30 \n",
+ " \n",
+ " \n",
+ " 2 \n",
+ " 0 \n",
+ " 2 \n",
+ " 5560.941283 \n",
+ " 0.0 \n",
+ " 272.251908 \n",
+ " 30 \n",
+ " \n",
+ " \n",
+ " 3 \n",
+ " 0 \n",
+ " 3 \n",
+ " 5546.233027 \n",
+ " 0.0 \n",
+ " 265.027138 \n",
+ " 30 \n",
+ " \n",
+ " \n",
+ " 4 \n",
+ " 0 \n",
+ " 4 \n",
+ " 5442.605704 \n",
+ " 0.0 \n",
+ " 274.290276 \n",
+ " 30 \n",
+ " \n",
+ " \n",
+ "
\n",
+ "
, tau_g_excl_baseline=, beta_m=, beta_rf=, beta_om=, beta_orf=, eta_m=, eta_rf=, eta_om=, eta_orf=, gamma_c=, gamma_n=, xi_c=, xi_n=, alpha_m=, alpha_rf=, alpha_om=, alpha_orf=, ec_m=, ec_rf=, ec_om=, ec_orf=, slope_m=, slope_rf=, slope_om=, slope_orf=, sigma=, roi_m=, roi_rf=, mroi_m=, mroi_rf=, contribution_m=, contribution_rf=, contribution_om=, contribution_orf=, contribution_n=), media_effects_dist='log_normal', hill_before_adstock=False, max_lag=8, unique_sigma_for_each_geo=False, media_prior_type='roi', rf_prior_type=None, paid_media_prior_type=None, roi_calibration_period=array([[False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [False, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True],\n",
+ " [ True, True]]), rf_roi_calibration_period=None, organic_media_prior_type='contribution', organic_rf_prior_type='contribution', non_media_treatments_prior_type='contribution', non_media_baseline_values=None, knots=None, baseline_geo=None, holdout_id=None, control_population_scaling_id=None, non_media_population_scaling_id=None, adstock_decay_spec='geometric', saturation_spec='hill', enable_aks=False)"
+ ]
+ },
+ "execution_count": 5,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "namespace = {}\n",
+ "exec(code, namespace)\n",
+ "model_spec_calibrated = namespace[\"model_spec\"]\n",
+ "assert np.array_equal(np.asarray(model_spec_calibrated.roi_calibration_period), mask)\n",
+ "model_spec_calibrated"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "1df0ff91",
+ "metadata": {},
+ "source": [
+ "## Fit the national Meridian model - default vs calibrated priors\n",
+ "\n",
+ "The national data is the geo panel summed across geos. We build Meridian's\n",
+ "`InputData` with its DataFrame builder and fit twice: once with the default\n",
+ "`ModelSpec()` (Meridian's default `roi_m` prior, `LogNormal(0.2, 0.9)`), once\n",
+ "with the calibrated spec from the generated snippet."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 6,
+ "id": "8f903e93",
+ "metadata": {
+ "execution": {
+ "iopub.execute_input": "2026-08-16T19:53:02.437982Z",
+ "iopub.status.busy": "2026-08-16T19:53:02.437881Z",
+ "iopub.status.idle": "2026-08-16T19:54:21.815732Z",
+ "shell.execute_reply": "2026-08-16T19:54:21.815009Z"
+ }
+ },
+ "outputs": [
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/data/input_data_builder.py:706: UserWarning: The `population` argument is ignored in a nationally aggregated model. It will be reset to [1, 1, ..., 1]\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/data/input_data.py:575: UserWarning: Revenue from the `kpi` data is used when `kpi_type`=`revenue`. `revenue_per_kpi` is ignored.\n",
+ " warnings.warn(\n",
+ "WARNING: All log messages before absl::InitializeLog() is called are written to STDERR\n",
+ "I0000 00:00:1786909982.522983 32840500 service.cc:153] XLA service 0x912c26b00 initialized for platform Host (this does not guarantee that XLA will be used). Devices:\n",
+ "I0000 00:00:1786909982.522997 32840500 service.cc:161] StreamExecutor [0]: Host, Default Version (Driver: 0.0.0; Runtime: 0.0.0; Toolkit: 0.0.0; DNN: 0.0.0)\n",
+ "I0000 00:00:1786909982.538321 32840500 device_compiler.h:208] Compiled cluster using XLA! This line is logged at most once for the lifetime of the process.\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/model.py:103: UserWarning: In a nationally aggregated model, the `media_effects_dist` will be reset to `normal`.\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. tau_g_excl_baseline has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. eta_m has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. eta_rf has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. eta_om has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. eta_orf has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. xi_c has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. xi_n has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n"
+ ]
+ },
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "I0000 00:00:1786910002.193695 32840500 dump_mlir_util.cc:269] disabling MLIR crash reproducer, set env var `MLIR_CRASH_REPRODUCER_DIRECTORY` to enable.\n"
+ ]
+ },
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "W0000 00:00:1786910002.596307 32840500 assert_op.cc:39] Ignoring Assert operator mcmc_retry_init/assert_equal_1/Assert/AssertGuard/Assert\n"
+ ]
+ },
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "W0000 00:00:1786910010.912318 32840719 cpp_gen_intrinsics.cc:74] Empty bitcode string provided for eigen. Optimizations relying on this IR will be disabled.\n"
+ ]
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "WARNING:tensorflow:5 out of the last 5 calls to .f at 0x12c2f3740> triggered tf.function retracing. Tracing is expensive and the excessive number of tracings could be due to (1) creating @tf.function repeatedly in a loop, (2) passing tensors with different shapes, (3) passing Python objects instead of tensors. For (1), please define your @tf.function outside of the loop. For (2), @tf.function has reduce_retracing=True option that can avoid unnecessary retracing. For (3), please refer to https://www.tensorflow.org/guide/function#controlling_retracing and https://www.tensorflow.org/api_docs/python/tf/function for more details.\n"
+ ]
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "WARNING:tensorflow:6 out of the last 6 calls to .f at 0x12c2f3ec0> triggered tf.function retracing. Tracing is expensive and the excessive number of tracings could be due to (1) creating @tf.function repeatedly in a loop, (2) passing tensors with different shapes, (3) passing Python objects instead of tensors. For (1), please define your @tf.function outside of the loop. For (2), @tf.function has reduce_retracing=True option that can avoid unnecessary retracing. For (3), please refer to https://www.tensorflow.org/guide/function#controlling_retracing and https://www.tensorflow.org/api_docs/python/tf/function for more details.\n"
+ ]
+ },
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/arviz/data/inference_data.py:157: UserWarning: trace group is not defined in the InferenceData scheme\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/arviz/data/inference_data.py:1647: UserWarning: trace group is not defined in the InferenceData scheme\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/model.py:103: UserWarning: In a nationally aggregated model, the `media_effects_dist` will be reset to `normal`.\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. tau_g_excl_baseline has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. eta_m has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. eta_rf has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. eta_om has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. eta_orf has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. xi_c has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/meridian/model/prior_distribution.py:1391: UserWarning: Hierarchical distribution parameters must be deterministically zero for national models. xi_n has been automatically set to Deterministic(0).\n",
+ " warnings.warn(\n"
+ ]
+ },
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "W0000 00:00:1786910040.502335 32840500 assert_op.cc:39] Ignoring Assert operator mcmc_retry_init/assert_equal_1/Assert/AssertGuard/Assert\n"
+ ]
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "both fits done\n"
+ ]
+ },
+ {
+ "name": "stderr",
+ "output_type": "stream",
+ "text": [
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/arviz/data/inference_data.py:157: UserWarning: trace group is not defined in the InferenceData scheme\n",
+ " warnings.warn(\n",
+ "~/.venvs/diffdiff-mmm-meridian/lib/python3.12/site-packages/arviz/data/inference_data.py:1647: UserWarning: trace group is not defined in the InferenceData scheme\n",
+ " warnings.warn(\n"
+ ]
+ }
+ ],
+ "source": [
+ "nat = panel.groupby(\"week\").agg(\n",
+ " sales=(\"sales\", \"sum\"),\n",
+ " search_spend=(\"search_spend\", \"sum\"),\n",
+ " tv_spend=(\"tv_spend\", \"sum\"),\n",
+ ").reset_index()\n",
+ "nat[\"time\"] = nat[\"week\"].map(dict(enumerate(times)).get)\n",
+ "nat[\"geo\"] = \"national\"\n",
+ "nat[\"population\"] = 12_000_000.0\n",
+ "\n",
+ "from meridian.data import data_frame_input_data_builder as dfb\n",
+ "from meridian.model import model as m_model, spec as m_spec\n",
+ "\n",
+ "builder = (\n",
+ " dfb.DataFrameInputDataBuilder(kpi_type=\"revenue\")\n",
+ " .with_kpi(nat, kpi_col=\"sales\")\n",
+ " .with_population(nat)\n",
+ " .with_media(\n",
+ " nat,\n",
+ " media_cols=[\"search_spend\", \"tv_spend\"],\n",
+ " media_spend_cols=[\"search_spend\", \"tv_spend\"],\n",
+ " media_channels=[\"search\", \"tv\"],\n",
+ " )\n",
+ ")\n",
+ "input_data = builder.build()\n",
+ "\n",
+ "\n",
+ "def fit_meridian(model_spec):\n",
+ " mer = m_model.Meridian(input_data=input_data, model_spec=model_spec)\n",
+ " mer.sample_posterior(n_chains=4, n_adapt=500, n_burnin=500, n_keep=500, seed=SEED)\n",
+ " return mer\n",
+ "\n",
+ "\n",
+ "mer_default = fit_meridian(m_spec.ModelSpec())\n",
+ "mer_calibrated = fit_meridian(model_spec_calibrated)\n",
+ "print(\"both fits done\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "309f825d",
+ "metadata": {},
+ "source": [
+ "## Read out the ROI posterior\n",
+ "\n",
+ "Meridian's `Analyzer.roi()` returns the posterior ROI per channel in KPI units."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "id": "62cf511d",
+ "metadata": {
+ "execution": {
+ "iopub.execute_input": "2026-08-16T19:54:21.817640Z",
+ "iopub.status.busy": "2026-08-16T19:54:21.817543Z",
+ "iopub.status.idle": "2026-08-16T19:54:22.579406Z",
+ "shell.execute_reply": "2026-08-16T19:54:22.578982Z"
+ }
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ " default prior: ROI mean 3.29, 90% interval (3.02, 3.57)\n",
+ " calibrated prior: ROI mean 2.54, 90% interval (2.47, 2.62)\n",
+ " truth: 2.50\n",
+ " DiD measurement: 2.49 ± 0.05\n"
+ ]
+ }
+ ],
+ "source": [
+ "import arviz as az\n",
+ "from meridian.analysis import analyzer as m_analyzer\n",
+ "\n",
+ "\n",
+ "def roi_posterior(mer):\n",
+ " roi = np.asarray(m_analyzer.Analyzer(mer).roi()) # (chains, draws, channels)\n",
+ " return roi[..., 0] # search is channel 0\n",
+ "\n",
+ "\n",
+ "roi_default = roi_posterior(mer_default)\n",
+ "roi_calibrated = roi_posterior(mer_calibrated)\n",
+ "\n",
+ "for name, roi in [(\"default prior\", roi_default), (\"calibrated prior\", roi_calibrated)]:\n",
+ " flat = roi.ravel()\n",
+ " lo, hi = np.quantile(flat, [0.05, 0.95])\n",
+ " print(f\"{name:>17}: ROI mean {flat.mean():.2f}, 90% interval ({lo:.2f}, {hi:.2f})\")\n",
+ "print(f\"{'truth':>17}: {BETA_SEARCH:.2f}\")\n",
+ "print(f\"{'DiD measurement':>17}: {prior.roi_mean:.2f} ± {prior.roi_sd:.2f}\")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "id": "1ca40c00",
+ "metadata": {
+ "execution": {
+ "iopub.execute_input": "2026-08-16T19:54:22.580665Z",
+ "iopub.status.busy": "2026-08-16T19:54:22.580579Z",
+ "iopub.status.idle": "2026-08-16T19:54:22.711622Z",
+ "shell.execute_reply": "2026-08-16T19:54:22.711064Z"
+ }
+ },
+ "outputs": [
+ {
+ "data": {
+ "image/png": "iVBORw0KGgoAAAANSUhEUgAAAxYAAAG4CAYAAADYN3EQAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvJkbTWQAAAAlwSFlzAAAPYQAAD2EBqD+naQAAbDJJREFUeJzt3QecE9X6//Fn6b1KlSIC0hUpKoKKIKBYwI4XFIELig1EUbAhFxUQC3YUFVAU9SJYUEBAQEVBaYoFFATpYKHXZTf/1/fc/+SX3c3uJptdspt83q9XCDuZJCdnJpPzzHnOmQSfz+czAAAAAIhAvkieDAAAAAAEFgAAAACyBT0WAAAAACJGYAEAAAAgYgQWAAAAACJGYAEAAAAgYgQWAAAAACJGYAEAAAAgYgQWAAAAACJGYIE8ISEhwW677TbLbTZs2ODK9sQTT0S7KLnaww8/7OoJabVt29bd4vl7FKvbbMGCBa7Ode+58cYb7aSTTsrhEsbWceOvv/6KdlEQ4v4NEFggR7333nvuwDN9+vQ0j5122mnusfnz56d5rEaNGnb22WfnWLl+/vln96OlwCAUn376qVsfse2xxx6zDz74INrFQIgNGt0mT54cdJ3WrVu7xxs3bkx9xjFvP/FupUqVsvPOO88++eSTdJ/z008/WY8ePezEE0+0woULW9WqVa179+5ueWoTJ050r7t06VLLjTim4XgjsECOatOmjbv/6quvUizfu3ev/fjjj1agQAFbtGhRisc2bdrkbt5zcyqwGD58eFiBhdZH1jzwwAN26NChXF990fgR/uyzz9wN4StSpIi9/fbbaZbre/3111+7x3PbNhs/frytWbMm28uE9HXo0MHefPNNe+ONN+yee+6xtWvX2qWXXmqzZ89Os+60adOsWbNmNm/ePOvVq5e9+OKL1qdPH3cCTMuDnSSL12Paueee647rugc8Bfz/A3KAzvTUqlUrTWDxzTffmM/ns6uvvjrNY97fORlY4Pg4cOCAFS9e3AWQuuH/HDx40IoVK2aFChXKtmo5duyYJScnZ+tr5madO3e2jz76yKXKnHDCCf7lCjYqVapkdevWtV27duWqbVawYEGLZ9HYR0855RTXA+G58sorrWHDhvbMM89Yp06d/MvXrVtn119/vZ188sn2xRdfWIUKFfyPDRgwwM455xz3+A8//ODWiVeHDx922y9fvnzZGrx7vxfI2+ixQI5TgLBixYoUZ6zVS9GoUSO76KKLbPHixe6HJvAxdS0rlSE1nXlRaoO6p/X8WbNmpXj8jz/+sFtuucXq1atnRYsWtfLly7vgJbBnQl3XWibnn3++v4s8vTxR5US/8MIL7v+BXeqpvfLKK1a7dm1XtpYtW9p3332XZp3Vq1fbVVddZeXKlXMH5BYtWriGUShUR2PHjnWfW89Vw+mmm25K0XAaNmyYO9jrbFugfv36uR+C77//PkUqybvvvmv33XefVa5c2R3QL7vsMtdblNqSJUvswgsvtNKlS7uGlVIJUvc0efnQ6g3617/+ZWXLlvUHh8HGWHj5/v/973/dj7y2V6tWrWzVqlXu8Zdfftnq1KnjPqvy2YP1LoVTLp2l1LYsU6aMW19nI9VQDCyPftgmTZrk38ZaPz3h1KHKr/122bJl7uyeyqrnpJevv3PnTneWVNtYn19pgypXeuN7tF94+57qPzNvvfWW+47otZs3b+4aUYHSGweQejuqvlW2YPT6gY22d955x71XyZIlXTpKkyZNXMMukBp2uoWqS5cu7jNrHwqkwOKaa66x/PnzB32e0qdUFu1z+i5269Yt4m22efNm69q1q9sHKlasaHfeeacdOXIkzXsHq1ttQ6V+6nilMqlsU6dOTfNc7zuT2XEwPc8995xbX59F308df1L3+GzZssV69+7t9j3v9V9//fUU6xw9etQeeughV059l/SZ1ehOndaa2T6q46G2kxrw+tzaZ+6///405d69e3eG391wNWjQwAWiqfe1MWPGuNfVsTwwqBCtr2OSjhGPP/542O/ppUzpu6bjtra1vgc33HBD0OBXPSWqey8V69Zbb3X1EOi3335zQZKOPfouV6tWze3Le/bsCemYFsq29o5z+v6q51npYdp/lHWQ3hgLfR+975fqTUGd3iuQylGiRAm3DXSCQMcFpZsh7+MUInKcGpfqhlYj0PsxVuNPP6S66SCotKhTTz3V/1j9+vXdgTd1T4a6qRU46CD07LPPuoPqxo0b/euqMa8UCB1cdZDVD9tLL73k3lc/ZjogqpFwxx13uOeroaAfGfHuU9OPwNatW23OnDnucwSjH+d9+/a5dXWg1Q/PFVdcYb///rv/DKXycxUs6cA8ZMgQ92OsMShqjLz//vt2+eWXZ1iPem39OOlHVeVfv369Pf/88y5oU53pfXTg//jjj12jVA101ZO6+5V+MWLEiDSNwEcffdSV995773WNWf34X3DBBbZy5Ur3oyCff/65CwD1Q+EFLhMmTLB27drZl19+aWeccUaK11TQpjPF6oJXr1RG9HwFVvrRlJEjR9oll1zi0hX0w6ptrR9d1ad+AFUWT7jlUgNGvWd6j+XLl9urr77qGoCjR492j2vb/vvf/3bPUyAmaghlJpQ6lL///tuVV/umfmj1Yx6MAnDtrwqE1IhUmfVDrR9iNSx05jSQPrPOIKrMahyooZyRhQsXumBI+5DWVz0rOPv222/DHo+gs7d9+/Z139/A5+p7+Ouvv7r9UfTdue6666x9+/b++v7ll1/cfhv4efS4hJqiqO+zgospU6ZY//793TIFz/quafvqzHKw7fXggw+6/UHb+88//3QNbh0X9F1S4zUr20xl17FI9aqGoPanwP01IwqwFJCqYaVGuxpx+h7NmDHDLr744rCPg8HoGKCy6cSG6lz7jOpHx2WdCJAdO3bYWWed5Q9g1LieOXOmO56oITlw4EC3nv6v+tU21fbXse+1115zgaT2o6ZNm2a6j+q9FYzouKXlCrbUyNTxS9sonO9uuPSbo+NK6u+33lvlULmC0T6ixzMan5EZ1av2MQXpSonT75NOiHmNdNFjSr3VcUT7tbeevlfesV77iepbwevtt9/uggs13rXP6DihACyjY1qo29qj3w+dnLr77rvde6bX4+T9RunkmraX3kf7t8qd+vul3it9BrURFHzq+4wY4ANy2E8//aTWpW/EiBHu78TERF/x4sV9kyZNcn9XqlTJ98ILL7j/792715c/f35f3759U7yGnl+oUCHf2rVr/cu+//57t/y5557zLzt48GCa9//mm2/cem+88YZ/2X//+1+3bP78+SF9hltvvdWtn9r69evd8vLly/v++ecf//IPP/zQLf/444/9y9q3b+9r0qSJ7/Dhw/5lycnJvrPPPttXt27dDN//yy+/dK/31ltvpVg+a9asNMtXrVrl6urf//63b9euXb4TTzzR16JFC1fvHn1uPU+Pqc497733nlv+zDPP+MunsnXq1Mn9P7Cea9Wq5evQoYN/2bBhw9xzr7vuujTl9x4LpL8LFy7s6tDz8ssvu+WVK1dOUa6hQ4e65d66WSlX7969U7z/5Zdf7rZbIO2XPXv29IUi1DqU8847zy0bN25cmtfRY7p5xo4d69adPHmyf9nRo0d9rVq18pUoUcL/Xt6+V6pUKd/OnTtDKrPW123p0qX+ZX/88YevSJEirj48qoOaNWtmuh13797tnnvvvfemWO+OO+5wdbl//37394ABA1w5jx07lmH59J7B3je9utf3eMaMGb6EhATfxo0b3WODBw/2nXzyye7/qtdGjRr5n7dhwwZ3fHn00UdTvJ6+MwUKFEixPCvbTNvec+DAAV+dOnXSHGeC1W3q45a2d+PGjX3t2rXL0nEwmC5duqSoi2D69Onjq1Kliu+vv/5Ksbxbt26+0qVL+8up7XjkyJEU6+hYo2N54Pcso3303HPP9ZUsWdLtf4ECv8/hfHfTo+frc/3555+uDNr3L7zwQrd8zJgxKfZlLVM9ZeSyyy5z63nfwwkTJri/v/vuuwyf563XvHlzt309jz/+uFuu3wxRGbWNO3bs6EtKSvKv9/zzz7v1Xn/9dff3ihUr/N+BjKR3TAt1W3vfNX2nUu+n3mPe/q3PVbFiRbfvHjp0yL+evqNa76GHHvIvU5m0bMiQIRmWH3kPqVDIceoJ0Jk0b+yEziiqe9ab9Un3XvqKxl4kJSUFHV+hszeBZ5jUw6GuZPUKeALPECcmJrozjkqn0VkSnenKKddee61LLfB4Z7y8sv3zzz/u7KXOvOnsnnLCdVP5dMZGXdqpu4oD6Yy1zkBpEKL3XN10tl7dyYEpCDpzrLNdOqun19Z66goPNsZB3fA66+nR2cwqVaq4weqis+4qm85oqqze+2r76QytuvUD09jk5ptvDrne9BqBaSFnnnmmu9cZ2MByecu9+syOcmkb6bk6OxeJzOrQozO1OpOXGT1PZx91NtijM5Q627x//37X4xBIdZU6bSMjSjfTfhM4A5vO+qtnS9+9cGif9HoMvN4pvYZ6RLy0INH3T9tGPRcZUU9FqL0Vno4dO7oz4DrLrzLoPrDuAulMv/YLfQ8Dv0eqb/WypU7lCWebaZtr23t09tU7S5yZwOOWzqTrjLr2z2DHrFCOg8FoGyhdK1iKpqju1HOqQc36f2D96DiiMnnlUYqZd8Za9anjm84+K7UqWJlT76PqJdJ3VL2Q2v8CBUszjfS7q94Uvb96OVRGpYqqV3TQoEH+dXRclsDvcjDe41k9bmifCBxnox4JHZu948XcuXNdb4R6DNQL61HPkLaz11ui757oextuWlg429rTs2fPFPtpMJoZS7226k0LHHuhXjdlIQTr6fF6GhE7SIVCjtMPhYIHr7GnIEIHeDX4RY8ppUe8ACNYYJH6B0jUmA/MT1VKgrpf1fWuhnpgKo6Xd5oTUpfNCzK8simtRWVRCoZuweiArDSpYNSIVvlVb+k9N9DgwYNdA0tpCUpJ0hiGYNSYSr2ttF28xp3e1/tRSY/KFRhUKWUhq/Xm/VhWr1496HKvPrNSroy2kX6wsyqzOvRo24YyYFVpEXrNwEZFYKqeHg8UTn0HK683uFWNEzX41MgON7BSIKH0M6WKqGGk9AelSXnU0FDan9KKVA8KBtS4VwpWpNRIU9qQ0hGV8qGxEl5qT2rab/Q9DFYH3mtldZtpm6duFGvMQCiUvvLII4+4gDlwXEawRnYox8FglKqnbaM6Ulm1DVRP3lg2bXul0Gh8gW6ZHWd0suLJJ5904yR0Eiej/TH1Mi8ICjX1LrPvrgIbNcY9agB7xwxR8Kt0H62jwErHRO3vgd8xL2DwAoz0hBqApCf1vqcTQwpKveOF9/1Ove9oP9SAce9x1akCo6eeesqNmVKwpXQ6pewFfvZgwt3W3vtlJr2yiwKL1BO1KKBSyjJiC4EFjgsFCspfVd6/N77Co/+rIaxAQAce5ScHm3EjvYGYgcGDck0VVOhsj87M6gCrH2flSKc+g52dMiub997KTw0c0BrIC7SC0fMVVOgHJJjUZ6z1w+01vr3B0FnhlVuDGlPnTQf+MAbK7KxWKPUWan2GU65Q9p+cFE69RPt107uYYbAeDe3PGnugAdEKLHSv4ERn1j3ad9Vo1tlV5XHrpu+pgpLUg9KzQg3kcePGudx0jSNKL5DWfqPPpvcPtj9Esi9nlQIyNQhVdxrvokamAhzVT7CpdLO6HyswVa6+ghgN9tYZa72fBmGrh9P7Tqlhml7A7o2D0zbWmB/1SunYre2rcumkTrDB95HWY2afWePZAnvyVH7l+nvUePX2Rw0U1oBiBRqavEPPFf1WqO6DjcsJpMcVcEZyMiK7KLDTdvjwww/d9Mfq1dQ20IQoGTXYw9nWOfldUI9g6hMoyPsILHDcr2ehwCJwYJjSMnSA0eA1DSTUgT+rNJOKDpQ64Ho0aDD1bBrhXgU60qtGe4GSGgyBDa5QKfVBZxt1djGzA7x+NPRjox8+1bPOzilFw/sBDeQFH4E/1Opd8X5UvJQLvVZWyp1TcqpcWdnOmdVhuGrWrOkaL9qOgT+6OjPsPR6J1OUVDbRW6o4XoOqMcOrvTLDeEq/Rp4a9GnIaTKsZi5S2kboxqDOuSr3QTZ9NvRiaZUc9eBkF1aEeX3RWW8eQjAb0ar/R9tHZV/XSZBdtEw1g12sH7kOhXK9CDXyljSjo0nHQo8Aiuyk1TWmbuunsvY4JGig9dOhQt+11Fl7BY2bfKR1ndUxTalng59UkCuEcD1Vn2UHH+8AeG52cymwijKefftpNLqBJM7zPoIkjNMhdv1PBes0VBKpnQc+P5PungMaj9MZt27b5f/e877f2ncATbNpemrAj9bbR7Gq66bNo4hL9RijIVg9Yese0cLZ1OALLrkk0AmlZpMcu5A2EijgulNeqH0+dcVfPRGCPhX5MdeEhTemqPOxIrl+hxkzqM3ea8SX1mVYv9ztY4ymYcNdPTWf0NNOPGlL6EQnWNZ0RpY3oM2hmjtSU2xxYLnWN6wdGXdxaX3WtPFblz6amC0YFdv2rwaDyKWXFC/rUGNOMHfoBDLfcOSWnyqXtHO42zqwOw6UGxvbt2116UeA21n6sM+qa4jUSGscUmD+t1CGd8VRqjBcMqG6VShZ49lafKb2LgyntSQ07Nbi0PQKvGSDKhw+kgMkLvAJTf8KdbtajxpNmR1LDNjAFKzU1pPUZdYY+9XFCf6cuZzjbTDPHBU4R601bmhmVR+UPPEap8ZrdFzVL/dkU6KlnR59bqUwqh8ZCKNAJ1uAP/E55+0lgHeqkkPatUKhhqx4aTW2q2awi7UHU8UANZO+WXo9VYArOXXfd5WYm077vUe+LTtxoP05dX0q30lgPBeBaL6u0TwSmjmm2J32/veOFyq9to/05sC40TkTfSW+WMI3x0PMCKcDQdyvwOxXsmBbOtg73d16/dQpsAsugHkLVdeoZzhCb6LHAcaEDpaaf0xkfBRKBg0dFjV+vlyGSwEJnnDTFnrq19eOiHzqd6U89DaPSZ3Rw1dlNHaxVJp1hSW8Mg1dedTUr9UPPVXpVOBQ46bPp4K8zujobpVx0lVGDKr1rTASjxqR+7NTNrZQSNQLV+6GzXxrYren81Cuhg7fOAKvHQmeGRWeS9Xm9PPdAGvSqMmmAqsqiqVJ19ljlE/1IaRC4fvQ0x7nWUxqAgkMNdFWPgVLcjrecKpe2s/YXBWfexR29gePpyawOszK4UwGotqGuoaDB7WqwqqdPr53V3G6P8tq1DwdONyuBV5bXvq2cfJ3N1XpqJKsBpLP8wQbnnn766e51tS8q5UYnCgJpyks1zPQdU4qGej4UKGm/DJzmOdzpZgMpj163jChg0plcnaHXeyiVR/WpM8EKmlT3SlcMl7a1xokptUvbTCk1Og6FMn2mGlva3zTeRD0/ym3XsUL7UGZpOeHQMUMpajqjrdQ1HStUZr2/t0+NGjXKfX+0z+sz6Riq7aZtru+F/u8dZ9Vbof1Dz1f9qTGp9YMF+sGo4azvjfYV1bu+a9omGuCrY1xO0/dLaWD6DdB+4I1/UGqepv3VcVpTr3rlUsNeJ2c0UUEo01CnRz0P2s91skhn8fX9Uz0oHc4LurR/6vuofULLvfX0G+oF7ZoMROlcGl+k76WCDO1zXtCQ2TEt1G0dDv0mqT51LNRvliZR8Kab1XFM13ZBHIj2tFSIH96UoZpeNbVp06a5xzT9YLApKfWYpnxNTdM2Bk6lpykPe/Xq5TvhhBPc1JyajnT16tVp1pPx48e7KfQ0/WRmU8+qTLfffruvQoUKbmpL76vjTacYOG1hYJk1XWKgdevW+W644QY3nWrBggXdVKWXXHKJb+rUqb5QvPLKK266wqJFi7q60vS199xzj2/r1q2ujC1btvRVq1bNTZ0YSFOfqjzvvvtuimkCp0yZ4raLpgjUa1588cVppn/0pja84oor3BSPmiJW9XnNNdf45s2bl2ZqSE3rGOp0s6m3aXr1GTi9aHaVy5v+MXC6W+0rmgZTdaHHMpp6Npw6TD3taUZTl8qOHTv8+7GmntR2VnlDqauMeHWuqWw1Xa/q7PTTTw+673/22Wdu2ki9f7169dxzgm3H1NNmPvbYY2ke0/6t6TNVR3q9GjVq+G666Sbftm3bIp5uNiPp1fv777/va9OmjZuKU7f69eu7elmzZk2mz01vm2mbayrSYsWKue2mKXa96aAzm272tdde828PlUXbOtTvjAQ7vqWmqZy1b3vfldq1a7upeffs2ZNm39N7VK9e3R2jdKzSVNk69gROCavtrPf19iFNKZr6s2W2j/74449u6tgyZcq4aYu1nz344INZ+u6mJ706k4cffjjosf+HH35w02ZrOlavDvS3piVOLdzpZhcuXOjr16+fr2zZsu43qnv37r6///47zfqaXlb7gt5f0/j279/f/b55fv/9dzcNr7aj6q5cuXK+888/3zd37twUr5PRMS2UbZ3Rdy31dLMe/c5on9C+oXLpM27evDnFOiqHvnuIPQn6J9rBDYDjS7noyvPVGebAKTJBHWaVzkrqjKTO7gabuQiIZ96F4zQrlVKGgFjFGAsAQER0fkqpIkp/IKgAgPjFGAsAQJZosoWPPvrI5WprWuPAgbAAgPhDYAEAyBLNHqMBx7qq83333ecfgAoAiE+MsQAAAAAQMcZYAAAAAIgYgQUAAACA+B5jkZyc7K52qov7BLtsPQAAAIDIZv7bt2+fu8iiLlAbs4GFgorq1atHuxgAAABATNu0aZNVq1YtdgML9VR4H7RUqVLRLg4QNQM/fdh2HdpjZYuWtrGdH3bL6tevb9u2bbMqVarY6tWr2ToAACBse/fudSfyvXZ3zAYWXvqTggoCC8SzRzrfY0m+ZMufkM//XfC6K3XP9wMAAEQilGEHeTqwAPA/VUtVpioAAEBUMSsUAAAAgIgRWAAAAACIGKlQQAz46o9v7cixo1a4QCFrU/OMaBcHAHK1pKQkS0xMjHYxgFyhYMGClj9//mx5LQILIAZM/n66/XNot5UrWobAAgAymI9/+/bttnv3buoICFCmTBmrXLlyxNeFI7AAAABxwQsqKlasaMWKFePiuoh7Pp/PDh48aDt37nR1oSnqI0FgAQAA4iL9yQsqypcvH+3iALlG0aJF3b2CC30/IkmLYvA2AACIed6YCvVUAEjJ+15EOvaIwAIAAMSNSHPIgViUkE3fCwILAAAAABEjsAAAAMhj2rZtawMHDgzrOR988IHVqVPH5dCH+9zMznbrtY+HBQsWuPfLKzN7nXTSSTZ27FiLFwzeBgAAcW39hEnH7b1q9epp0XLTTTdZr1697I477rCSJUvmyHts2LDBatWqZStWrLCmTZtm++ufffbZtm3bNitdurTlBd99950VL17c4gWBBQAAQIzbv3+/m/WnU6dOVrVqVcuN055q5q4CBTJumhYqVMhdbyG3O3r0qCtrhQoVsuV18gpSoYAYUKZIKXdxPN0DAGLLgQMH7IYbbrASJUq46ww8+eSTadY5cuSI3X333XbiiSe6M+RnnnmmSxsS3Xs9FO3atXOpRFr2999/23XXXeeeo1mBmjRpYlOmTMk0lUc9EQ8//HDQsqq3Qk4//XT3PkrZyiilaebMmda8eXMrXLiwffXVV+5zqEdF054WKVLE2rRp4876ZyUVauLEie7CbzNmzLB69eq5z3jVVVe56zZMmjTJfbayZcu691NQkx59Vn3ml19+2apXr+5e55prrrE9e/b417nxxhuta9eu9uijj7rATe8XrP42btxoXbp0cduyVKlS7nV27NiR5r1effVVV5eqg7yEHgsgBozqODSi54+ZvDTTdQb3aBHRewAAsmbw4MG2cOFC+/DDD12D+7777rPly5enSDW67bbb7Oeff7Z33nnHNWynT59uF154oa1atcqlD61Zs8Y1dt9//333d7ly5ezPP/90jfp7773XNXI/+eQTu/7666127dp2xhlnZKms3377rXvu3LlzrVGjRpmebR8yZIg98cQTdvLJJ7tG/j333OPKqIZ/zZo17fHHH3e9LGvXrnVlDpeCiGeffdbVy759++yKK66wyy+/3AUcn376qf3+++925ZVXWuvWre3aa69N93X0/u+99559/PHHtnfvXuvTp4/dcsst9tZbb/nXmTdvnqvHOXPmBH2N5ORkf1Ch7Xns2DG79dZb3ft6QaD3XqqDadOmRXRNiWggsAAAAMjFKUyvvfaaTZ482dq3b++WqdFdrVq1FGfBJ0yY4O69NCf1XsyaNcstf+yxx1xAImqce6lE6qnQep7bb7/dZs+e7RrQWQ0svNQfXYQwlJSl//znP9ahQwd/z8xLL73kehouuugit2z8+PGuoa46UIAVLl2XQa+pYEnUY/Hmm2+6XgI18Bs2bGjnn3++zZ8/P8PA4vDhw/bGG2+4OpPnnnvOLr74Ytd75H3O4sWLu56G9IIpBR4K9NavX+96PkSvqQBMvTItW7b0pz9peaRpVNFAYAEAAJBLrVu3zjU0ldrkUXDgpdqIGqtK5TnllFNSPFdpRRldZVzPUdChQGLLli3uffSc43kRwRYtWqT4rAoE1HvgKViwoAtyfvnllyy9vj6LF1RIpUqVXHqSgorAZRp/kpEaNWr4gwpp1aqV64FQT5AXWDRp0iTDHhp9BgUUXlAhCmzUe6LHvMBCPTV5MagQAgsAAIA83quhlJlly5alSZ0JbECnNmbMGHvmmWfcGAA1inXGXdPQKsDw5MuXzw2sDhTp1ZkD5fSMSQpMAml8RrBlChJyy2cpnodnkWLwNhADXvnuLXtq0Xh3DwCIHTrbrobwkiVL/Mt27dplv/76q/9vDZRW74POuus6FYG3jNKRFi1a5HL+e/ToYaeddpob5xD4uqIz55re1aPxBUrlSY93xj6jwdAZfVY9X+UKDGKUJqQz+9GkNLOtW7f6/168eLELugJ7jjLToEED27Rpk7t5NC5GA9Gj/fmyCz0WQAxYvu1H++fQbjczFAAgdqjHQQOFNb5AaU0aK3H//fe7Rq1HKVDdu3d3M0cp51+BhgZmK6f/1FNPdWMBgqlbt65NnTrVvv76azdw+qmnnnJjDwIbuZpFSmMeLr30Upey89BDD2U4oFjlK1q0qBvfoXEgmtUo1GtO6Ex9//793WdVupfSjzR4WwOwVQfRpM/Rs2dPN9BcwZVmktKMTuFMfXvBBRe4niFtK/USafC2BoCfd955KVLC8jJ6LAAAAHIxpSydc845rnGvxqmmYNVsToE0SFuBxV133eXOomvqU53pV+M8PQ888IA1a9bMzbqkaWHVSNbzAg0dOtQ1fC+55BIXoOjxwDELqek6FJqFSVOzaiC5ekTCMWrUKDdLk2anUtk0Q5IGlCvwiSb1/mhGqc6dO1vHjh1dwPbiiy+G9RoJCQluZi99lnPPPddtS/USvfvuuxYrEnypE+fyEEWMioI1j7Cm9wLi1c0fDfX3WIy7bKRbpjNFGoynwWabN2/O8PlMNwsg1mlWH6Xw5MVrAyC6dG2JDz74wFauXBmX34+9YbS3o9pjofy7Bx980H0IdZspAh4xYkSaQUIAAAAAcreoBhajR492cws///zzbpot/a1cOs0NDAAAAKRH17rQGJRgN02jizgbvK3BQsq98wYVaV5hXUpeV20EAAAA0qOL0R06dCjoY1m5SndGqVC6IZcHFrqk/CuvvOKmNtOMBt9//7199dVXblYCAAAAID2BF6xD7hDVwGLIkCFuQEj9+vXd1GUac/Hoo4+6abiC0dUgdfPouQAAAADifIyFLiH/1ltv2dtvv23Lly+3SZMmufmBdR/MyJEj3ah07xZ4SXQAAAAAcdpjoQugqNeiW7du7m9dNOSPP/5wAYQuQpKa5lIeNGhQih4LggvArHWNFnbg6EErXqgY1QEAAOIvsNCVFAOvHClKiUpOTg66fuHChd0NQErXN72SKgEAAPEbWOgKkhpToatCNmrUyFasWOEGbvfu3TuaxQIAAACQl8ZY6HoVV111ld1yyy3WoEEDu/vuu+2mm25yF8kDAABA5DRVatOmTf1/33jjjda1a1f/323btrWBAwfmuc+RXRYsWGAJCQm2e/duywtOOukkGzt2rOVGUe2xKFmypKuY3Fo5AAAg9o2ZvPS4vdfgHi0s2p555hnz+XzHLRj44IMPbOXKlZZb6fIH27ZtcxMD5QXfffedFS9e3HKjqAYWALLHwE8ftl2H9ljZoqVtbGcu4gMASF92NKCPHj1qhQoVytXVrOBJlzIoUCDj5q4+R+XKlS23O/r/67xChQrZ8joxlwoFIHscPnbEDh077O4BALFFk9o8/vjjVqdOHTeJjcamaoyq595773UXGi5WrJidfPLJ9uCDD1piYmK6r5c6FUqOHTtmt912mws6TjjhBPcagb0aSr9RqvoNN9xgpUqVsn79+mX63hMnTrThw4e7CyAr1Ug3LROlHf373/92jWS9Xrt27dx6gUaNGmWVKlVyGS59+vSxw4cPh5TSNHPmTGvevLmrK114WddAu+OOO6xixYpWpEgRa9OmjTvrn5VUKJW/TJkyNmPGDKtXr5773Err14REkyZNcvVUtmxZ934KajJL63r55ZfdDKd6nWuuucb27NmTZjtpW1etWtW9n7ctArN9Nm7caF26dLESJUq4utTr7NixI8176UrltWrVcnWQU+ixAAAAyMU03f748ePt6aefdo1ipe2sXr3a/7ga3mrwqvG5atUq69u3r1t2zz33hPweahSr8f7tt9/a0qVLXeCgAEav5dG1xh566CEbNmxYSO997bXX2o8//mizZs2yuXPnpugtufrqq61o0aIuCNAyNbDbt29vv/76q5UrV85d60wN4hdeeMF95jfffNOeffZZF7xkRpcyUFm1rhr5Ksv777/vPmPNmjVdkNapUydbu3ate69wKYhQWd555x3bt2+fXXHFFXb55Ze7gOPTTz+133//3a688kpr3bq1q4P06P31OT/++GN3CQXVv8Yd6xpvnnnz5rlgYc6cOekGnV5QsXDhQhcg3nrrre59FTAFvpfqYNq0aW4G1pxCYAEAAJBLqeGqMRHPP/+8/xpftWvXdo1tzwMPPOD/v85mazIcNXrDCSx01lyBi87c68y4ggT9HRhYqFfhrrvuSvG8jN5bgYMavEpFCkw1Ui+CApidO3f6LyOgQEBjMaZOneqCGp2RV0NbN3nkkUdccJJZr4X85z//sQ4dOrj/HzhwwF566SUX/Fx00UVumYI0NdRfe+01d021cKlHRq+p7SDqsVDgs2PHDvd5GzZsaOeff77Nnz8/w8BCn+WNN96wE0880T+p0cUXX2xPPvmkv740lkI9DemlLinw0LZav369/9puek3NtqpemZYtW/rTn7Q80jSqzJAKBQAAkEv98ssvLpVHZ/PT8+6777qz42qMqmGrxr7SY8Jx1llnuaDC06pVK/vtt99SpPO0aNEiW95bKU/79++38uXLu+d4NzWO161b5//cZ555ZornqUyhCCynXk+BgMroKViwoJ1xxhnuPbJCaUteUCFK11JQVaJEiRTLFDhlRD1CXlDhfT71QKxZs8a/TBePzmg8hD6DAorAC0YrsFHvSeDnU09NTgcVQo8FAABALqWz/hn55ptvrHv37m4sg9J7lFakHgOd9c5uqWciyup7K6ioUqVKilQdjxrE2V3O7KbAJJACsmDLktO54HM0PsvxmkWKHgsAAIBcqm7dui64UMpLMF9//bU7G33//fe7M/Va/48//gj7fZYsWZLi78WLF7vXyigfP5T31tn21IOYmzVrZtu3b3cpUhqQHnjTwHHR9c2ClSlc6llQGRYtWuRfph4MpQnpzH40bdy40bZu3Zri8+XLl88/SDsUqqdNmza5m+fnn392A9Gj8fnosQAAAMilNIOPZl7SmAU1kJXS8+eff9pPP/3kxh+oMa8GqnoKlE//ySef2PTp08N+H73GoEGD3IWKly9f7vL9M+t5COW9lSKkFCddx6JatWpuYPcFF1zg0n4045EGUmtWKTWw9XwNglaQMmDAADcrkv6vz6wBzfrMoQzeTn2mvn///m4shQZqK/1I76kB2N74jWhu2549e7rxJRq8rZmkNKNTOFPfqi6VLqWeI41L0eBtDQA/77zzgqau5TR6LAAAAHIxTeGqQdOakUlnqDUg2Mvfv+yyy+zOO+90U8VqSlH1Imj9cGka2UOHDrmxB5pVSA17b0rZ9ITy3pod6cILL3SDmZXjP2XKFJcmpNmTzj33XOvVq5cLLLp16+Z6OzQ2QfQZ9VoKqDR1rB5TgJAVmrZW5bj++utdb4lmSJo9e7abMSqa6tSp42aU6ty5s3Xs2NFOPfVUe/HFF8N6DdXlhx9+6D6L6lOBhoIvjX2JhgTf8br0Yg5QdKd8Ps35q6m4gHh180dD7Z9Du61c0TI27rKRbpnODG3ZssUNDNu8eXPEV53NDVeLBYCs0gw8OnOe0/P4A3nxiuQZfT/CaW+TCgXEgL4t/mVHk45aofy5+yqoAAAgdpEKBcSA5lWbWKvqzd09AADIGl3rInAK3MDbY489RrVmgh4LAAAAwMxdjE5jTYLJylW6M0qF0i3WEFgAAAAAZikuWIfwEVgAMeD3f/6wY8lJViBffju5XM1oFwcAAMQhAgsgBjz+1bg0s0IBAAAcTwzeBgAAABAxAgsAAAAAESOwAAAAABAxAgsAAAAAESOwAAAAyKNuvPFG69q1a46/j6650LRpU4uWhIQE++CDDyyaFixY4Mqxe/fusJ/72muvWceOHS2ndevWzZ588kmLFgILAACAPOqZZ56xiRMnWm50PIOBtm3b2sCBA3Pl6x0+fNgefPBBGzZsmH/Z+PHj7ZxzzrGyZcu62wUXXGDffvttSIFN6tv27dv96zzwwAP26KOP2p49eywaCCwAAADyqNKlS1uZMmWiXYw8wefz2bFjx477+06dOtVKlSplrVu3ThEkXHfddTZ//nz75ptvrHr16q5HY8uWLZm+3po1a2zbtm3+W8WKFf2PNW7c2GrXrm2TJ0+2aCCwAAAAyMXUMG3SpIkVLVrUypcv785uHzhwIGgqlM6033777e5su86EV6pUyZ0d1/q9evWykiVLWp06dWzmzJn+56jHI3Vwop4GnQ1Pz3fffWcdOnSwE044wQU35513ni1fvtz/+EknneTuL7/8cvc63t/y4YcfWrNmzaxIkSJ28skn2/Dhw1M0+H/77Tc799xz3eMNGza0OXPmZFg/qoOFCxe63hvvLP6GDRv8Z/j1WZs3b26FCxe2r776Kmj6mOpLdZfR63mWLVtmLVq0sGLFitnZZ5/tGvoZeeedd+zSSy9Nseytt96yW265xaWX1a9f31599VVLTk62efPmWWYUSFSuXNl/y5cvZXNe76X3jAYukAcAAOKWGoiBqSTHgxqDS5cuDWldnZHWme3HH3/cNdL37dtnX375pTv7np5JkybZPffc41Jr3n33Xevfv79Nnz7dPf++++6zp59+2q6//nrbuHGjaxxnhcrRs2dPe+6551xZlNffuXNnFxQoeFHgoQbwhAkT7MILL7T8+fO756nsN9xwgz377LMuFWjdunXWr18/95hShdS4vuKKK1xAtGTJEpfSk1lKkgKAX3/91Z2t/89//uOWVahQwR8MDBkyxJ544gkXxCjYykxmr3f//fe7z6tlN998s/Xu3dsWLVqU7uspmFF9Z+TgwYOWmJho5cqVy7R8CkaOHDniyqexL4E9IXLGGWe4dCito2DqeCKwAGLA0xcNM5/5LMHSP7sEAEhLQUUo6SfRosBCZ/PV2K5Zs6Zbpt6LjJx22mku116GDh1qo0aNcj0Lffv2dcseeughe+mll+yHH36ws846K0vlateuXYq/X3nlFdfroTP9l1xyiWt0i5YpkPKod0INfQUlosb+iBEjXCCkwGLu3Lm2evVqmz17tlWtWtWt89hjj9lFF12UblnUY1KoUCEXJAW+l0fBgXpXQpXZ66nRrh4a0We5+OKL3TiKIkWKpFlXA70VHHmfJT333nuvW0e9UempUqWKjRs3zgXDChrUy6FeFgVg6gHy6HWOHj3q9m1vnzleCCyAGFC0YNqDGQAgc8EajrnpPRUktG/f3gUTnTp1cnn4V111VYZn3k899VT//9VToPSpwGBEvQGyc+fOLH+GHTt2uOBF6UZ6naSkJHfWXb0gGfn+++/d2X01zj16rhrmev4vv/zixhsENsRbtWplkVBDPDsF1q8a+6I6qFGjRpp1Dx065O6DBR0eBX5KXVJdZrRevXr13M2jNCz1+KgH6s033/QvV8qcqD6PNwILAAAQt0JNSYoWBQYaY/D111/bZ5995lKPlIqjs9S1atUK+pyCBQum+FtjBAKXeWMnlHYkytFPnVqltJyMqMfh77//dmlDOiuulBsFADpTnpH9+/e7Xgv1wKSWUaM6EsWLF0/xd1Y+b6CM6jI1BXVaZ9euXRaMUrQUWKinJjBgCZXSnpRqFeiff/5x916v0fHE4G0AAIBcTA1T5dGrQb5ixQqXpqMxE9lFDVCNmfAGhMvKlSszfI56He644w43rqJRo0YusPjrr7/SNMDVGxFIKTsa7KwB5KlvavA3aNDANm3a5FLAPIsXL870M6hOUr9XRp838PWDfd5wXi+zcjVs2NB+/vnnNI9p3IzSwGbNmpXlXhWV2+s18fz4449WrVo1l/52vNFjAcSAGWvm2sHEw1asYBG7pF76+ZkAgLxFPROaKUgpUBoMrb///PNP1wDPLmeeeaYbT6CB3QoW9B6ZXRujbt26Lv1GDeK9e/fa4MGD/Sk4Hs0EpbIrKFLgofQtje/QGAylDSmlS8GE0qPUGH7kkUfcGINTTjnF9YiMGTPGvbZ6aDKj91K5NcC6RIkSGQ6C1vgQvfYbb7zhelk0Nave//TTT8/S62WmU6dOrlchcBD66NGjXV28/fbb7r28CQT0Xrp542M0/kfllLFjx7peKgVySh3TGIvPP//c9WQF0gD543ExvlzXY6GKDHahj1tvvTWaxQLynBlr5tnUnz5x9wCA2KHrH3zxxReuZ0ANbo1r0IxEGQ1mDpcazWpcf/rpp24sxpQpU9xsQ5ldSVrpPeqB0IxHCkgCr6cgKqfSuDRmwmu0q5E9Y8YM1xhu2bKlGzyuMQLeIGMFGuqN0dgEpfn8+9//TjEeIz133323SxtT74B6JDIa66Ey6IJ1GjCuMqi3RjNVZfX1MtOnTx9Xt4EXrdPgeaWNKbhSj4N3U2qUR70qge+r9e+66y63jTR4XAGZUqg0BsejgENTBXsD9Y+3BF9G85XlMEXcgd1MihY1al8XC/HmEs6IoliN3NeG0hcPiFc3fzTU/jm028oVLWPjLhvplqkbVGc6TjzxRNu8eXOGzx8zOfMc48E9snfwGwAcT2pwrV+/3p3xzalcfiA9V199tQvC1AuRkxSwKDBL3YsRyfcjnPZ2VHssFAEGXuBDEayuFuhN4QUAAADkdWPGjPGnOOUkjWvRAP9oyTVjLNS9o264QYMGZXilRwAAACAvOemkk9wV0XOaUseiKdcEFsoH00VEdBn19OhiILoFds0AAAAAiL5cM92sBgFpIFJGVyYcOXKky/HybhoMBAAAACD6ckVg8ccff7hR7Zl132jAiwaOeDfNcwwAAAAg+nJFKtSECRPcFGUXX3xxhutpDmTdAAAAAOQuUe+x0CXQFVjoQigFCuSKOAcAAABAmKLeklcKlC7+0bt372gXBcizapWtbuWLlbVShXN+KjsAAIBcGVjokuNRvEYfEBPuPeeWaBcBAJALaMp+XSCta9euFssmTpxoAwcOdDOKIveIeioUgJynK2tndAMA5F6ail8Bg266AFqlSpWsQ4cO9vrrr7uU8kDbtm1zs2x6vOfpVrx4catbt657vWXLllledu2119qvv/5quUHbtm1dkBPKet620NWtTznlFDfjabAT7JMmTbKWLVtasWLFrGTJku7i0bqQdKAFCxa418pNwRWBBQAAQC534YUXuqBhw4YNNnPmTDv//PNtwIABdskll9ixY8f861WuXDnNRDcay6rn/vTTT/bCCy/Y/v377cwzz7Q33njD8qLExEQrWrSom/gnr+nbt6/bFmvWrHGznT700EM2bty4FOvcfffddtNNN7ng6YcffrBvv/3W2rRpY126dLHnn3/ecjMCCwAhodcDAKJHwYKChhNPPNGaNWtm9913n3344YcuyFBakEdnsHXR4UBlypRxz9XVn5WCPnXqVOvevbvddttttmvXrnTfU2fCdSmAChUqWKlSpaxdu3b2/fffu8f+/PNP95qPPfaYf/2vv/7aChUqZPPmzXN/P/zww9a0aVN7+eWX3bXHdPb9mmuucZcMCPTqq69agwYN3Fn8+vXr24svvuh/TIGUPtO7777rztprnbfeest9Zn0uj/de6sWpUaOGlShRwm655RZLSkqyxx9/3JVVgcijjz4a8mcMfN0333zT1Z+uo9atWzfbt2+fe1y9PwsXLrRnnnnG3xuhMqdHdaCy1KxZ03r16mWnnnqqzZkzx//44sWL7cknn7QxY8a4AKNOnTqublRu9YoMGjQoV19ugcACiAGjv3zR7p/7uLsHAMQHNYJPO+00mzZtWtjPvfPOO13jOLBRm9rVV19tO3fudMGLUqcU0LRv397++ecf1xBXI14N76VLl7rXuv76612wonU8a9eutffee88+/vhjmzVrlq1YscI1+D0KEnTWXg3nX375xQUqDz74oEsFCjRkyBDXQ6N1OnXqFLS869atc2XV+0yZMsVdfFmXMti8ebNr/I8ePdoeeOABW7JkSUifMfB1FawpFUk3vdaoUaPcYwooWrVq5e+J2LZtW0gXcFb605dffmmrV692wZhH5VZQpB6L1O666y7XW/P+++9bbhX1wdsAIrd+1yb759BuK1f0/87eAABCM2PNXJux5n9n2TObgS/1ZBk6oaNjcGYuqdfeLql3QbZvEp3hV7pMVp4n6Z1d/+qrr1wKjhrdXmrVE0884RrY6vHo16+fde7c2TWo1fvRokULN4ZDYwYCHT582KVcqadFnnvuOdfY11l5nbkfNmyY+/8VV1zhHq9Vq5b9/PPPrpdDlyLw6Gy9t056NN5EwY7GJDRs2NCliynl6NNPP7V8+fJZvXr1XHAxf/58lwoWymf0Xlc9JHpdUQClXhkFQ+rBUGDg9URkRr0x6qE5evSoCxLUA3PHHXf4H9e4kdq1a6cINjxVq1Z1vSq5ZWxJMAQWAAAgrh1MPOxOzmRG03qntvfI/pCeq/fICTrzrfSbrDxP0nuu0oE0FqN8+fIplh86dMidwfeoId64cWP773//6874px7fobQkL6gQnd1XQ10NfjXU9Vp9+vRxAYpHY0bUYA+kwCUzSlXyGv+iQe758+d3QUXgMgUS4XzG1K9bpUoV/2uEq3v37nb//fe7FDQFVWeffba7BcrLs6USWAAAgLhWrGCRkHp8g10rSMtCea7eIycoNUhn+bPyPEnvuWpwqwGtmYdSCxzboAb41q1bXbCg3o8mTZqEXAa9h4wfP971IARSQBBIvSGZ0YxZgbxZtFIv82bSCvUzZvQa4SpdurQbNyFKEdP/zzrrLLvggv/1ZmmmKPWkqEcjda+F6nnv3r1undyKwAIAAMQ1pShlNU0pmtcR+vzzz23VqlVuvES4xo4d69JqvAZtahprsH37ditQoIA7Yx+MGr89evRwsxcpzUiDoFWewNmadBFkNYiVxuMNTvbSktR7oOW///67O5N/vIXyGUOhAECDxMOlsRQaN6JB2hp7ooBFA8OfffZZlwp2++23p1hfvUMKcq688krLrQgsAAAAcrkjR464RrAasDt27HADlDWeQdPN3nDDDRk+VzMf6bl6DeXnq9GqcQQa+xB4Zj6QAg6lLelCe5pVSWfJFSB88skndvnll7vUJKX0aIYnNYTVSNZYht69e6e43oLGEGishBrFOtuu8QSaGcobjzB8+HC3TGfyNaWuyqjB4EoV0gxIOSmUzxgKBSUaEL5hwwZXD+XKlUuRfpURDdIeMWKEG5B91VVXufIo2Bg8eLAL3FQ2jcWYPHmyGyiugDCUweHRQmABAACQyymQUNqOzq6XLVvWzQalBr0a7Zk1YjWtqdfI13gHXRNBg5Z1xj49OnuuQEHBg57vTS977rnnup4GpQ+pkauB0Or5EE3JqnK99NJL1r9/f7dMqT4adK2B3pppSYFQ4HSy6uXQwGdNr6rGtFKelE4VygXnIpXZZwyVehy0HRo2bOjGZ6xfvz7kHhAFIQoMNbuW6knbUvWqaWhVT5rFSmlh2lYKBi+99FLLzRJ8eXiEiCJfRbiKlr2dGohHN3801D8r1LjL/jcjR7Vq1WzLli3uR2TAqJRzmueEwT1CO7MDANGg2YnU4NOYAjWwkfPUWFZjeOXKlVR3Hv5+hNPe5joWAAAAACJGYAEAAAAgYoyxAGKALrykOdJzajpDAACykgqlG+IHgQUQA3Liaq4AAADhIBUKAAAAQMQILAAAQNzIw5NhArn+e0FgAcSAQ4mH7WDiIXcPAEhLVyyWgwcPUj1AKt73wvueZBVjLIAYcOfM4WmuYwEA+D+6yJiuMr1z5073ty7KpgukAfHeU3Hw4EH3vdD3Q9+TSBBYAACAuKCrKosXXAD4HwUV3vcjEgQWAAAgLqiHokqVKlaxYkVLTEyMdnGAXEHpT5H2VHgILAAAQFxRIyq7GlIA/g+DtwEAAABEjMACAAAAQMQILAAAAABEjMACAAAAQMQILAAAAABEjMACAAAAQMSYbhaIAfe0udmOJSdZgXxMnwgAAKKDwAKIASeXqxntIgAAgDgX9VSoLVu2WI8ePax8+fJWtGhRa9KkiS1dujTaxQIAAACQV3osdu3aZa1bt7bzzz/fZs6caRUqVLDffvvNypYtG81iAQAAAMhLgcXo0aOtevXqNmHCBP+yWrVqRbNIQJ60bOsqO5p01ArlL2TNqzaJdnEAAEAcimoq1EcffWQtWrSwq6++2ipWrGinn366jR8/Pt31jxw5Ynv37k1xA2A2funb9vTXr7p7AACAuAssfv/9d3vppZesbt26Nnv2bOvfv7/dcccdNmnSpKDrjxw50kqXLu2/qbcDAAAAQJwHFsnJydasWTN77LHHXG9Fv379rG/fvjZu3Lig6w8dOtT27Nnjv23atOm4lxkAAABALgssqlSpYg0bNkyxrEGDBrZx48ag6xcuXNhKlSqV4gYAAAAgzgMLzQi1Zs2aFMt+/fVXq1mTOfkBAACAvCSqgcWdd95pixcvdqlQa9eutbffftteeeUVu/XWW6NZLAAAAAB5KbBo2bKlTZ8+3aZMmWKNGze2ESNG2NixY6179+7RLBYAAACAvHQdC7nkkkvcDQAAAEDeFdUeCwAAAABxGlgcOHAgZ0oCIMuKFChsRQsUcfcAAAB5IhWqUqVKds0111jv3r2tTZs2OVMqAGEZ2/lhagwAAOStHovJkyfbP//8Y+3atbNTTjnFRo0aZVu3bs2Z0gEAAACIzcCia9eu9sEHH9iWLVvs5ptvdlPE6roTGoA9bdo0O3bsWM6UFAAAAEDsDd6uUKGCDRo0yH744Qd76qmnbO7cuXbVVVdZ1apV7aGHHrKDBw9mb0kBAAAAxN50szt27LBJkybZxIkT7Y8//nBBRZ8+fWzz5s02evRod+G7zz77LHtLCyCoN1e+bweOHrTihYrZ9U2vpJYAAEDuDyyU7jRhwgSbPXu2NWzY0G655Rbr0aOHlSlTxr/O2WefbQ0aNMjusgJIx6KNS+2fQ7utXNEyBBYAACBvBBa9evWybt262aJFi9yVs4NROtT999+fHeUDAAAAEIuBxbZt26xYsWIZrlO0aFEbNmxYJOUCAAAAEMuBRcmSJV1wUbFixRTL//77b7csKSkpO8sHII8YM3lppusM7tHiuJQFAADkgVmhfD5f0OVHjhyxQoUKZUeZAAAAAMRqj8Wzzz7r7hMSEuzVV1+1EiVK+B9TL8UXX3xh9evXz5lSAgAAAIiNwOLpp5/291iMGzfO8ufP739MPRUnnXSSWw4AAAAg/oQcWKxfv97dn3/++W7K2bJly+ZkuQAAAADE8uDt+fPn50xJAAAAAMR2YDFo0CAbMWKEFS9e3P0/I0899VR2lQ1AiJpVaWz7jx60EoUyngoaAAAgqoHFihUrLDEx0f//9GhgN4Djr1/L7lQ7AADI/YFFYPoTqVAAAAAAIr6ORWp79+61Dz74wFavXh3pSwEAAACIl8Dimmuuseeff979/9ChQ9aiRQu3rEmTJvb+++/nRBkBAAAAxFpgoQvhnXPOOe7/06dPd9e12L17t7uA3iOPPJITZQSQiSGfjbSbPxrq7gEAAPJEYLFnzx4rV66c+/+sWbPsyiuvtGLFitnFF19sv/32W06UEUAmdh/ea/8c2u3uAQAA8kRgUb16dfvmm2/swIEDLrDo2LGjW75r1y4rUqRITpQRAAAAQKxdIG/gwIHWvXt3K1GihNWsWdPatm3rT5HSOAsAAAAA8SfswOKWW26xM844wzZt2mQdOnSwfPn+1+lx8sknM8YCAAAAiFNhBxaimaB0C6QxFgAAAADiU9iBRVJSkk2cONHmzZtnO3futOTk5BSPf/7559lZPgAAAACxGFgMGDDABRbqoWjcuLElJCTkTMkAAAAAxG5g8c4779h7771nnTt3jvjNH374YRs+fHiKZfXq1eMq3gAAAECsBxaFChWyOnXqZFsBGjVqZHPnzv2/AhXI0rAPAAAAAFEUdiv+rrvusmeeecaef/75bEmDUiBRuXLliF8HiFdjJi+1MglNrZQds3z7Cri/Zf/BxBT3AAAAuSqw+Oqrr2z+/Pk2c+ZM19tQsGDBFI9PmzYtrNfT1bqrVq3qLq7XqlUrGzlypNWoUSPcYgFxrbyvdrSLAAAA4lzYgUWZMmXs8ssvz5Y3P/PMM91AcI2r2LZtmxtvcc4559iPP/5oJUuWTLP+kSNH3M2zd+/ebCkHAAAAgOMcWEyYMMGyy0UXXeT//6mnnuoCDV3NW4PD+/Tpk2Z99WakHuwNAAAAIPr+d9nsMB07dswNuH755Zdt3759btnWrVtt//79ERVGvSGnnHKKrV27NujjQ4cOtT179vhvuvo3ALNDttsO2S53DwAAkCd6LP744w+78MILbePGjS4tqUOHDi5tafTo0e7vcePGZbkwCkzWrVtn119/fdDHCxcu7G4AUlpTYKYlJhy0gr5i1vTYdVQPAADI/T0WukBeixYtbNeuXVa0aFH/co270NW4w3H33XfbwoULbcOGDfb111+718ifP79ddx0NIwAAACCmeyy+/PJLFwToehaBTjrpJNuyZUtYr7V582YXRPz9999WoUIFa9OmjS1evNj9HwAAAEAMBxbJycmWlJQUNEgINpNTZlfxBgAAABCHqVAdO3a0sWPH+v/WRfI0NmLYsGHWuXPn7C4fAAAAgFjssXjyySetU6dO1rBhQzt8+LD961//che5O+GEE2zKlCk5U0oAAAAAsRVYVKtWzb7//nuXxvTDDz+43gpdc6J79+4pBnMDAAAAiB8FsvSkAgWsR48e2V8aAAAAALEbWHz00Uchv+Bll10WSXkAAAAAxGpg0bVr1xR/a8C2z+dLs0yCzRgFAAAAILblC3WKWe/22WefWdOmTW3mzJm2e/dud9P/mzVrZrNmzcr5EgNIo+GxLnZaYjd3DwAAkCfGWAwcONDGjRvnLmbn0SxRxYoVs379+tkvv/yS3WUEkIlCVow6AgAAees6FuvWrbMyZcqkWV66dGnbsGFDdpULAAAAQCwHFi1btrRBgwbZjh07/Mv0/8GDB9sZZ5yR3eUDAAAAEIupUK+//rpdfvnlVqNGDatevbpbtmnTJqtbt6598MEHOVFGAJnYmbDakhMSLZ+voFX01ae+AABA7g8s6tSp4y6MN2fOHFu9erVb1qBBA7vgggv8M0MBOL625l9hiQkHraCvmFU8RmABAADyyAXyFEB07NjR3QAAAAAg7DEWAAAAAJAagQUAAACAiBFYAAAAADi+gcWxY8fsjTfeSDHVLAAAAACEFVgUKFDAbr75Zjt8+DA1BwAAACDrqVC6CN7KlSvDfRoAAACAGBb2dLO33HKLu/K2LorXvHlzK168eIrHTz311OwsHwAAAIBYDCy6devm7u+4444U17Xw+XzuPikpKXtLCCBTRXylLb8VsoK+otQWAADIG4HF+vXrc6YkALKsflJnag8AAOStwKJmzZo5UxIAAAAA8RNYyLp162zs2LH2yy+/uL8bNmxoAwYMsNq1a2d3+QAAAADE4qxQs2fPdoHEt99+6wZq67ZkyRJr1KiRzZkzJ2dKCQAAACC2eiyGDBlid955p40aNSrN8nvvvdc6dOiQneUDEIJ1+RfYMTtsBayI1U5qS50BAIDc32Oh9Kc+ffqkWd67d2/7+eefs6tcAMKwL2Gb7c23xd0DAADkicCiQoUKQS+Qp2UVK1bMrnIBAAAAiOVUqL59+1q/fv3s999/t7PPPtstW7RokY0ePdpdOA8AAABA/Ak7sHjwwQetZMmS9uSTT9rQoUPdsqpVq9rDDz+c4qJ5AAAAAOJH2KlQurq2Bm9v3rzZ9uzZ4276v6ab1WNZpcHgev7AgQOz/BoAAAAA8tB1LDzqucgO3333nb388stu6loAAAAAMRpYNGvWzObNm2dly5a1008/PcOeieXLl4dVgP3791v37t1t/Pjx9sgjj4T1XAAAAAB5KLDo0qWLFS5c2P2/a9eu2VqAW2+91S6++GK74IILMg0sjhw54m6evXv3ZmtZAAAAAORgYDFs2DB3n5SUZOeff75LWSpTpoxF6p133nE9HEqFCsXIkSNt+PDhEb8vAAAAgCgO3s6fP7917NjRdu3aFfEbb9q0yQ34fuutt6xIkSIhPUezUHkDxnXTawAwq5BczyolNXL3AAAAeWLwduPGjd01LGrVqhXRGy9btsx27tzpxm941CPyxRdf2PPPP+9SnhTIBFI6lpeSBeD/nJj8f98jAACAPBFYaBzE3XffbSNGjLDmzZtb8eLFUzxeqlSpkF6nffv2tmrVqhTLevXqZfXr17d77703TVABAAAAIPcKO7Do3Lmzu7/ssstSzA7l8/nc3+p1CHWqWvV+BFKQUr58+TTLAQAAAMRYYDF//vycKQkAAACA+AkszjvvvJwpiZktWLAgx14biGUrC0yxxISDVtBXzJoeuy7axQEAAHEorFmhPF9++aX16NHDzj77bNuyZYtb9uabb9pXX32V3eUDAAAAEIuBxfvvv2+dOnWyokWLumtQeBes0/Svjz32WE6UEQAAAECsBRaaFWrcuHE2fvx4K1iwoH9569atXaABAAAAIP6EHVisWbPGzj333DTLS5cubbt3786ucgEAAACI5cCicuXKtnbt2jTLNb7i5JNPzq5yAQAAAIjlwKJv3742YMAAW7JkibtuxdatW+2tt95yF83r379/zpQSAAAAQGxNNztkyBBLTk52V84+ePCgS4sqXLiwCyxuv/32nCklAAAAgNgKLNRLcf/999vgwYNdStT+/futYcOGVqJEiZwpIQAAAIDYS4Xq3bu37du3zwoVKuQCijPOOMMFFQcOHHCPAQAAAIg/YQcWkyZNskOHDqVZrmVvvPFGdpULQBhOTmprpxzr5O4BAABydSrU3r17zefzuZt6LIoUKeJ/LCkpyT799FOrWLFiTpUTQAZK+apQPwAAIG8EFmXKlHHjK3Q75ZRT0jyu5cOHD8/u8gEAAACIpcBi/vz5rreiXbt29v7771u5cuX8j2m8Rc2aNa1q1ao5VU4AAAAAsRBYnHfeee5+/fr1VqNGDddDASB32JuwzXyWZAmWn7QoAACQNwZv//LLL7Zo0SL/3y+88II1bdrU/vWvf9muXbuyu3wAQvB7/gX2a4HZ7h4AACBPBBa6foUGcsuqVats0KBB1rlzZ9eTof8DAAAAiD9hXyBPAYSuXyEaa3HppZfaY489ZsuXL3cBBgAAAID4E3aPhQZqHzx40P1/7ty51rFjR/d/Deb2ejIAAAAAxJeweyzatGnjUp5at25t3377rb377rtu+a+//mrVqlXLiTICAAAAiLUei+eff94KFChgU6dOtZdeeslOPPFEt3zmzJl24YUX5kQZAQAAAMRaj4Wmmp0xY0aa5U8//XR2lQkAAABArAcWkpSUZB988IGbelYaNWpkl112meXPnz+7ywcAAAAgFgOLtWvXutmftmzZYvXq1XPLRo4cadWrV7dPPvnEateunRPlBAAAABBLYyzuuOMOFzxs2rTJTTGr28aNG61WrVruMQAAAADxJ+wei4ULF9rixYvd9LKe8uXL26hRo9xMUQCOv6bHrqPaAQBA3uqxKFy4sO3bty/N8v3797trXAAAAACIP2EHFpdccon169fPlixZYj6fz93Ug3HzzTe7AdwAAAAA4k/YgcWzzz7rxli0atXKihQp4m5KgapTp44988wzOVNKAAAAALE1xqJMmTL24Ycf2m+//eamm01ISLAGDRq4wAJAdGzJt9yS7Kjlt0J2YnIzNgMAAMj9PRaeunXr2qWXXupSo7IaVOjK3aeeeqqVKlXK3dQLoit4AwjPn/nW2I78P7l7AACAPBNYvPbaa9a4cWN/KpT+/+qrr4b9OtWqVXOzSS1btsyWLl1q7dq1sy5duthPP/2UlWIBAAAAyCupUA899JA99dRTdvvtt7seBvnmm2/szjvvdNez+M9//hPya6nHI9Cjjz7qejE0GFxX8wYAAAAQo4GFGv7jx4+36677v3nzNRuUUpoUbIQTWARKSkqy//73v3bgwAF/wAIAAAAgRgOLxMREa9GiRZrlzZs3t2PHjoVdgFWrVrlA4vDhw1aiRAmbPn26NWzYMOi6R44ccTfP3r17w34/AAAAALlgjMX111/vei1Se+WVV6x79+5hF6BevXq2cuVKd12M/v37W8+ePe3nn38Ouu7IkSOtdOnS/lv16tXDfj8AAAAAuaDHwhu8/dlnn9lZZ53l/lZQoPEVN9xwgw0aNMi/nsZiZEZX6/ZmlVKvx3fffeeuh/Hyyy+nWXfo0KEpXl89FgQXAAAAQB4MLH788Udr1ux/8+SvW7fO3Z9wwgnupsc8ur5FViQnJ6dIdwpUuHBhdwMAAACQxwOL+fPnZ9ubqwfioosusho1ati+ffvs7bfftgULFtjs2bOz7T0AAAAA5NJUqOyyc+dOlz61bds2N2ZCM0spqOjQoUM0iwXkOSV9VeyY77AVsCLRLgoAAIhTUQ0sNFYDQORqJ7WlGgEAQN678jYAAAAABCKwAAAAABAxAgsAAAAAeXuMBYDssTr/p5aYcMgK+opa/aTOVCsAADjuCCyAGHA4YY8lJhy0JDsa7aIAAIA4RSoUAAAAgIgRWAAAAACIGKlQAADEmPUTJmX4eK1ePY9bWQDED3osAAAAAESMwAIAAABAxAgsAAAAAESMMRZALjdm8tJoFwEAACBTBBYAAMTZwGsGdwPICQQWQAyomnS6JSckWj5fwWgXBQAAxCkCCyAGVPTVN/NFuxQAACCeMXgbAAAAQMTosQAAIA/JbHwEAEQLgQUQA47aQftfLlSCFbJi0S4OAACIQwQWQAz4ucCHlphw0Ar6ilnTY9dFuzgAACAOMcYCAAAAQMQILAAAAABEjMACAAAAQMQYYwEAwHHCjE4AYhk9FgAAAAAiRmABAAAAIGIEFgAAAAAiRmABAAAAIGIEFgAAAAAixqxQQAyod+wiM/OZWUK0iwIAAOIUgQUQA4pamWgXAQAAxLmoBhYjR460adOm2erVq61o0aJ29tln2+jRo61evXrRLBYAAHEtlOtt1OrV87iUBUDeEdUxFgsXLrRbb73VFi9ebHPmzLHExETr2LGjHThwIJrFAgAAAJCXeixmzZqV4u+JEydaxYoVbdmyZXbuuedGrVxAXvN3wjpLtmOWzwpYeV/taBcHQBzIrFeDHg0g/uSqMRZ79uxx9+XKlYt2UYA8ZVP+by0x4aAV9BWz8scILADEBoIXIG/JNYFFcnKyDRw40Fq3bm2NGzcOus6RI0fczbN3797jWEIAAAAAuf46Fhpr8eOPP9o777yT4WDv0qVL+2/Vq1c/rmUEAAAAkIsDi9tuu81mzJhh8+fPt2rVqqW73tChQ126lHfbtGnTcS0nAAAAgFyYCuXz+ez222+36dOn24IFC6xWrVoZrl+4cGF3A5A3jZm8NMPHB/docdzKAgAAYiiwUPrT22+/bR9++KGVLFnStm/f7pYrzUnXtQAAAHkT18IA4k9UA4uXXnrJ3bdt2zbF8gkTJtiNN94YpVIBAIDcEnwAyDuingoFAAAAIO/LFYO3AQAAAORtueY6FgCyThfGC7wHEB2k9gCIZwQWQAxolNQl2kUAAABxjlQoAAAAABGjxwIAgBCQ5gQAGaPHAgAAAEDE6LEAYsCGfF/ZsYQjVsBX2E5KbhPt4gAAgDhEYAHEgN35NlliwsH/zQqVHO3SAACAeEQqFAAAAICIEVgAAAAAiBiBBQAAAICIEVgAAAAAiBiBBQAAAICIEVgAAAAAiBiBBQAAAICIEVgAAAAAiBgXyANiQPnk2v4rbwMAAEQDgQUQA6onnxHtIgAAgDhHYAEAAPKk9RMmZbpOrV49j0tZABBYAADiRGaNUBqgABAZBm8DAAAAiBipUEAMWFVgqh21g1bIilmTY1dFuzgAACAOEVgAMSDJEi05IdGSfInRLgoQ0/n6AID0kQoFAAAAIGIEFgAAAAAiRmABAAAAIGIEFgAAAAAiRmABAAAAIGIEFgAAAAAiRmABAAAAIG8HFl988YVdeumlVrVqVUtISLAPPvggmsUBAAAAkBcDiwMHDthpp51mL7zwQjSLAeR5JyW1ttrH2rl7AACAuLvy9kUXXeRuACJTxleDKgQAAPEbWITryJEj7ubZu3dvVMsDAAAAIA8O3h45cqSVLl3af6tevXq0iwQAAAAgr/VYDB061AYNGpSix4LgAjA7YH+ZLyHJEnz5rbidkGerZMzkpRk+PrhHi+NWFgCxYf2ESRk+XqtXz+NWFiDW5anAonDhwu4GIKXfCsyxxISDVtBXzJoeu47qAQAAx12eCiwAAACOZ4+G0KsB5IHAYv/+/bZ27Vr/3+vXr7eVK1dauXLlrEYNZrkBAAAA8oqoBhZLly61888/3/+3N36iZ8+eNnHixCiWDAAAAECeCSzatm1rPp8vmkUAAOQBDMAFgNwvT003CwAAACB3IrAAAAAAEDFmhQJy+bUbgFgXyqw8x+M1AACRoccCAAAAQMQILAAAAABEjFQoIAY0OXZltIsAAADiHIEFEAPyW6FoFwEAYhbTHQOhIRUKAAAAQMTosQAAAMhh9HogHhBYADFge75VlmSJlt8KWuXkJtEuDgAAiEMEFkAM2J7vR0tMOGgFfcUILADgOOM6KsD/MMYCAAAAQMQILAAAAABEjMACAAAAQMQILAAAAABEjMACAAAAQMQILAAAAABEjOlmAQA5Os1mrV49qWEAiAMEFgCAHMUc/wAQHwgsgBhQ3HeCJfoOWUErGu2iAAByAL2DyAsILIAYUDepQ7SLAADI5cEHaYnIaQQWAAAAUUbKIGIBgQUAIF00dgAAoSKwAHLQmMlLqV/kagQOAIDsQmABxIDf8s+xRPvf4O1YHm8RSqA2uEeL41IWAACQEoEFEAMOJPxliQkHraCvWLSLAgAA4hRX3gYAAAAQMQILAAAAABEjsAAAAAAQG2MsXnjhBRszZoxt377dTjvtNHvuuefsjDPOiHaxACBPY8YnAEBcBRbvvvuuDRo0yMaNG2dnnnmmjR071jp16mRr1qyxihUrRrt4AAAAMSE7TjZw9W7k6sDiqaeesr59+1qvXr3c3wowPvnkE3v99ddtyJAh0S4egBibkja3TEdLbwKAvCiUYxfBR/yKamBx9OhRW7ZsmQ0dOtS/LF++fHbBBRfYN998E82iASHhAngAAGTviRMCk7wrqoHFX3/9ZUlJSVapUqUUy/X36tWr06x/5MgRd/Ps2bPH3e/du/c4lBbx5pl3l1tekZj/iCUmHDXz5bfDSfvdMp8v2X9/+ND/lsHshxfHZVoNJ/X4V4aPb5j8dkTPl32HDrE5ACCHjtPIPl472+fzZbpugi+UtXLI1q1b7cQTT7Svv/7aWrVq5V9+zz332MKFC23JkiUp1n/44Ydt+PDhUSgpAAAAEL82bdpk1apVy709FieccILlz5/fduzYkWK5/q5cuXKa9ZUypYHenuTkZPvnn3+sfPnylpCQYNGK4qpXr+4qu1SpUlEpQyyhPqnL3Ij9krrMjdgvqcvcin0ztupSfRD79u2zqlWrZrpuVAOLQoUKWfPmzW3evHnWtWtXf7Cgv2+77bY06xcuXNjdApUpU8ZyA21sAgvqMzdi36QucyP2S+oyN2K/pD5zq1JRbmeWLl06b8wKpR6Inj17WosWLdy1KzTd7IEDB/yzRAEAAADI/aIeWFx77bX2559/2kMPPeQukNe0aVObNWtWmgHdAAAAAHKvqAcWorSnYKlPeYFSs4YNG5YmRQvUZ7Sxb1KXuRH7JXWZG7FfUp+5VeE81s6M6qxQAAAAAGJDvmgXAAAAAEDeR2ABAAAAIGIEFgAAAAAiRmARYOTIkdayZUsrWbKkVaxY0V1bY82aNZlW4n//+1+rX7++FSlSxJo0aWKffvppisc1jEWzXlWpUsWKFi1qF1xwgf32228Wy7JSl+PHj7dzzjnHypYt626qp2+//TbFOjfeeKO7GGLg7cILL7RYlpW6nDhxYpp60v4Z7/tlVuuzbdu2aepTt4svvjiu982XXnrJTj31VP/86q1atbKZM2dm+ByOl9lTlxwvs2+/5HiZfXXJsTJ0o0aNcr8TAwcOjKljJoFFgIULF9qtt95qixcvtjlz5lhiYqJ17NjRXVcjPV9//bVdd9111qdPH1uxYoVrpOj2448/+td5/PHH7dlnn7Vx48bZkiVLrHjx4tapUyc7fPiwxaqs1OWCBQtcXc6fP9+++eYbd6VJPWfLli0p1lNjbdu2bf7blClTLJZlpS5FPwKB9fTHH3+keDwe98us1ue0adNS1KW+3/nz57err746rvfNatWquR/HZcuW2dKlS61du3bWpUsX++mnn4Kuz/Ey++qS42X21aVwvMyeuuRYGZrvvvvOXn75ZRe0ZSRPHjM1KxSC27lzp2bM8i1cuDDdKrrmmmt8F198cYplZ555pu+mm25y/09OTvZVrlzZN2bMGP/ju3fv9hUuXNg3ZcqUuKn6UOoytWPHjvlKlizpmzRpkn9Zz549fV26dPHFs1DqcsKECb7SpUun+zj7ZXj1mdrTTz/t9s39+/f7l7Fv/k/ZsmV9r776atB643iZfXWZGsfLrNclx8uc2y85Vqa1b98+X926dX1z5szxnXfeeb4BAwb40pMXj5n0WGRgz5497r5cuXLprqMz6+p2CqRIUctl/fr17sJ/gevosuhnnnmmf514EEpdpnbw4EF3Njn1c3SmTiks9erVs/79+9vff/9t8STUuty/f7/VrFnT9fykPsPEfhl+fQZ67bXXrFu3bu7MUKB43jeTkpLsnXfecT0/SpcIhuNl9tVlahwvI6tLjpc5s19yrExLPeZKo03ddoyVY2auuEBebpScnOzy3lq3bm2NGzdOdz1t0NRXCdffWu497i1Lb51YF2pdpnbvvfda1apVU3xhlGpyxRVXWK1atWzdunV233332UUXXeS+QEpNiXWh1qUatq+//rrrZlXD+YknnrCzzz7bBRfq2ma/DK8+A2ncj7qh9YMZKF73zVWrVrlGhrrdS5QoYdOnT7eGDRsGXZfjZfbVZWocL7Nelxwvc2a/5FiZlgKz5cuXu1SoUOTFYyaBRQYRpRoPX3311fHdIjEoK3WpnE59AXUGOHDQsc4SezSISQ3n2rVru/Xat29vsS7UutSPQOAZJQUVDRo0cDmdI0aMOA4ljd19UwGF9r0zzjgjxfJ43TfVKFu5cqULYKdOnWo9e/Z041hCbRAj8rrkeBnZfsnxMmf2S46VKW3atMkGDBjgxvalnkwllpAKFcRtt91mM2bMcIOIdXY3I5UrV7YdO3akWKa/tdx73FuW3jqxLJy69Ojsun4oP/vss0wHNp188sl2wgkn2Nq1ay3WZaUuPQULFrTTTz/dX0/xvl9mtT6VAqCAVwPpMhMv+2ahQoWsTp061rx5czfj1mmnnWbPPPNM0HU5XmZfXXo4XmZfXXo4XkZelxwr09IA+J07d1qzZs2sQIEC7qYATQOv9X+lmsXCMZPAItWUXWpsqJvv888/dykNmdGZjnnz5qVYpmjUO1us19DGDVxn7969buR+qDmK8VKX3uwGOqM+a9Ysa9GiRabrb9682eWxa5q1WJXVugykA5a6s716itf9MtL61LR/R44csR49emS6bjzsm+mll6mOguF4mX11KRwvs68uA3G8jLwuOVampZ5r/Q6r98e7qZ3TvXt39/9gKbN58pgZlSHjuVT//v3dTDoLFizwbdu2zX87ePCgf53rr7/eN2TIEP/fixYt8hUoUMD3xBNP+H755RffsGHDfAULFvStWrXKv86oUaN8ZcqU8X344Ye+H374wc1qVKtWLd+hQ4d8sSordal6KlSokG/q1KkpnqMZFET3d999t++bb77xrV+/3jd37lxfs2bN3OwKhw8f9sWqrNTl8OHDfbNnz/atW7fOt2zZMl+3bt18RYoU8f30009xvV9mtT49bdq08V177bVplsfrvqk60mxa+szah/R3QkKC77PPPnOPc7zMubrkeJl9dcnxMvvq0sOxMjSpZ4WKhWMmgUVgZZgFvWkqusCdQNNKBnrvvfd8p5xyimsUN2rUyPfJJ5+keFzTgT344IO+SpUquSnA2rdv71uzZo0vlmWlLmvWrBn0OfoiiRp+HTt29FWoUMF9sbR+3759fdu3b/fFsqzU5cCBA301atRw+6T2u86dO/uWL1/ui/f9MpLv+erVq9163g9qoHjdN3v37u0+q/YzfXbtQ4H1w/Ey5+qS42X21SXHy+z9jnOszHpgEQvHzAT9E52+EgAAAACxgjEWAAAAACJGYAEAAAAgYgQWAAAAACJGYAEAAAAgYgQWAAAAACJGYAEAAAAgYgQWAAAAACJGYAEAAAAgYgQWABBlEydOtDJlylhecf3119tjjz2WLa+1YcMGS0hIsJUrV1qsatu2rQ0cODDL+8PDDz9sTZs2jbgcP//8s1WrVs0OHDgQ8WsBQDAEFgCQB3kNcu9Wrlw5O++88+zLL79Ms+4///zjGrY1a9a0QoUKWdWqVa137962cePGFOvdeOON1rVr1wzf9/vvv7dPP/3U7rjjjmz/TMhZDRs2tLPOOsueeuopqhpAjiCwAIAcdPTo0Ryt37lz59q2bdvsiy++cAHDJZdcYjt27EgRVKgxqfXGjRtna9eutXfeecfdt2zZ0n7//few3u+5556zq6++2kqUKGHxLKe3a3ZLTEx097169bKXXnrJjh07Fu0iAYhBBBYA4srUqVOtSZMmVrRoUStfvrxdcMEFKVJDXn31VWvQoIEVKVLE6tevby+++GKK59977712yimnWLFixezkk0+2Bx980N9oC0xb0evUqlXLvY7s3r3bbrrpJqtUqZJb1rhxY5sxY0aK1549e7Z7bzXaL7zwQhcwZEafoXLlyu717rvvPtu7d68tWbLE//j9999vW7dudYHFRRddZDVq1LBzzz3XvVfBggXt1ltvDbnukpKSXP1deumlKZarjurWres+lz7fVVdd5X9s1qxZ1qZNG5fao7Iq8Fm3bl2G7/Pjjz+6sqoe9HpKvfrrr79C3oaBFixY4Hp0PvnkEzv11FNdGRVo6T0CffXVV3bOOee416xevbrrkQl8zZNOOslGjBhhN9xwg5UqVcr69esX9P30HK2jslepUsWefPLJNOvs2rXLrVO2bFm3H+mz/vbbbxaq7777zjp06GAnnHCClS5d2vVULV++PMU6+swKIC677DIrXry4Pfroo265nqdgc+HChSG/HwCEisACQNxQQ/26665zaUC//PKLa3ReccUV5vP53ONvvfWWPfTQQ64Rpsc1jkCBw6RJk/yvUbJkSZcDr3z1Z555xsaPH29PP/10ivdRb8D7779v06ZNc2MHkpOTXeNx0aJFNnnyZPfcUaNGWf78+f3POXjwoD3xxBP25ptvut4HpSndfffdIX+2Q4cO2RtvvOH+r3Qn0fuqd6J79+4u+AikBvQtt9ziAgw1NEPxww8/2J49e6xFixb+ZUuXLnWN8P/85z+2Zs0aF0gocAlsaA8aNMitN2/ePMuXL59dfvnlrmzBKABr166dnX766e45ej31wFxzzTUhbcP0DB482DXy1SivUKGCC468gFCBjgK5K6+80n3Gd9991wUat912W4rX0PY57bTTbMWKFW6/SO991Gj/8MMP7bPPPnPlS93oV8qZPttHH31k33zzjSt7586dUwSoGdm3b5/17NnTlXHx4sUuqNPztTyQglzV9apVq1x9efuGAt9gKXMAEDEfAMSJZcuWqfXp27BhQ9DHa9eu7Xv77bdTLBsxYoSvVatW6b7mmDFjfM2bN/f/PWzYMF/BggV9O3fu9C+bPXu2L1++fL41a9YEfY0JEya4cq1du9a/7IUXXvBVqlQp3fddv369e07RokV9xYsX9yUkJLi/VZajR4+6dbZv3+6WPf3000FfY9q0ae7xJUuWuL979uzp69KlS7rvOX36dF/+/Pl9ycnJ/mXvv/++r1SpUr69e/f6QvHnn3+691y1alWKz7FixQp/fXfs2DHFczZt2uTWUf1ltg1Tmz9/vlv/nXfe8S/7+++/Xb29++677u8+ffr4+vXrl+J5X375pdtmhw4dcn/XrFnT17Vr1wzfa9++fb5ChQr53nvvvTTvNWDAAPf3r7/+6sqzaNEi/zp//fWXW8d7nvaH0qVLp9inTjvttHTfNykpyVeyZEnfxx9/7F+m9xg4cGDQ9S+//HLfjTfemOFnAYCsoMcCQNzQ2eb27du7NBqNE1Bvg9JSvDPrOnPdp08fl8bi3R555JEUqTs6m926dWvXA6DHH3jggTSDoDVIWmfFPeq10Gw8SqFKj1Jiateu7f9baTQ7d+7M9DOpPDqDrh6SOnXquN4UpTgFyuxsfji9IoULF3ZpNh6l1ujzKi1MKUvq9VHvi0cpPuph0ONKIVJKkaSus8DB4fPnz0+xDZSSJtoOGW3DjLRq1cr/fw10r1evnuvx8N5T9Rb4np06dXK9KuvXr/c/L7CnJhiVT2MvzjzzzDTv5dF7FihQIMU6SucKLE9m1IPTt29f11OhVCjV6/79+9PUaXrlVW9V4DYCgOxSINteCQByOaUezZkzx77++muXpqKByBqDoDEJatiLGqqBjT7veaK0FaUVDR8+3DU81ahTqlHqPHrltKduyGUmdTCgxnsoAYHGA6iBqZsG5Cr1ReMHFAAouNHYhvQarFqu91FAEgrl9KtBqsazl26l1DCl+ijlR3WqVDKl4CjlSO+tlCMFHqpXDS5XY13jQdIb/KwGsp4zevToNI8p2MpoG2pMS1boPTX+JdhMVxqTkt52jRalQf39998uFU91q22twCl1naZXXqW+BQaxAJBd6LEAEFfUkFaPg4IDnelXA3n69OlukLAavpolSQ3twJvXYFVjVg05NWR1NliN+T/++CPT99Sg4c2bN9uvv/6ao59Ng6Z1NtwbcK7xDBqb8Pbbb9v27dvT9D5oPQVIOqseCu9aChojEkjvqQHUjz/+uBujoKlwP//8c9f41bgL9eqol0ED0zPrXWjWrJn99NNPrmcj9XbwGsrpbcOMaCyCR2XQtlB5vPfUZ0r9frp5AVQo1FhXgBg4eN57L4/eUwFg4DpePWk62FBorI6CII2raNSokQssAge3Z0aBp8awAEB2I7AAEDfUmNOAbA2cVdqIBlf/+eef/gamGqojR460Z5991jUGNeh1woQJ/nn/FUjoeeqlUNqL1susQSuatUcDmjU4WGfblV4zc+ZMNzA5O6nBrQanBoZ7qS76vErbUsqS3nPTpk1ucLgCCg0WfuGFF0J+ffWAqBGuQcMezWylelC6l4IsDSBXr4RSezTrkdJ8XnnlFTegXcGGBnJnRLNU6Yy60qfU66F61gBzTZOqWaky24bp0eByDR5Xo1qDp9X74l2zQzN9KWjUYG19DqVvafB16sHbmVEKlVLpNIBbn9V7LwV4Hu1DXbp0calMqkelYfXo0cNOPPFEtzwUeg0N8lePk+pDvWih9IqJgr4tW7a4QBAAshuBBYC4oVx0Nap1plfjHXQmXWlMmrFJ/v3vf7tpYhVMKIdfAYFy770eC03deeedd7oGp87eqzGa3uxAqWkMhK4boQazzkzfc889rqGcE2kyChief/5597ca9jpbf/7557t0H51VVy+G7tVw19iHcKiONI7Co3QnNe41k5Ma97pWxpQpU9yZdDWoFYQtW7bMpT+p7saMGZPh66vXSGfkVTcdO3Z020EX99P76PUy24bpUbA1YMAAa968ueu9+fjjj/29EepR0kxOCiY15azO5iulS2UJlz6fXkPpXGq8a6pdvWcg7V9apql3lcKklDdddDB1Olx6XnvtNdcToiBP41oUTFasWDGk52rbqF7V8wYA2S1BI7iz/VUBADFJKVTqjdCg8cAB0bmVxn4oqFJDXMFJPNMYDPV2KDVOqWQAkN3osQAAhEwpN0p3CienH7mDUsd0EUWCCgA5hVmhAABhadu2LTWWB3kD0gEgp5AKBQAAACBipEIBAAAAiBiBBQAAAICIEVgAAAAAiBiBBQAAAICIEVgAAAAAiBiBBQAAAICIEVgAAAAAiBiBBQAAAICIEVgAAAAAsEj9P+exS1nTkNUgAAAAAElFTkSuQmCC",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "fig, ax = plt.subplots(figsize=(8, 4.5))\n",
+ "bins = np.linspace(2.0, 4.0, 80)\n",
+ "ax.hist(roi_default.ravel(), bins=bins, density=True, alpha=0.55,\n",
+ " label=\"default roi_m prior\", color=\"#c44e52\")\n",
+ "ax.hist(roi_calibrated.ravel(), bins=bins, density=True, alpha=0.65,\n",
+ " label=\"calibrated roi_m prior\", color=\"#4c72b0\")\n",
+ "ax.axvline(BETA_SEARCH, color=\"black\", lw=2, label=\"simulated truth (2.5)\")\n",
+ "ax.axvline(prior.roi_mean, color=\"#55a868\", lw=2, ls=\"--\", label=\"DiD experiment ROI\")\n",
+ "ax.set_xlabel(\"search ROI (sales per dollar)\")\n",
+ "ax.set_ylabel(\"posterior density\")\n",
+ "ax.set_title(\"What the experiment prior buys: Meridian search-ROI posterior\")\n",
+ "ax.legend()\n",
+ "plt.tight_layout()\n",
+ "plt.show()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "204a3490",
+ "metadata": {},
+ "source": [
+ "With the default prior the national model overcredits search - the launches rode\n",
+ "the demand ramp, and nationally there is no holdout to say otherwise. The\n",
+ "calibrated prior carries the geo experiment's information into the model: the\n",
+ "posterior centers on the truth and tightens."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "ed4f381a",
+ "metadata": {},
+ "source": [
+ "## Acceptance checks\n",
+ "\n",
+ "Sampler health first - the whole-model maximum R-hat (Meridian's own health-check\n",
+ "convention treats R-hat < 1.2 as converged) plus a tighter bar on the reported\n",
+ "`roi_m` for the experiment channel - then the calibration claims with slack\n",
+ "margins against draw noise."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "id": "fb43c031",
+ "metadata": {
+ "execution": {
+ "iopub.execute_input": "2026-08-16T19:54:22.712911Z",
+ "iopub.status.busy": "2026-08-16T19:54:22.712807Z",
+ "iopub.status.idle": "2026-08-16T19:54:22.823473Z",
+ "shell.execute_reply": "2026-08-16T19:54:22.822974Z"
+ }
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "default: whole-model max Rhat=1.0038, Rhat(roi_m, search)=0.9998\n"
+ ]
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "calibrated: whole-model max Rhat=1.0175, Rhat(roi_m, search)=1.0045\n",
+ "error: 0.79 -> 0.04; 90% width: 0.55 -> 0.15\n",
+ "all acceptance checks passed\n"
+ ]
+ }
+ ],
+ "source": [
+ "for name, mer in [(\"default\", mer_default), (\"calibrated\", mer_calibrated)]:\n",
+ " post = mer.inference_data.posterior\n",
+ " with warnings.catch_warnings():\n",
+ " # constant/deterministic posterior variables have undefined R-hat;\n",
+ " # arviz emits RuntimeWarnings computing them - excluded below.\n",
+ " warnings.simplefilter(\"ignore\", RuntimeWarning)\n",
+ " rh_all = az.rhat(post)\n",
+ " finite_rhats = np.concatenate(\n",
+ " [np.asarray(rh_all[v]).ravel() for v in rh_all.data_vars]\n",
+ " )\n",
+ " max_rhat_model = float(finite_rhats[np.isfinite(finite_rhats)].max())\n",
+ " rhat_roi = float(az.rhat(post[\"roi_m\"]).sel(media_channel=\"search\")[\"roi_m\"])\n",
+ " print(f\"{name}: whole-model max Rhat={max_rhat_model:.4f}, \"\n",
+ " f\"Rhat(roi_m, search)={rhat_roi:.4f}\")\n",
+ " # Meridian's own health-check convention treats Rhat < 1.2 as converged;\n",
+ " # the reported ROI parameter is held to a tighter bar.\n",
+ " assert np.isfinite(max_rhat_model) and max_rhat_model < 1.2, (\n",
+ " f\"{name} fit: some model parameter has not converged \"\n",
+ " f\"(max Rhat {max_rhat_model:.4f})\"\n",
+ " )\n",
+ " assert np.isfinite(rhat_roi) and rhat_roi < 1.05, (\n",
+ " f\"{name} fit: roi_m chains disagree (Rhat {rhat_roi:.4f}) - sampling is \"\n",
+ " \"unhealthy here\"\n",
+ " )\n",
+ "\n",
+ "flat_d, flat_c = roi_default.ravel(), roi_calibrated.ravel()\n",
+ "err_default = abs(flat_d.mean() - BETA_SEARCH)\n",
+ "err_calibrated = abs(flat_c.mean() - BETA_SEARCH)\n",
+ "width_default = float(np.quantile(flat_d, 0.95) - np.quantile(flat_d, 0.05))\n",
+ "width_calibrated = float(np.quantile(flat_c, 0.95) - np.quantile(flat_c, 0.05))\n",
+ "print(f\"error: {err_default:.2f} -> {err_calibrated:.2f}; \"\n",
+ " f\"90% width: {width_default:.2f} -> {width_calibrated:.2f}\")\n",
+ "\n",
+ "assert err_calibrated < err_default, (\n",
+ " \"calibrated ROI mean is not closer to the simulated truth - the experiment \"\n",
+ " \"prior did not correct the bias in this environment\"\n",
+ ")\n",
+ "assert width_calibrated < 0.8 * width_default, (\n",
+ " \"calibrated 90% interval is not clearly narrower - the experiment prior did \"\n",
+ " \"not add information in this environment\"\n",
+ ")\n",
+ "print(\"all acceptance checks passed\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "4aa3ca94",
+ "metadata": {},
+ "source": [
+ "## Where to go next\n",
+ "\n",
+ "- **PyMC-Marketing**: notebook 29 runs the same hand-off as a *lift test* for a\n",
+ " spend-boost experiment - the other classic geo design - including the grain\n",
+ " reconciliation a non-total export needs.\n",
+ "- **API reference**: the `diff_diff.mmm` module documents `to_meridian_roi_prior`\n",
+ " (spend-weighted pooling of multiple experiments, `se_widening`, `mroi_m`),\n",
+ " `meridian_calibration_mask`, and `MeridianROIPrior.to_code` - including why\n",
+ " channel order must match your `InputData` exactly, and Meridian's own guidance\n",
+ " on when to prefer `full_model_window=True` over a calibration mask.\n",
+ "- **Estimator choices**: the DiD side works with any diff-diff estimator that\n",
+ " supports `aggregate('total')` (CallawaySantAnna, EfficientDiD, ImputationDiD,\n",
+ " TwoStageDiD)."
+ ]
+ }
+ ],
+ "metadata": {
+ "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.13"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/docs/tutorials/README.md b/docs/tutorials/README.md
index 4adbb65b..9ee2ebb1 100644
--- a/docs/tutorials/README.md
+++ b/docs/tutorials/README.md
@@ -159,6 +159,20 @@ A loyalty program looks dead on mean DiD ($0.22, p = 0.90) but in truth lifted t
- Covariate-composition confounding fixed with `covariates=` (quantile-regression conditioning, qte `xformla` parity), and the `practitioner_next_steps()` close
- Companion drift-test file (`tests/test_t27_cic_distributional_effects_drift.py`); fully self-contained (runs live, no committed data files)
+### 29. MMM Calibration with PyMC-Marketing Lift Tests (`29_mmm_calibration_pymc.ipynb`)
+Executed end-to-end against the real framework: a staggered geo spend-boost experiment measured with `CallawaySantAnna`, exported via `to_pymc_marketing_lift_test`, and fed to `MMM.add_lift_test_measurements` - the demand-chasing-biased MMM's search-ROI posterior moves to the simulated truth and tightens:
+- Grain reconciliation for a national-weekly lift row (`aggregate('simple')` + `scale=G_treated`), the `aggregate('total')` campaign headline, and the double-count guardrail
+- With/without-calibration money plot plus in-notebook sampler-health and posterior acceptance asserts
+- Requires `pip install "pymc-marketing==1.0.0"` on Python 3.12 (this tutorial only); cannot share an environment with the Meridian tutorial (conflicting `arviz` pins); fully self-contained data
+- Companion drift-test file (`tests/test_t29_mmm_calibration_pymc_drift.py`)
+
+### 30. MMM Calibration with Meridian ROI Priors (`30_mmm_calibration_meridian.ipynb`)
+Executed end-to-end against the real framework: a geo-holdout launch experiment (the `roi_m` estimand exactly), exported via `to_meridian_roi_prior` from the `aggregate('total')` container, with `meridian_calibration_mask` and the generated `to_code()` snippet run verbatim:
+- A national Meridian model that cannot see the geo holdout overcredits the launch; the experiment prior centers its ROI posterior on the truth and tightens it
+- Window/spend/mass alignment asserts, default-vs-calibrated money plot, and sampler-health checks
+- Requires `pip install "google-meridian==1.8.0"` on Python 3.12 (this tutorial only); cannot share an environment with the PyMC-Marketing tutorial (conflicting `arviz` pins); fully self-contained data
+- Companion drift-test file (`tests/test_t30_mmm_calibration_meridian_drift.py`)
+
## Running the Notebooks
1. Install diff-diff with dependencies:
diff --git a/docs/tutorials/index.rst b/docs/tutorials/index.rst
index f484960d..e91c47da 100644
--- a/docs/tutorials/index.rst
+++ b/docs/tutorials/index.rst
@@ -1,5 +1,5 @@
.. meta::
- :description: Hands-on diff-diff tutorials — 28 Jupyter notebooks covering basic 2x2 DiD, staggered adoption, synthetic DiD, power analysis, and business applications.
+ :description: Hands-on diff-diff tutorials — 31 Jupyter notebooks covering basic 2x2 DiD, staggered adoption, synthetic DiD, power analysis, and business applications.
:keywords: DiD tutorial, difference-in-differences examples, causal inference notebooks
Tutorials
@@ -79,6 +79,20 @@ Practitioner walkthroughs built around marketing and policy scenarios.
See which quantiles moved when the average hides the action, with
Changes-in-Changes.
+ .. grid-item-card:: MMM Calibration: PyMC-Marketing
+ :link: 29_mmm_calibration_pymc
+ :link-type: doc
+
+ Turn a geo spend-boost experiment into a PyMC-Marketing lift test and
+ watch the MMM's ROI posterior correct - executed end-to-end.
+
+ .. grid-item-card:: MMM Calibration: Meridian
+ :link: 30_mmm_calibration_meridian
+ :link-type: doc
+
+ Export a geo-holdout launch as Meridian's ROI prior and calibration
+ mask, then fit the real model - executed end-to-end.
+
.. toctree::
:maxdepth: 1
:caption: Business Applications
@@ -93,6 +107,8 @@ Practitioner walkthroughs built around marketing and policy scenarios.
Spillover DiD (TVA) <23_spillover_tva>
Composition Drift & Calibration <26_composition_drift_calibration>
Distributional Effects (CiC) <27_cic_distributional_effects>
+ MMM Calibration: PyMC-Marketing <29_mmm_calibration_pymc>
+ MMM Calibration: Meridian <30_mmm_calibration_meridian>
Fundamentals
------------
diff --git a/tests/test_guides.py b/tests/test_guides.py
index a2ba3694..ef6eb953 100644
--- a/tests/test_guides.py
+++ b/tests/test_guides.py
@@ -39,6 +39,25 @@ def test_content_stability_practitioner_workflow():
assert "8-step" in get_llm_guide("practitioner").lower()
+def test_practitioner_step8_mmm_handoff_route_qualified():
+ """The Step-8 MMM hand-off must stay ROUTE-QUALIFIED: both exporters named,
+ the total-container route AND the scaled simple/group route both present
+ inside the Step-8 slice - an unqualified total-first claim (or a total-route
+ deletion) must fail here, not ship silently."""
+ guide = get_llm_guide("practitioner")
+ start = guide.index("## Step 8")
+ end = guide.index("## Common Pitfalls")
+ step8 = guide[start:end]
+ for needle in (
+ "to_pymc_marketing_lift_test",
+ "to_meridian_roi_prior",
+ "aggregate('total')",
+ "scale=",
+ "SAME\n observation",
+ ):
+ assert needle in step8, f"Step-8 MMM hand-off lost {needle!r}"
+
+
def test_content_stability_self_reference_after_rewrite():
assert "get_llm_guide" in get_llm_guide("concise")
diff --git a/tests/test_mmm_interop_meridian.py b/tests/test_mmm_interop_meridian.py
new file mode 100644
index 00000000..953106cb
--- /dev/null
+++ b/tests/test_mmm_interop_meridian.py
@@ -0,0 +1,148 @@
+"""Meridian interop smoke tests (schema/defaults-drift canary).
+
+These tests exercise :func:`diff_diff.mmm.to_meridian_roi_prior`,
+:meth:`MeridianROIPrior.to_code`, and :func:`diff_diff.mmm.meridian_calibration_mask`
+against the REAL installed google-meridian API - no sampling anywhere. Skipped
+unless google-meridian >= 1.8 is installed (NOT a diff-diff dependency); the
+dedicated ``mmm-interop.yml`` CI job installs it and runs this file, with an
+import-canary step so an install failure cannot be silently skipped into a green
+run.
+"""
+
+import numpy as np
+import pandas as pd
+import pytest
+
+pytest.importorskip(
+ "meridian",
+ minversion="1.8",
+ reason="google-meridian>=1.8 required (installed only in the mmm-interop CI job)",
+)
+
+import tensorflow_probability as tfp # noqa: E402
+from meridian.data import data_frame_input_data_builder as m_dfb # noqa: E402
+from meridian.model import model as m_model # noqa: E402
+from meridian.model import prior_distribution, spec # noqa: E402
+
+from diff_diff.mmm import ( # noqa: E402
+ _MERIDIAN_PARAM_DEFAULTS,
+ meridian_calibration_mask,
+ to_meridian_roi_prior,
+)
+
+TIMES = [str(d.date()) for d in pd.date_range("2024-01-01", periods=8, freq="W-MON")]
+CHANNELS = ["search", "tv"]
+
+
+def _prior():
+ return to_meridian_roi_prior(
+ incremental_outcome=5000.0, incremental_outcome_se=800.0, spend=2000.0
+ )
+
+
+def _mask():
+ return meridian_calibration_mask(
+ media_times=TIMES,
+ media_channels=CHANNELS,
+ channel="search",
+ window=(TIMES[4], TIMES[-1]),
+ )
+
+
+class TestMeridianDefaultsPin:
+ def test_param_defaults_match_live_meridian(self):
+ # The multi-channel to_code() snippet keeps non-experiment channels on what
+ # it BELIEVES are Meridian's default LogNormal params - if Meridian changes
+ # its defaults, this canary must fail rather than silently rot the snippet.
+ # Meridian stores float32 tensors: compare after float32 rounding.
+ live = prior_distribution.PriorDistribution()
+ for param, (mu, sigma) in _MERIDIAN_PARAM_DEFAULTS.items():
+ dist = getattr(live, param)
+ live_mu = np.float32(np.asarray(dist.parameters["loc"]))
+ live_sigma = np.float32(np.asarray(dist.parameters["scale"]))
+ assert (live_mu, live_sigma) == (np.float32(mu), np.float32(sigma)), (
+ f"Meridian default drift for {param}: "
+ f"ours=({mu}, {sigma}) live=({live_mu}, {live_sigma})"
+ )
+
+
+class TestToCodeExecutes:
+ def test_multi_channel_snippet_execs_with_value_retention(self):
+ prior = _prior()
+ mask = _mask()
+ code = prior.to_code(channel="search", media_channels=CHANNELS, roi_calibration_period=mask)
+ ns: dict = {}
+ exec(code, ns) # noqa: S102 - executing our own generated snippet is the test
+ model_spec = ns["model_spec"]
+ assert isinstance(model_spec, spec.ModelSpec)
+ # Value retention: the exec'd prior's roi_m LogNormal carries our mu/sigma
+ # for the experiment channel and Meridian's defaults for the other channel.
+ dist = model_spec.prior.roi_m
+ loc = np.asarray(dist.parameters["loc"], dtype=float)
+ scale = np.asarray(dist.parameters["scale"], dtype=float)
+ default_mu, default_sigma = _MERIDIAN_PARAM_DEFAULTS["roi_m"]
+ np.testing.assert_allclose(loc, [prior.mu, default_mu], rtol=1e-6)
+ np.testing.assert_allclose(scale, [prior.sigma, default_sigma], rtol=1e-6)
+ # And the serialized mask round-trips exactly.
+ assert np.array_equal(np.asarray(model_spec.roi_calibration_period), mask)
+
+ def test_single_channel_snippet_execs(self):
+ prior = _prior()
+ code = prior.to_code(single_channel=True, full_model_window=True)
+ ns: dict = {}
+ exec(code, ns) # noqa: S102
+ model_spec = ns["model_spec"]
+ assert isinstance(model_spec, spec.ModelSpec)
+ dist = model_spec.prior.roi_m
+ assert float(np.asarray(dist.parameters["loc"])) == pytest.approx(prior.mu)
+ assert float(np.asarray(dist.parameters["scale"])) == pytest.approx(prior.sigma)
+ assert model_spec.roi_calibration_period is None
+
+
+class TestModelAcceptsExports:
+ def test_mask_and_prior_accepted_by_meridian(self):
+ prior = _prior()
+ mask = _mask()
+ roi_prior = tfp.distributions.LogNormal(
+ [prior.mu, _MERIDIAN_PARAM_DEFAULTS["roi_m"][0]],
+ [prior.sigma, _MERIDIAN_PARAM_DEFAULTS["roi_m"][1]],
+ name="roi_m",
+ )
+ model_spec = spec.ModelSpec(
+ prior=prior_distribution.PriorDistribution(roi_m=roi_prior),
+ media_prior_type="roi",
+ roi_calibration_period=mask,
+ )
+ # Build InputData through the PUBLIC DataFrameInputDataBuilder (the
+ # route notebook 30 teaches) rather than Meridian's internal test-utils
+ # fixtures, so upstream fixture churn cannot break this canary.
+ rng = np.random.default_rng(1)
+ rows = []
+ for geo in ["g0", "g1", "g2", "g3"]:
+ for t in TIMES:
+ rows.append(
+ {
+ "geo": geo,
+ "time": t,
+ "sales": 1000.0 + 50.0 * rng.random(),
+ "population": 1_000_000.0,
+ "search_spend": 100.0 + 10.0 * rng.random(),
+ "tv_spend": 200.0 + 10.0 * rng.random(),
+ }
+ )
+ frame = pd.DataFrame(rows)
+ input_data = (
+ m_dfb.DataFrameInputDataBuilder(kpi_type="revenue")
+ .with_kpi(frame, kpi_col="sales")
+ .with_population(frame)
+ .with_media(
+ frame,
+ media_cols=["search_spend", "tv_spend"],
+ media_spend_cols=["search_spend", "tv_spend"],
+ media_channels=CHANNELS,
+ )
+ .build()
+ )
+ mer = m_model.Meridian(input_data=input_data, model_spec=model_spec)
+ assert mer is not None
+ assert np.array_equal(np.asarray(mer.model_spec.roi_calibration_period), mask)
diff --git a/tests/test_mmm_interop_pymc.py b/tests/test_mmm_interop_pymc.py
new file mode 100644
index 00000000..812f9720
--- /dev/null
+++ b/tests/test_mmm_interop_pymc.py
@@ -0,0 +1,174 @@
+"""PyMC-Marketing interop smoke tests (schema-drift canary).
+
+These tests exercise :func:`diff_diff.mmm.to_pymc_marketing_lift_test` against the
+REAL installed pymc-marketing API - no sampling anywhere. They are skipped unless
+pymc-marketing >= 1.0 is installed (it is NOT a diff-diff dependency); the dedicated
+``mmm-interop.yml`` CI job installs it and runs this file, with an import-canary
+step so an install failure cannot be silently skipped into a green run.
+
+Value-retention oracle: pymc-marketing 1.0 stores an added lift test as a
+``lift_measurements`` potential inside the PyMC model graph and retains no
+DataFrame attribute, so retention is asserted FUNCTIONALLY - the measured values
+enter the likelihood (perturbing ``delta_y``/``sigma`` changes the model's logp at
+a fixed parameter point; an identical DataFrame reproduces it exactly).
+"""
+
+import numpy as np
+import pandas as pd
+import pytest
+
+pytest.importorskip(
+ "pymc_marketing",
+ minversion="1.0",
+ reason="pymc-marketing>=1.0 required (installed only in the mmm-interop CI job)",
+)
+
+from pymc_marketing.mmm import MMM, GeometricAdstock, LogisticSaturation # noqa: E402
+
+from diff_diff.mmm import to_pymc_marketing_lift_test # noqa: E402
+
+LIFT_COLUMNS = ["channel", "x", "delta_x", "delta_y", "sigma"]
+
+
+def _national_frame(n: int = 40) -> tuple[pd.DataFrame, pd.Series]:
+ rng = np.random.default_rng(7)
+ dates = pd.date_range("2024-01-01", periods=n, freq="W-MON")
+ search = (1500 + 450 * (np.arange(n) > n // 2)).astype(float)
+ tv = 3000 + 100 * rng.random(n)
+ y = pd.Series(45000 + 2.0 * search + 0.5 * tv + rng.normal(0, 200, n), name="y")
+ return pd.DataFrame({"date": dates, "search": search, "tv": tv}), y
+
+
+def _build_mmm(X: pd.DataFrame, y: pd.Series, **kwargs) -> MMM:
+ mmm = MMM(
+ date_column="date",
+ channel_columns=["search", "tv"],
+ adstock=GeometricAdstock(l_max=2),
+ saturation=LogisticSaturation(),
+ **kwargs,
+ )
+ mmm.build_model(X, y)
+ return mmm
+
+
+def _model_logp(mmm: MMM) -> float:
+ model = mmm.model
+ return float(model.compile_logp()(model.initial_point()))
+
+
+class TestLiftTestSchema:
+ def test_exporter_frame_accepted(self):
+ X, y = _national_frame()
+ mmm = _build_mmm(X, y)
+ df = to_pymc_marketing_lift_test(
+ channel="search", x=1500.0, delta_x=450.0, delta_y=930.0, sigma=110.0
+ )
+ assert list(df.columns) == LIFT_COLUMNS
+ mmm.add_lift_test_measurements(df)
+ assert "lift_measurements" in mmm.model.named_vars
+
+ def test_wrong_schema_rejected(self):
+ # Regression pin: the API validates the exporter's column names - a renamed
+ # column must FAIL, otherwise the schema contract is no longer exercised.
+ X, y = _national_frame()
+ mmm = _build_mmm(X, y)
+ df = to_pymc_marketing_lift_test(
+ channel="search", x=1500.0, delta_x=450.0, delta_y=930.0, sigma=110.0
+ ).rename(columns={"delta_y": "lift"})
+ with pytest.raises(KeyError, match="delta_y"):
+ mmm.add_lift_test_measurements(df)
+
+ def test_multiple_rows_accepted(self):
+ X, y = _national_frame()
+ mmm = _build_mmm(X, y)
+ df = to_pymc_marketing_lift_test(
+ channel=["search", "tv"],
+ x=[1500.0, 3000.0],
+ delta_x=[450.0, 300.0],
+ delta_y=[930.0, 240.0],
+ sigma=[110.0, 60.0],
+ )
+ mmm.add_lift_test_measurements(df)
+ assert "lift_measurements" in mmm.model.named_vars
+
+
+class TestLiftTestDims:
+ def test_dims_row_round_trips(self):
+ # dims= is a first-class exporter parameter emitting real DataFrame
+ # columns; verify a dims-carrying row against the real multidimensional API.
+ rng = np.random.default_rng(11)
+ n = 30
+ dates = pd.date_range("2024-01-01", periods=n, freq="W-MON")
+ frames = []
+ for geo in ["east", "west"]:
+ search = (700 + 200 * (np.arange(n) > n // 2)).astype(float)
+ tv = 1500 + 50 * rng.random(n)
+ frames.append(
+ pd.DataFrame(
+ {
+ "date": dates,
+ "geo": geo,
+ "search": search,
+ "tv": tv,
+ "y": 20000 + 2.0 * search + 0.5 * tv + rng.normal(0, 100, n),
+ }
+ )
+ )
+ panel = pd.concat(frames, ignore_index=True)
+ X = panel[["date", "geo", "search", "tv"]]
+ y = panel["y"]
+ mmm = MMM(
+ date_column="date",
+ channel_columns=["search", "tv"],
+ dims=("geo",),
+ adstock=GeometricAdstock(l_max=2),
+ saturation=LogisticSaturation(),
+ )
+ mmm.build_model(X, y)
+ df = to_pymc_marketing_lift_test(
+ channel="search",
+ x=700.0,
+ delta_x=200.0,
+ delta_y=410.0,
+ sigma=55.0,
+ dims={"geo": "east"},
+ )
+ assert list(df.columns) == ["channel", "geo", "x", "delta_x", "delta_y", "sigma"]
+ mmm.add_lift_test_measurements(df)
+ assert "lift_measurements" in mmm.model.named_vars
+
+
+class TestLiftTestValueRetention:
+ """Functional retention: the exported values actually enter the likelihood."""
+
+ def test_identical_frames_identical_logp(self):
+ X, y = _national_frame()
+ df = to_pymc_marketing_lift_test(
+ channel="search", x=1500.0, delta_x=450.0, delta_y=930.0, sigma=110.0
+ )
+ logps = []
+ for _ in range(2):
+ mmm = _build_mmm(X, y)
+ mmm.add_lift_test_measurements(df.copy())
+ logps.append(_model_logp(mmm))
+ assert logps[0] == pytest.approx(logps[1], rel=0, abs=0)
+
+ @pytest.mark.parametrize(
+ "column,changed",
+ [("x", 3000.0), ("delta_x", 900.0), ("delta_y", 1860.0), ("sigma", 220.0)],
+ )
+ def test_perturbed_values_change_logp(self, column, changed):
+ X, y = _national_frame()
+ base = to_pymc_marketing_lift_test(
+ channel="search", x=1500.0, delta_x=450.0, delta_y=930.0, sigma=110.0
+ )
+ perturbed = base.copy()
+ perturbed[column] = changed
+ mmm_a = _build_mmm(X, y)
+ mmm_a.add_lift_test_measurements(base)
+ mmm_b = _build_mmm(X, y)
+ mmm_b.add_lift_test_measurements(perturbed)
+ assert _model_logp(mmm_a) != _model_logp(mmm_b), (
+ f"changing {column} did not change the model logp - the lift value is "
+ "being ignored or remapped by the framework"
+ )
diff --git a/tests/test_openai_review.py b/tests/test_openai_review.py
index 468f5f2f..afac8902 100644
--- a/tests/test_openai_review.py
+++ b/tests/test_openai_review.py
@@ -2203,6 +2203,7 @@ def test_rust_test_yml_pr_filter_covers_ai_review_surfaces(self, workflow_paths)
".github/workflows/notebooks.yml",
".github/workflows/ci-gate.yml",
".github/workflows/release-build-check.yml",
+ ".github/workflows/mmm-interop.yml",
)
def test_rust_test_yml_push_filter_covers_guarded_workflows(self, workflow_paths):
@@ -2269,6 +2270,7 @@ class TestCiWorkflowLabelEventGuard:
EXPECTED_JOBS = {
"rust-test.yml": ("rust-tests", "python-tests", "python-fallback"),
"notebooks.yml": ("execute-notebooks", "interop-notebooks"),
+ "mmm-interop.yml": ("mmm-pymc-interop", "mmm-meridian-interop"),
"docs-tests.yml": ("doc-snippets", "sphinx-build", "docs-deps-py39-smoke"),
# release-build-check.yml is a single reusable-workflow caller job gated on
# ready-for-ci (it build-tests the PyPI release path on PRs); lock its guard too.
diff --git a/tests/test_t29_mmm_calibration_pymc_drift.py b/tests/test_t29_mmm_calibration_pymc_drift.py
new file mode 100644
index 00000000..73daa008
--- /dev/null
+++ b/tests/test_t29_mmm_calibration_pymc_drift.py
@@ -0,0 +1,250 @@
+"""Drift guards for docs/tutorials/29_mmm_calibration_pymc.ipynb.
+
+The notebook is CI-excluded from execution (it needs pymc-marketing, which
+conflicts with the Meridian tutorial's environment), so these tests protect its
+committed surface with diff-diff-only recomputation:
+
+- DGP constant + formula-line sync (the t26 pattern): the notebook's DGP must
+ equal the copy recomputed here, or stale committed outputs would stay green;
+- sampler-configuration needles: sampling cannot be weakened or de-seeded
+ silently;
+- integration-content presence: the framework fits and acceptance asserts cannot
+ be deleted while CI stays green;
+- rendered-number cross-check: the DiD-side numbers the notebook prints are
+ recomputed with diff-diff alone and matched against the committed outputs
+ (loose bands per the t27 normal-suite convention - never bit-exact floats).
+"""
+
+import re
+
+import numpy as np
+import pandas as pd
+import pytest
+
+from diff_diff import CallawaySantAnna
+from diff_diff.mmm import to_pymc_marketing_lift_test
+from tests._tutorial_drift import (
+ _read_notebook,
+ assert_quotes_in_rendered,
+ notebook_output_text,
+)
+
+NB = "docs/tutorials/29_mmm_calibration_pymc.ipynb"
+
+# ---- locked DGP constants (must equal the notebook's copy verbatim) ----
+SEED = 2026
+N_GEOS, N_WEEKS = 15, 104
+BETA_SEARCH, SEARCH_BASE, BOOST = 2.0, 100.0, 50.0
+FREEZE_WEEK = 52
+TREAT_COHORTS = {60: [0, 1, 2], 68: [3, 4, 5], 76: [6, 7, 8]}
+DEMAND_AMP, CHASE = 300.0, 0.5
+TV_BASE, TV_BETA, TV_ADSTOCK, TV_SAT = 200.0, 600.0, 0.5, 0.004
+BASE_SALES, GEO_SD, NOISE_SD = 3000.0, 100.0, 15.0
+
+DGP_CONSTANT_NEEDLES = [
+ "SEED = 2026",
+ "N_GEOS, N_WEEKS = 15, 104",
+ "BETA_SEARCH = 2.0",
+ "SEARCH_BASE = 100.0",
+ "BOOST = 50.0",
+ "FREEZE_WEEK = 52",
+ "TREAT_COHORTS = {60: [0, 1, 2], 68: [3, 4, 5], 76: [6, 7, 8]}",
+ "DEMAND_AMP = 300.0",
+ "CHASE = 0.5",
+ "TV_BASE, TV_BETA, TV_ADSTOCK, TV_SAT = 200.0, 600.0, 0.5, 0.004",
+ "BASE_SALES, GEO_SD, NOISE_SD = 3000.0, 100.0, 15.0",
+]
+
+# Mechanism lines: linear zero-carryover experiment channel, tv adstock, the
+# demand-chasing history, and the frozen baseline + constant boost.
+DGP_FORMULA_NEEDLES = [
+ "+ BETA_SEARCH * search",
+ "carry = tv[t] + TV_ADSTOCK * carry",
+ "search = SEARCH_BASE * (1 + CHASE * season[t])",
+ "search = SEARCH_BASE + (BOOST if boosted else 0.0)",
+ "DEMAND_AMP * season[t]",
+]
+
+SAMPLER_NEEDLES = [
+ "chains=4, draws=500, tune=3000",
+ "target_accept=0.999",
+ 'nuts_sampler_kwargs={"max_treedepth": 14}',
+ "random_seed=SEED",
+]
+
+INTEGRATION_NEEDLES = [
+ "add_lift_test_measurements(",
+ "mmm_plain = fit_mmm()",
+ "mmm_cal = fit_mmm(df_lift)",
+ "aggregation_result=simple,\n scale=G_TREATED,",
+ "assert div == 0",
+ "assert mtd == 0",
+ "assert max_rhat_model < 1.01",
+ "assert rhat_roi < 1.01",
+ "assert ess_bulk >= ess_floors[name]",
+ "assert ess_tail >= 300",
+ 'ess_floors = {"plain": 300, "calibrated": 800}',
+ "assert err_cal < err_plain",
+ "assert width_cal < 0.8 * width_plain",
+]
+
+
+def _notebook_source() -> str:
+ nb = _read_notebook(NB)
+ return "\n".join(
+ "".join(c["source"]) if isinstance(c["source"], list) else c["source"]
+ for c in nb["cells"]
+ if c["cell_type"] == "code"
+ )
+
+
+def _build_panel() -> pd.DataFrame:
+ rng = np.random.default_rng(SEED)
+ geo_start = {g: w for w, gg in TREAT_COHORTS.items() for g in gg}
+ season = np.sin(2 * np.pi * np.arange(N_WEEKS) / 52.0)
+ rows = []
+ for g in range(N_GEOS):
+ tv = TV_BASE * (1 + 0.1 * rng.random(N_WEEKS))
+ ad = np.zeros(N_WEEKS)
+ carry = 0.0
+ for t in range(N_WEEKS):
+ carry = tv[t] + TV_ADSTOCK * carry
+ ad[t] = carry
+ tvc = TV_BETA * (1 - np.exp(-TV_SAT * ad))
+ geo_level = rng.normal(0, GEO_SD)
+ for t in range(N_WEEKS):
+ if t < FREEZE_WEEK:
+ search = SEARCH_BASE * (1 + CHASE * season[t])
+ else:
+ boosted = g in geo_start and t >= geo_start[g]
+ search = SEARCH_BASE + (BOOST if boosted else 0.0)
+ sales = (
+ BASE_SALES
+ + geo_level
+ + DEMAND_AMP * season[t]
+ + BETA_SEARCH * search
+ + tvc[t]
+ + rng.normal(0, NOISE_SD)
+ )
+ rows.append((g, t, sales, search, tv[t], geo_start.get(g, 0)))
+ return pd.DataFrame(
+ rows, columns=["geo", "week", "sales", "search_spend", "tv_spend", "first_treat"]
+ )
+
+
+def _parse_number(output: str, label_regex: str) -> float:
+ match = re.search(label_regex, output)
+ assert match, f"notebook output lost the line matching {label_regex!r}"
+ return float(match.group(1).replace(",", ""))
+
+
+class TestNotebookSourceSync:
+ def test_dgp_constants_match(self):
+ src = _notebook_source()
+ missing = [n for n in DGP_CONSTANT_NEEDLES if n not in src]
+ assert not missing, f"notebook DGP constants drifted from the test copy: {missing}"
+
+ def test_dgp_formula_lines_present(self):
+ src = _notebook_source()
+ missing = [n for n in DGP_FORMULA_NEEDLES if n not in src]
+ assert not missing, f"notebook DGP mechanism lines drifted: {missing}"
+
+ def test_sampler_config_pinned(self):
+ src = _notebook_source()
+ missing = [n for n in SAMPLER_NEEDLES if n not in src]
+ assert not missing, f"notebook sampler configuration drifted: {missing}"
+
+ def test_integration_cells_present(self):
+ src = _notebook_source()
+ missing = [n for n in INTEGRATION_NEEDLES if n not in src]
+ assert not missing, f"notebook integration/acceptance cells missing: {missing}"
+
+ def test_no_kernelspec_metadata(self):
+ nb = _read_notebook(NB)
+ assert "kernelspec" not in nb["metadata"], (
+ "committed notebooks stay kernelspec-free; local execution names the "
+ "venv kernel explicitly"
+ )
+
+
+class TestRenderedOutputs:
+ def test_acceptance_checks_passed_in_committed_run(self):
+ assert_quotes_in_rendered(
+ NB, ["all acceptance checks passed", "both fits done"], surface="output"
+ )
+
+ def test_money_plot_output_present(self):
+ nb = _read_notebook(NB)
+ has_image = any(
+ "image/png" in (out.get("data") or {})
+ for c in nb["cells"]
+ if c["cell_type"] == "code"
+ for out in c.get("outputs", [])
+ )
+ assert has_image, "the committed money-plot image output is missing"
+
+
+@pytest.fixture(scope="module")
+def fitted():
+ panel = _build_panel()
+ res = CallawaySantAnna().fit(
+ panel, outcome="sales", unit="geo", time="week", first_treat="first_treat"
+ )
+ return panel, res
+
+
+class TestDidNumbersRecompute:
+ """Recompute the DiD/export side with diff-diff only; match printed numbers."""
+
+ def test_total_and_att_match_notebook(self, fitted):
+ panel, res = fitted
+ total = res.aggregate("total")
+ out = notebook_output_text(NB)
+ printed_att = _parse_number(out, r"per-geo-week ATT: ([\d.,-]+)")
+ printed_total = _parse_number(out, r"CS total incremental sales:\s+([\d.,-]+)")
+ np.testing.assert_allclose(printed_att, res.overall_att, rtol=2e-3)
+ np.testing.assert_allclose(printed_total, total.att[0], rtol=2e-3)
+ mass = int(((panel.first_treat > 0) & (panel.week >= panel.first_treat)).sum())
+ assert total.n[0] == mass
+
+ def test_lift_row_matches_notebook(self, fitted):
+ # Every lift-row FIELD is pinned individually (not just the ROI ratio,
+ # where the treated-geo multiplier would cancel): x/delta_x against the
+ # frozen-baseline national-weekly grain, delta_y/sigma against the
+ # aggregate('simple') container scaled by the treated-geo count.
+ panel, res = fitted
+ simple = res.aggregate("simple")
+ g_treated = int((panel.first_treat > 0).sum() / N_WEEKS)
+ assert g_treated == 9
+ df = to_pymc_marketing_lift_test(
+ channel="search",
+ x=SEARCH_BASE * N_GEOS,
+ delta_x=BOOST * g_treated,
+ aggregation_result=simple,
+ scale=g_treated,
+ )
+ assert float(df.x[0]) == SEARCH_BASE * N_GEOS == 1500.0
+ assert float(df.delta_x[0]) == BOOST * g_treated == 450.0
+ np.testing.assert_allclose(float(df.delta_y[0]), res.overall_att * g_treated, rtol=1e-12)
+ np.testing.assert_allclose(float(df.sigma[0]), res.overall_se * g_treated, rtol=1e-12)
+ out = notebook_output_text(NB)
+ printed_roi = _parse_number(out, r"DiD measurement:\s+([\d.,-]+)")
+ np.testing.assert_allclose(printed_roi, df.delta_y[0] / df.delta_x[0], atol=5e-3)
+
+ def test_notebook_states_linearity_scope_of_lift_row(self):
+ # The staggered-average -> single-transition construction is exact only
+ # under this simulation's linear homogeneous response; the notebook must
+ # keep saying so (round-4 review: an unscoped "generalizes" claim would
+ # transport a mixed-regime average onto one spend step of a saturating
+ # curve).
+ from tests._tutorial_drift import notebook_markdown
+
+ md = notebook_markdown(NB)
+ assert "linear and homogeneous" in md
+ assert "boost all test geos simultaneously" in md
+ # The nonlinear-case guidance must keep the estimator fit on the FULL
+ # panel (a post-only refit strips every cohort's pre-treatment baseline
+ # and leaves Callaway-Sant'Anna unidentified - round-5 review) and
+ # restrict only the post-fit aggregation.
+ assert "keep the FULL panel in the estimator fit" in md
+ assert "re-fit the DiD restricted" not in md
diff --git a/tests/test_t30_mmm_calibration_meridian_drift.py b/tests/test_t30_mmm_calibration_meridian_drift.py
new file mode 100644
index 00000000..304d22b7
--- /dev/null
+++ b/tests/test_t30_mmm_calibration_meridian_drift.py
@@ -0,0 +1,215 @@
+"""Drift guards for docs/tutorials/30_mmm_calibration_meridian.ipynb.
+
+Same protection pattern as the t29 twin (see its module docstring): the notebook
+is CI-excluded from execution (google-meridian conflicts with the PyMC tutorial's
+environment), so its committed surface is guarded here with diff-diff-only
+recomputation - DGP constant/formula sync, sampler-config needles,
+integration-content presence, and rendered-number cross-checks with loose bands.
+"""
+
+import re
+
+import numpy as np
+import pandas as pd
+import pytest
+
+from diff_diff import CallawaySantAnna
+from diff_diff.mmm import meridian_calibration_mask, to_meridian_roi_prior
+from tests._tutorial_drift import (
+ _read_notebook,
+ assert_quotes_in_rendered,
+ notebook_output_text,
+)
+
+NB = "docs/tutorials/30_mmm_calibration_meridian.ipynb"
+
+# ---- locked DGP constants (must equal the notebook's copy verbatim) ----
+SEED = 30301
+N_GEOS, N_WEEKS = 12, 104
+BETA_SEARCH, LAUNCH_SPEND = 2.5, 150.0
+LAUNCH_COHORTS = {30: [0, 1, 2], 44: [3, 4, 5], 58: [6, 7]}
+DEMAND_AMP = 80.0
+TV_BASE, TV_BETA, TV_ADSTOCK, TV_SAT = 250.0, 800.0, 0.4, 0.003
+BASE_SALES, GEO_SD, NOISE_SD, SHOCK_SD = 5000.0, 150.0, 20.0, 60.0
+
+DGP_CONSTANT_NEEDLES = [
+ "SEED = 30301",
+ "N_GEOS, N_WEEKS = 12, 104",
+ "BETA_SEARCH = 2.5",
+ "LAUNCH_SPEND = 150.0",
+ "LAUNCH_COHORTS = {30: [0, 1, 2], 44: [3, 4, 5], 58: [6, 7]}",
+ "DEMAND_AMP = 80.0",
+ "TV_BASE, TV_BETA, TV_ADSTOCK, TV_SAT = 250.0, 800.0, 0.4, 0.003",
+ "BASE_SALES, GEO_SD, NOISE_SD, SHOCK_SD = 5000.0, 150.0, 20.0, 60.0",
+]
+
+# Mechanism lines: linear zero-carryover experiment channel, launched-geos-only
+# spend, the demand ramp, and the common national shocks.
+DGP_FORMULA_NEEDLES = [
+ "+ BETA_SEARCH * search",
+ "carry = tv[t] + TV_ADSTOCK * carry",
+ "search = LAUNCH_SPEND if launched else 0.0",
+ "DEMAND_AMP * ramp[t] + shock[t]",
+ "shock = rng.normal(0, SHOCK_SD, N_WEEKS)",
+]
+
+SAMPLER_NEEDLES = [
+ "n_chains=4, n_adapt=500, n_burnin=500, n_keep=500, seed=SEED",
+]
+
+INTEGRATION_NEEDLES = [
+ "exec(code, namespace)",
+ "mer_default = fit_meridian(m_spec.ModelSpec())",
+ "mer_calibrated = fit_meridian(model_spec_calibrated)",
+ "to_meridian_roi_prior(aggregation_result=total, spend=total_spend)",
+ "assert np.isfinite(max_rhat_model) and max_rhat_model < 1.2",
+ "assert np.isfinite(rhat_roi) and rhat_roi < 1.05",
+ "assert err_calibrated < err_default",
+ "assert width_calibrated < 0.8 * width_default",
+ "assert total.n[0] == spend_positive",
+]
+
+
+def _notebook_source() -> str:
+ nb = _read_notebook(NB)
+ return "\n".join(
+ "".join(c["source"]) if isinstance(c["source"], list) else c["source"]
+ for c in nb["cells"]
+ if c["cell_type"] == "code"
+ )
+
+
+def _build_panel() -> pd.DataFrame:
+ rng = np.random.default_rng(SEED)
+ geo_start = {g: w for w, gg in LAUNCH_COHORTS.items() for g in gg}
+ ramp = np.clip((np.arange(N_WEEKS) - 26) / 52.0, 0, None)
+ shock = rng.normal(0, SHOCK_SD, N_WEEKS)
+ rows = []
+ for g in range(N_GEOS):
+ tv = TV_BASE * (1 + 0.1 * rng.random(N_WEEKS))
+ ad = np.zeros(N_WEEKS)
+ carry = 0.0
+ for t in range(N_WEEKS):
+ carry = tv[t] + TV_ADSTOCK * carry
+ ad[t] = carry
+ tvc = TV_BETA * (1 - np.exp(-TV_SAT * ad))
+ geo_level = rng.normal(0, GEO_SD)
+ for t in range(N_WEEKS):
+ launched = g in geo_start and t >= geo_start[g]
+ search = LAUNCH_SPEND if launched else 0.0
+ sales = (
+ BASE_SALES
+ + geo_level
+ + DEMAND_AMP * ramp[t]
+ + shock[t]
+ + BETA_SEARCH * search
+ + tvc[t]
+ + rng.normal(0, NOISE_SD)
+ )
+ rows.append((g, t, sales, search, tv[t], geo_start.get(g, 0)))
+ return pd.DataFrame(
+ rows, columns=["geo", "week", "sales", "search_spend", "tv_spend", "first_treat"]
+ )
+
+
+def _parse_number(output: str, label_regex: str) -> float:
+ match = re.search(label_regex, output)
+ assert match, f"notebook output lost the line matching {label_regex!r}"
+ return float(match.group(1).replace(",", ""))
+
+
+class TestNotebookSourceSync:
+ def test_dgp_constants_match(self):
+ src = _notebook_source()
+ missing = [n for n in DGP_CONSTANT_NEEDLES if n not in src]
+ assert not missing, f"notebook DGP constants drifted from the test copy: {missing}"
+
+ def test_dgp_formula_lines_present(self):
+ src = _notebook_source()
+ missing = [n for n in DGP_FORMULA_NEEDLES if n not in src]
+ assert not missing, f"notebook DGP mechanism lines drifted: {missing}"
+
+ def test_sampler_config_pinned(self):
+ src = _notebook_source()
+ missing = [n for n in SAMPLER_NEEDLES if n not in src]
+ assert not missing, f"notebook sampler configuration drifted: {missing}"
+
+ def test_integration_cells_present(self):
+ src = _notebook_source()
+ missing = [n for n in INTEGRATION_NEEDLES if n not in src]
+ assert not missing, f"notebook integration/acceptance cells missing: {missing}"
+
+ def test_no_kernelspec_metadata(self):
+ nb = _read_notebook(NB)
+ assert "kernelspec" not in nb["metadata"], (
+ "committed notebooks stay kernelspec-free; local execution names the "
+ "venv kernel explicitly"
+ )
+
+
+class TestRenderedOutputs:
+ def test_acceptance_checks_passed_in_committed_run(self):
+ assert_quotes_in_rendered(
+ NB, ["all acceptance checks passed", "both fits done"], surface="output"
+ )
+
+ def test_money_plot_output_present(self):
+ nb = _read_notebook(NB)
+ has_image = any(
+ "image/png" in (out.get("data") or {})
+ for c in nb["cells"]
+ if c["cell_type"] == "code"
+ for out in c.get("outputs", [])
+ )
+ assert has_image, "the committed money-plot image output is missing"
+
+ def test_generated_snippet_shown(self):
+ assert_quotes_in_rendered(
+ NB,
+ ["Generated by diff_diff.mmm.to_meridian_roi_prior"],
+ surface="output",
+ )
+
+
+@pytest.fixture(scope="module")
+def fitted():
+ panel = _build_panel()
+ res = CallawaySantAnna().fit(
+ panel, outcome="sales", unit="geo", time="week", first_treat="first_treat"
+ )
+ return panel, res
+
+
+class TestDidNumbersRecompute:
+ """Recompute the DiD/export side with diff-diff only; match printed numbers."""
+
+ def test_total_and_prior_match_notebook(self, fitted):
+ panel, res = fitted
+ total = res.aggregate("total")
+ total_spend = float(panel.search_spend.sum())
+ prior = to_meridian_roi_prior(aggregation_result=total, spend=total_spend)
+ out = notebook_output_text(NB)
+ printed_total = _parse_number(out, r"CS total incremental sales:\s+([\d.,-]+)")
+ printed_roi = _parse_number(out, r"experiment ROI: ([\d.,-]+)")
+ np.testing.assert_allclose(printed_total, total.att[0], rtol=2e-3)
+ np.testing.assert_allclose(printed_roi, prior.roi_mean, atol=5e-3)
+ # launch-design alignment: the total's mass is the spend-positive count
+ assert total.n[0] == int((panel.search_spend > 0).sum())
+
+ def test_mask_window_matches_notebook(self, fitted):
+ panel, _ = fitted
+ geo_start = {g: w for w, gg in LAUNCH_COHORTS.items() for g in gg}
+ first_launch = min(geo_start.values())
+ times = [str(d.date()) for d in pd.date_range("2024-01-01", periods=N_WEEKS, freq="W-MON")]
+ mask = meridian_calibration_mask(
+ media_times=times,
+ media_channels=["search", "tv"],
+ channel="search",
+ window=(times[first_launch], times[-1]),
+ )
+ out = notebook_output_text(NB)
+ printed_true = _parse_number(out, r"search weeks in window: (\d+)")
+ assert int(printed_true) == int(mask[:, 0].sum())
+ # every spend-positive week is inside the window (the time-axis
+ # containment the TIME x CHANNEL mask expresses)
+ assert panel.loc[panel.search_spend > 0, "week"].min() >= first_launch