-
Notifications
You must be signed in to change notification settings - Fork 51
Expand file tree
/
Copy pathmmm.py
More file actions
1620 lines (1533 loc) · 82.6 KB
/
Copy pathmmm.py
File metadata and controls
1620 lines (1533 loc) · 82.6 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
"""Assemble Marketing Mix Model (MMM) calibration inputs from experiment results.
MMM practitioners calibrate their models against experimental evidence. The two
dominant Python MMM frameworks consume that evidence in different shapes:
- **PyMC-Marketing** (and prophetverse, which mirrors the same schema) ingests a
*lift-test* DataFrame with columns ``channel``, optional model dims (e.g. ``geo``),
``x`` (baseline channel spend), ``delta_x`` (spend change during the test),
``delta_y`` (measured incremental outcome), and ``sigma`` (the experiment's standard
error), which it scores against the model's saturation curve via
``MMM.add_lift_test_measurements``. See
https://www.pymc-marketing.io/en/stable/notebooks/mmm/mmm_lift_test.html.
- **Google Meridian** has no experiment-ingestion API; calibration means setting a
lognormal prior per channel - ``roi_m`` for the return of a channel's full spend
(a full-holdout/zero-spend estimand), ``mroi_m`` for a marginal return. Google's
documented workflow maps the experiment's ROI point estimate to the prior mean and
its standard error to the prior standard deviation, converted to lognormal
``(mu, sigma)`` via the closed form used by
``meridian.model.prior_distribution.lognormal_dist_from_mean_std``:
``mu = ln(m) - ln((s/m)^2 + 1)/2``, ``sigma = sqrt(ln((s/m)^2 + 1))``. See
https://developers.google.com/meridian/docs/advanced-modeling/set-custom-priors-past-experiments.
**Design: explicit numbers in, or the pinned aggregation contract in - validated out.**
Reconciling an experiment's estimate to a calibration input requires context diff-diff
cannot always see - the target MMM's row granularity (per-geo vs national), its time
window, and the outcome's scale (additive levels vs a log/rate/share). The default
route therefore stays fully explicit: the CALLER supplies the already-scoped
incremental outcome and its standard error (the numbers they read off a fitted
result's ``summary()``, aggregated to the population and window their MMM row
represents). Alternatively, both exporters accept ``aggregation_result=`` - the pinned
:class:`~diff_diff.aggregation.AggregationResult` container returned by
``results.aggregate('simple')`` / ``results.aggregate('group')`` (with ``scale=``,
deriving ``effect = att * scale`` / ``se_out = se * scale`` per row) or by
``results.aggregate('total')`` (whose single row is the estimator-owned total
incremental outcome and takes NO scale - for overall-total exports this route
supersedes ``scale="auto"``: its finite-masked support eliminates the documented
Imputation/TwoStage raw-support overcount for that use; the total is meaningful
only for outcomes additive in levels, the same caveat every scale route carries).
The derivation is fail-closed: raw results objects are rejected (this module never
calls ``aggregate()`` itself), and ``scale="auto"`` (reading the row count off the
container) is honored only for the audited producers whose ``n`` matches the ATT's
averaged support on unweighted, fully identified fits. diff-diff does only what it can verify:
assemble the exact target schema, enforce the sign/positivity/monotonicity guards each
consumer requires, convert to the lognormal parameterization (with parity to Google's
closed form), pool multiple experiments, and emit ready-to-paste snippets. It is pure
numpy/pandas and never imports an MMM package.
"""
from __future__ import annotations
import ast
import math
import warnings
from dataclasses import dataclass, field
from typing import (
Any,
Dict,
List,
Literal,
Mapping,
Optional,
Sequence,
Tuple,
Union,
)
import numpy as np
import pandas as pd
# Runtime import (no cycle: aggregation.py depends only on results_base), so the
# exporters' public annotations resolve under typing.get_type_hints().
from diff_diff.aggregation import AggregationResult
from diff_diff.results_base import BaseResults
__all__ = [
"MeridianROIPrior",
"meridian_calibration_mask",
"to_meridian_roi_prior",
"to_pymc_marketing_lift_test",
]
_WRONG_SIGN_POLICIES = ("raise", "drop", "keep")
_LIFT_TEST_RESERVED = frozenset({"channel", "x", "delta_x", "delta_y", "sigma"})
# Container-mode routing (audited against each producer's aggregation source):
# scale="auto" reads per-row n off the container ONLY for these producers, whose n
# matches the treated observations (unit-periods) the ATT averages over on
# unweighted, fully identified fits, at the 'simple' and 'group' levels
# ('total' containers arrive already scaled and take neither scale nor "auto").
# Routing keys on AggregationResult.estimator - NEVER on n_kind
# alone: CallawaySantAnna repeated-cross-section fits report n_kind="obs" with
# n = treated+control observations, so n_kind is only a drift sanity guard here.
_SCALE_AUTO_ESTIMATORS = frozenset({"ImputationDiD", "TwoStageDiD"})
_CONTAINER_LEVELS = ("simple", "group", "total")
# The audited aggregate('total') producers: only their total containers are
# admitted (the "already scaled" claim is provenance-gated like "auto" -
# StackedDiD is staged out and unknown provenance could be any hand-built
# container). A new adopter extends this set in the PR that ships its total.
_TOTAL_ESTIMATORS = frozenset(
{"CallawaySantAnna", "DMLDiD", "EfficientDiD", "ImputationDiD", "TwoStageDiD"}
)
# Why scale="auto" is refused, per audited non-allowlisted producer. The generic
# closing prescription in the error is an EXAMPLE for the unweighted additive case,
# never a universal formula.
_SCALE_HINTS = {
"CallawaySantAnna": (
"CallawaySantAnna's n is not a treated unit-period count: the 'simple' "
"container reports treated+control units (treated+control OBSERVATIONS on "
"repeated-cross-section fits, where n_kind='obs' still does not mean "
"treated unit-periods), and the 'group' container reports contributing "
"(g, t) cells (n_kind='cells'). Where supported, use "
"results.aggregate('total') instead - its single row is the "
"estimator-owned total incremental outcome and needs no scale; on fits "
"it does not support it raises with the reason (repeated-cross-section, "
"declared survey_design, cluster-mass fits with incomplete treated "
"cells, or pre-upgrade results) - there, pass a numeric scale (a "
"caller-defined estimand, not the estimator-owned complete-case total)."
),
"DMLDiD": (
"DMLDiD's n is not a treated unit-period count: the 'simple' container "
"reports treated+control units and the 'group' container reports "
"contributing (g, t) cells (n_kind='cells'). Use "
"results.aggregate('total') instead - its single row is the "
"estimator-owned total incremental outcome over the per-cell "
"complete-case treated units and needs no scale. "
"PANEL fits support the total route (bootstrapped fits replay) except "
"declared-survey_design fits and bare-cluster fits whose kept cells "
"have incomplete treated support (both raise with the reason - the "
"same CS exclusions); repeated-cross-section fits "
"(panel=False) fail aggregate('total') closed on EVERY adopter "
"(CallawaySantAnna's RCS fits included) - estimator-owned totals need "
"per-unit tracking; on RCS fits pass a caller-defined numeric scale, "
"or use a panel fit when genuine longitudinal data exists. If you need a different "
"estimand (e.g. full-cohort exposure on a panel with missing cells), "
"pass a numeric scale."
),
"EfficientDiD": (
"EfficientDiD's n counts disjoint treated+control units ('simple') or "
"contributing (g, t) cells ('group'), not treated unit-periods. Where "
"supported, use results.aggregate('total') instead - its single row is "
"the estimator-owned total incremental outcome and needs no scale; its "
"only unsupported routing is a fit declaring a survey_design, which "
"raises with the reason - there, pass a numeric scale (a caller-defined "
"estimand, not the estimator-owned complete-case total)."
),
"StackedDiD": (
"StackedDiD's n is a deduplicated distinct-treated-unit count - units, not "
"treated unit-periods, so no per-observation mass exists - and under "
"weighting='population'/'sample_share' the ATT is a weighted estimand for "
"which raw treated exposure is not the right multiplier; derive the scoped "
"total per your weighting choice."
),
}
# 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; 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) 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.
_MERIDIAN_SINGLE_CHANNEL_TEMPLATE = """\
# Generated by diff_diff.mmm.to_meridian_roi_prior (template pinned to Meridian 1.7.0)
# TensorFlow-substrate snippet; for JAX-backed Meridian use
# `import tensorflow_probability.substrates.jax as tfp` instead.
# SINGLE-CHANNEL MODEL ONLY: a scalar {param} prior broadcasts to every media channel.
# For a multi-channel model, regenerate with to_code(channel=..., media_channels=[...]).
import tensorflow_probability as tfp
from meridian.model import prior_distribution, spec
{mask_prelude}roi_prior = tfp.distributions.LogNormal({mu!r}, {sigma!r}, name="{param}")
prior = prior_distribution.PriorDistribution({param}=roi_prior)
model_spec = spec.ModelSpec(
prior=prior,
media_prior_type="{prior_type}",
# {window_note}
roi_calibration_period={calibration_period},
)
"""
_MERIDIAN_MULTI_CHANNEL_TEMPLATE = """\
# Generated by diff_diff.mmm.to_meridian_roi_prior (template pinned to Meridian 1.7.0)
# TensorFlow-substrate snippet; for JAX-backed Meridian use
# `import tensorflow_probability.substrates.jax as tfp` instead.
# Channel order MUST match your Meridian InputData media channel order exactly:
# {channels}
# {channel!r} carries the experiment-informed prior; the other channels keep
# Meridian's default {param} prior LogNormal({default_mu}, {default_sigma}).
import tensorflow_probability as tfp
from meridian.model import prior_distribution, spec
{mask_prelude}mu = {mu_vector!r}
sigma = {sigma_vector!r}
roi_prior = tfp.distributions.LogNormal(mu, sigma, name="{param}")
prior = prior_distribution.PriorDistribution({param}=roi_prior)
model_spec = spec.ModelSpec(
prior=prior,
media_prior_type="{prior_type}",
# {window_note}
roi_calibration_period={calibration_period},
)
"""
def _is_sequence(value: Any) -> bool:
"""True for list-like containers (not strings, not mappings)."""
return isinstance(value, (list, tuple, np.ndarray, pd.Series))
def _broadcast(name: str, value: Any, n: int) -> List[Any]:
"""Broadcast a scalar to n rows, or validate a sequence's length."""
if _is_sequence(value):
values = list(value)
if len(values) != n:
raise ValueError(
f"{name} has length {len(values)} but {n} experiment(s) were given; "
f"pass one value per experiment or a single scalar"
)
return values
return [value] * n
def _seq_len(*values: Any) -> int:
"""Number of experiments implied by the first sequence argument (else 1)."""
for v in values:
if _is_sequence(v):
n = len(list(v))
if n == 0:
raise ValueError("empty sequence given; provide at least one experiment")
return n
return 1
def _finite_positive(name: str, value: Any, index: int) -> float:
v = float(value)
if not math.isfinite(v) or v <= 0:
raise ValueError(f"{name} must be finite and > 0; got {value!r} for experiment[{index}]")
return v
def _normalize_dims(
dims: Optional[Union[Mapping[str, str], Sequence[Mapping[str, str]]]], n: int
) -> Tuple[List[str], List[Mapping[str, str]]]:
"""Broadcast dims mappings and pin a single shared key set / column order."""
if dims is None:
return [], []
if isinstance(dims, Mapping):
rows: List[Mapping[str, str]] = [dims] * n
else:
rows = list(_broadcast("dims", dims, n))
for i, row in enumerate(rows):
if not isinstance(row, Mapping):
raise TypeError(
f"dims must be a mapping or a sequence of mappings (e.g. "
f"{{'geo': 'US-CA'}}); dims[{i}] is {type(row).__name__}"
)
dim_cols = list(rows[0].keys())
key_set = set(dim_cols)
reserved = key_set & _LIFT_TEST_RESERVED
if reserved:
raise ValueError(
f"dims keys {sorted(reserved)} collide with reserved lift-test columns "
f"{sorted(_LIFT_TEST_RESERVED)}; rename the model dimension"
)
for i, row in enumerate(rows):
if set(row.keys()) != key_set:
raise ValueError(
f"dims mappings must share one key set across rows; dims[{i}] has keys "
f"{sorted(row.keys())}, expected {sorted(key_set)}"
)
return dim_cols, rows
def _extract_aggregation_rows(
aggregation_result: Any,
scale: Optional[Union[float, Sequence[float], str]],
*,
effect_name: str,
se_name: str,
) -> Tuple[List[float], List[float], List[Any]]:
"""Derive per-row ``(effect, se)`` from a pinned ``AggregationResult``.
Every per-row value (``label``, ``target``, ``att``, ``se``, ``n``) is read
from ONE ``to_dataframe()`` call, so rows arrive in the order ``summary()``
prints (sorted by label when sortable, producer order otherwise) - the
alignment order for every per-row sequence kwarg, ``scale`` included.
Fail-closed throughout: anything the container cannot verify raises with the
remedy inline.
"""
if not isinstance(aggregation_result, AggregationResult):
if isinstance(aggregation_result, BaseResults):
raise TypeError(
f"aggregation_result must be an AggregationResult - the container "
f"returned by res.aggregate('simple'), res.aggregate('group'), or "
f"res.aggregate('total') on "
f"estimators that produce one; got "
f"{type(aggregation_result).__name__}. EventStudyResults and "
f"estimators whose aggregate() does not return that container "
f"(WooldridgeDiDResults returns the results object itself; "
f"HeterogeneousAdoptionDiD event-study results support only "
f"'event_study') have no container-mode route - pass the explicit "
f"{effect_name}/{se_name} arguments instead"
)
raise TypeError(
f"aggregation_result must be an AggregationResult (the return value of "
f"res.aggregate('simple'), res.aggregate('group'), or "
f"res.aggregate('total')); got "
f"{type(aggregation_result).__name__}. Only AggregationResult is "
f"supported."
)
level = aggregation_result.level
if level not in _CONTAINER_LEVELS:
raise ValueError(
f"aggregation_result.level must be 'simple', 'group', or 'total'; got "
f"{level!r}. "
f"Other levels ('calendar', 'dose', estimator-specific extras) have no "
f"defined MMM-experiment mapping - re-aggregate at a supported level."
)
frame = aggregation_result.to_dataframe()
if len(frame) == 0:
raise ValueError(
"aggregation_result has no rows (the aggregation selected no cells); "
"nothing to export"
)
if level == "total":
# A total container's row is ALREADY the total incremental outcome
# (the estimator applied its own finite-masked aggregation mass), so
# the simple/group target+scale machinery below does not apply: the
# row relays as-is and any scale would double-count. Because "already
# scaled" is a PRODUCER claim, admission is provenance-gated exactly
# like scale="auto": only the audited total adopters are trusted
# (StackedDiD is staged out - its total estimand is undefined under
# weighting= variants - and unknown/missing provenance could be any
# hand-built container whose att is NOT a total).
estimator = aggregation_result.estimator
if estimator not in _TOTAL_ESTIMATORS:
raise ValueError(
f"a level='total' container is accepted only from the audited "
f"total adopters {sorted(_TOTAL_ESTIMATORS)}; got estimator "
f"provenance {estimator!r}. Their aggregate('total') is the "
f"only producer whose single row is verifiably the "
f"estimator-owned total (StackedDiD totals are staged out - "
f"see DEFERRED.md); for other sources pass the explicit "
f"{effect_name}/{se_name} numbers, or a 'simple'/'group' "
f"container with a numeric scale."
)
if len(frame) != 1:
raise ValueError(
f"aggregation_result has level='total' but {len(frame)} rows; "
f"producers emit a single 'total' row, so this container is "
f"out of contract - re-aggregate with results.aggregate('total')"
)
total_target = frame["target"].tolist()[0]
if total_target != "total":
raise ValueError(
f"a 'total' container's single row must carry target='total'; "
f"got {total_target!r}. The container is rejected whole - "
f"re-aggregate with results.aggregate('total')"
)
total_label = frame["label"].tolist()[0]
n_kind = aggregation_result.n_kind
weight_arr = aggregation_result.weight
weight_ok = weight_arr is not None and len(weight_arr) == 1 and float(weight_arr[0]) == 1.0
if total_label != "total" or n_kind != "obs" or not weight_ok:
raise ValueError(
f"a 'total' container must carry the producer contract "
f"label='total', n_kind='obs', weight=[1.0]; got "
f"label={total_label!r}, n_kind={n_kind!r}, "
f"weight={None if weight_arr is None else list(weight_arr)!r}. "
f"The container schema has drifted from the audited contract - "
f"re-aggregate with results.aggregate('total')"
)
labels = frame["label"].tolist()
atts = [float(v) for v in frame["att"].tolist()]
ses = [float(v) for v in frame["se"].tolist()]
for i, (att_i, se_i) in enumerate(zip(atts, ses)):
if not (math.isfinite(att_i) and math.isfinite(se_i) and se_i > 0):
raise ValueError(
f"aggregation_result row [{i}] (label {labels[i]!r}) has "
f"att={att_i!r}, se={se_i!r}; the fit carries no usable point "
f"estimate/SE for this row, so it cannot calibrate an MMM"
)
if scale is not None:
raise ValueError(
"scale is not accepted with a level='total' container: this "
"container's rows are already totals; scale would double-count "
"- pass the container alone (aggregate('total') already "
"applied the estimator's own finite-masked aggregation mass)"
)
return atts, ses, labels
if level == "simple":
target_list = list(frame["target"].tolist())
offending = sorted({t for t in target_list if target_list.count(t) > 1})
if offending:
raise ValueError(
f"aggregation_result has level='simple' but multiple rows for "
f"target(s) {offending!r}; producers emit a single 'overall' row "
f"per target for 'simple', so this container is out of contract - "
f"re-aggregate, or use level='group' for per-cohort rows"
)
bad_targets = [(i, t) for i, t in enumerate(frame["target"].tolist()) if t != "att"]
if bad_targets:
idxs = [i for i, _ in bad_targets]
targets = sorted({t for _, t in bad_targets})
raise ValueError(
f"aggregation_result row(s) {idxs} have target {targets!r} where 'att' "
f"is required; only ATT rows map to an MMM incremental outcome "
f"(ContinuousDiD's 'acrt' dose derivative, HeterogeneousAdoptionDiD's "
f"WAS estimands, and ChaisemartinDHaultfoeuille's estimand relays are "
f"not per-unit-period ATTs). The container is rejected whole - no "
f"silent row filtering."
)
labels = frame["label"].tolist()
atts = [float(v) for v in frame["att"].tolist()]
ses = [float(v) for v in frame["se"].tolist()]
ns = [float(v) for v in frame["n"].tolist()]
for i, (att_i, se_i) in enumerate(zip(atts, ses)):
if not (math.isfinite(att_i) and math.isfinite(se_i) and se_i > 0):
raise ValueError(
f"aggregation_result row [{i}] (label {labels[i]!r}) has "
f"att={att_i!r}, se={se_i!r}; the fit carries no usable point "
f"estimate/SE for this row, so it cannot calibrate an MMM"
)
n_rows = len(frame)
if isinstance(scale, str):
if scale != "auto":
raise ValueError(
f"scale must be a number, a sequence of numbers, or the string "
f"'auto'; got {scale!r}"
)
estimator = aggregation_result.estimator
if estimator not in _SCALE_AUTO_ESTIMATORS:
if estimator is not None and estimator in _SCALE_HINTS:
hint = _SCALE_HINTS[estimator]
else:
hint = (
f"estimator provenance is {estimator!r}, whose n semantics "
f"diff-diff has not audited for this derivation."
)
raise ValueError(
f"scale='auto' is not available for this container: {hint} Pass a "
f"numeric scale=<the factor converting this row's per-observation "
f"ATT to your MMM row's total, e.g. treated units x treated "
f"periods for an unweighted additive fit>, in to_dataframe() "
f"order; auto-derivation is audited only for ImputationDiD and "
f"TwoStageDiD, whose n is the treated observations the ATT "
f"averages over."
)
n_kind = aggregation_result.n_kind
if n_kind != "obs":
raise ValueError(
f"aggregation_result from {estimator} reports n_kind={n_kind!r}, "
f"expected 'obs'; the container schema has drifted from the "
f"audited contract - pass a numeric scale explicitly"
)
for i, n_i in enumerate(ns):
if not (math.isfinite(n_i) and n_i > 0):
raise ValueError(
f"aggregation_result row [{i}] (label {labels[i]!r}) has "
f"n={n_i!r}; cannot auto-derive scale - pass a numeric scale"
)
scales = ns
elif scale is not None:
# bool is an int subclass, so float(True) == 1.0 would silently scale
# by one - a plausible typo for scale="auto" - and must fail closed.
scale_values = _broadcast("scale", scale, n_rows)
if isinstance(scale, (bool, np.bool_)) or any(
isinstance(v, (bool, np.bool_)) for v in scale_values
):
raise ValueError(
"scale must be a number, a sequence of numbers, or the string "
"'auto'; got a boolean (did you mean scale='auto'?)"
)
scales = [_finite_positive("scale", v, i) for i, v in enumerate(scale_values)]
else:
raise ValueError(
f"scale is required with aggregation_result: pass a numeric "
f"scale=<the factor converting each row's per-observation ATT to your "
f"MMM row's total; e.g. treated units x treated periods for an "
f"unweighted additive fit> (scalar or one value per row in "
f"to_dataframe() order), or scale='auto' (ImputationDiD/TwoStageDiD "
f"fits only - see the {effect_name} docstring for the assumptions "
f"'auto' acknowledges). For an overall total with no scale at all, "
f"pass results.aggregate('total') where the estimator supports it - "
f"it supersedes scale='auto' for total-report use"
)
effects = [att_i * s_i for att_i, s_i in zip(atts, scales)]
ses_out = [se_i * s_i for se_i, s_i in zip(ses, scales)]
return effects, ses_out, labels
def to_pymc_marketing_lift_test(
*,
channel: Union[str, Sequence[str]],
x: Union[float, Sequence[float]],
delta_x: Union[float, Sequence[float]],
delta_y: Optional[Union[float, Sequence[float]]] = None,
sigma: Optional[Union[float, Sequence[float]]] = None,
aggregation_result: Optional[AggregationResult] = None,
scale: Optional[Union[float, Sequence[float], Literal["auto"]]] = None,
dims: Optional[Union[Mapping[str, str], Sequence[Mapping[str, str]]]] = None,
on_wrong_sign: str = "raise",
) -> pd.DataFrame:
"""Assemble a PyMC-Marketing lift-test DataFrame from scoped experiment results.
Produces one row per experiment with columns ``channel``, any ``dims`` columns,
``x``, ``delta_x``, ``delta_y``, ``sigma`` - the exact schema consumed by
``pymc_marketing.mmm.MMM.add_lift_test_measurements`` (prophetverse's lift-test
API accepts the same shape). All values stay in original data units
(PyMC-Marketing rescales internally).
**Two input routes.** Either the caller supplies the scoped effect explicitly -
``delta_y`` and ``sigma``, the measured incremental outcome and its standard
error already aggregated to the population and time window ONE target-MMM row
represents - or passes ``aggregation_result=`` - the pinned container returned by
``results.aggregate('simple')`` / ``results.aggregate('group')`` (together with
``scale=``, this function deriving ``delta_y = att * scale`` and
``sigma = se * scale`` per container row) or by ``results.aggregate('total')``
(a single already-scaled total row; NO scale accepted). Rescaling is performed ONLY under that
explicit contract - reconciliation context the container cannot carry (the MMM's
row granularity, the outcome's scale) remains the caller's acknowledgement, see
``scale``. Either way, PyMC-Marketing scores one row's ``delta_y`` against
``saturation(x + delta_x) - saturation(x)``, so ``x``, ``delta_x``, ``delta_y``,
and ``sigma`` must all describe the SAME observation (same channel, same
population, same period span, additive-level outcome).
Parameters
----------
channel : str or sequence of str
MMM channel name(s); must match the target model's ``channel_columns``.
x : float or sequence of float
Baseline channel spend for the row's observation, in original spend units.
delta_x : float or sequence of float
Spend change during the experiment (nonzero; negative for go-dark/holdout
tests, with ``x + delta_x >= 0``), in the same units as ``x``.
delta_y : float or sequence of float, optional
Measured incremental outcome for the SAME observation as ``x``/``delta_x``,
in original outcome units. Finite; must be nonzero and share ``delta_x``'s
sign (see ``on_wrong_sign``). Required together with ``sigma`` unless
``aggregation_result`` is given (the two routes are mutually exclusive).
sigma : float or sequence of float, optional
Standard error of ``delta_y`` (finite, positive). Required together with
``delta_y`` unless ``aggregation_result`` is given.
aggregation_result : AggregationResult, optional
The pinned container returned by ``results.aggregate('simple')`` (one
experiment row), ``results.aggregate('group')`` (one row per cohort), or
``results.aggregate('total')`` (one already-scaled total row).
Mutually exclusive with ``delta_y``/``sigma``; requires ``scale`` for
'simple'/'group' containers and FORBIDS it for 'total'. Rows are
consumed in ``aggregation_result.to_dataframe()`` order - the order
``summary()`` prints - and every per-row sequence kwarg (``scale``, ``x``,
``delta_x``, ``dims``, ...) aligns to that order. Raw results objects and
``EventStudyResults`` are rejected; this function never calls
``aggregate()`` itself. On a bootstrapped fit, ``aggregate('simple')``
and, where supported, ``aggregate('total')`` relay the stored percentile
SE, which is used as-is. Group-container
cautions: (1) cohort rows come from ONE fit sharing controls and windows,
yet each emitted lift row is scored by PyMC-Marketing as an independent
observation with only its marginal ``sigma`` - the omitted cross-cohort
covariance can misstate the joint evidence (overstating it when the net
covariance is positive, as shared controls typically induce), and its
direction cannot be determined from marginal SEs alone; (2) a scalar
``x`` or
``delta_x`` replicates to every cohort row, asserting the channel's FULL
baseline spend / spend change produced only that cohort's ``delta_y`` -
pass per-cohort values; (3) the emitted frame carries NO cohort-label
column, so pass per-row ``dims`` to keep multi-cohort rows
distinguishable.
scale : float, sequence of float, or "auto", optional
Converts each container row's per-observation ATT to the row's total
incremental outcome: ``delta_y = att * scale``, ``sigma = se * scale``.
Required with a 'simple'/'group' ``aggregation_result``; forbidden with a
'total' container, whose row is already the estimator-owned total
(finite, positive; scalar or one value
per row in ``to_dataframe()`` order; e.g. treated units x treated periods
for an unweighted additive fit - an example, not a universal formula).
``scale="auto"`` derives ``scale`` from the container's per-row ``n`` and is
accepted ONLY for ImputationDiD and TwoStageDiD fits, whose ``n`` is the
treated observations the ATT averages over. Passing ``"auto"`` acknowledges
three assumptions the container cannot verify: the outcome is in additive
levels (not log/rate/share); the fit is unweighted (on survey-weighted fits
``att`` is a weighted average while ``n`` is a raw count); and every
treated observation's effect is identified (ImputationDiD averages finite
tau-hat only while ``n`` stays the raw count at both levels; TwoStageDiD's
'simple' ATT support excludes rows with non-finite first-stage residuals
while ``n`` reports the pre-filter count, though its 'group' counts are
post-filter - the affected fits warn at fit time on the degenerate
branch, and ``"auto"`` there overcounts). An
``att * scale`` (or ``se * scale``) that overflows the float range
surfaces via this function's ordinary finiteness errors - ``sigma`` is
validated before ``delta_y``, so a row where both overflow reports
``sigma``.
dims : mapping or sequence of mappings, optional
Extra model-dimension columns, e.g. ``{"geo": "US-CA"}`` for a geo-level MMM
built with ``MMM(dims=("geo",))``. Values must match the target model's
coordinate values exactly and identify the population ``delta_y`` describes
(do not label a multi-geo average with one specific geo). All rows must share
one key set; column order follows the first mapping.
on_wrong_sign : {"raise", "drop", "keep"}, default "raise"
Policy for rows PyMC-Marketing cannot use: ``sign(delta_y)`` contradicting
``sign(delta_x)`` (rejected upstream with ``NonMonotonicError``), or
``delta_y == 0`` (degenerate for its strictly-positive Gamma lift likelihood,
which its own monotonicity check does not catch). ``"raise"`` (default) errors
with guidance; ``"drop"`` warns and removes such rows (raises if that would
empty the frame); ``"keep"`` warns and emits them anyway (for non-PyMC
consumers - the frame is not valid PyMC-Marketing input).
Returns
-------
pd.DataFrame
Columns ``[channel, *dims, x, delta_x, delta_y, sigma]``, one row per
experiment.
Raises
------
ValueError
On invalid inputs (non-positive ``sigma``, non-finite values, ``delta_x == 0``,
``x + delta_x < 0``, broadcasting length mismatches, empty inputs,
heterogeneous ``dims`` key sets) or wrong-signed/zero rows under the default
policy.
"""
if on_wrong_sign not in _WRONG_SIGN_POLICIES:
raise ValueError(
f"on_wrong_sign must be one of {_WRONG_SIGN_POLICIES}; got {on_wrong_sign!r}"
)
if aggregation_result is not None:
if delta_y is not None or sigma is not None:
raise ValueError(
"pass either aggregation_result= or delta_y=/sigma=, not both; "
"aggregation mode derives delta_y and sigma from the container "
"(delta_y = att * scale, sigma = se * scale for 'simple'/'group' "
"containers; a 'total' container's row relays as-is)"
)
delta_y, sigma, _ = _extract_aggregation_rows(
aggregation_result, scale, effect_name="delta_y", se_name="sigma"
)
n = len(delta_y)
else:
if scale is not None:
raise ValueError(
"scale only applies with aggregation_result=; without a container, "
"pass the already-scaled delta_y and sigma directly"
)
if delta_y is None or sigma is None:
missing = (
"delta_y and sigma are"
if delta_y is None and sigma is None
else ("delta_y is" if delta_y is None else "sigma is")
)
raise ValueError(
f"{missing} required when aggregation_result is not given; pass "
f"both delta_y and sigma, pass "
f"aggregation_result=res.aggregate('simple'|'group') with scale= "
f"to derive them from a fitted result, or pass "
f"aggregation_result=res.aggregate('total') (no scale) where the "
f"estimator supports it"
)
n = _seq_len(channel, x, delta_x, delta_y, sigma, dims)
channels = _broadcast("channel", channel, n)
xs = _broadcast("x", x, n)
delta_xs = _broadcast("delta_x", delta_x, n)
delta_ys = _broadcast("delta_y", delta_y, n)
sigmas = _broadcast("sigma", sigma, n)
dim_cols, dim_rows = _normalize_dims(dims, n)
rows: List[Dict[str, Any]] = []
wrong_sign_rows: List[int] = []
zero_lift_rows: List[int] = []
for i in range(n):
x_i = float(xs[i])
dx_i = float(delta_xs[i])
dy_i = float(delta_ys[i])
sig_i = _finite_positive("sigma", sigmas[i], i)
if not math.isfinite(x_i) or x_i < 0:
raise ValueError(f"x must be finite and >= 0; got {xs[i]!r} for experiment[{i}]")
if not math.isfinite(dx_i) or dx_i == 0:
raise ValueError(
f"delta_x must be finite and nonzero (the experiment changed spend); "
f"got {delta_xs[i]!r} for experiment[{i}]"
)
post_spend = x_i + dx_i
if not math.isfinite(post_spend) or post_spend < 0:
raise ValueError(
f"x + delta_x must be finite and >= 0 (post-test spend cannot be "
f"negative or overflow; the saturation curve is evaluated at "
f"x + delta_x); got {x_i!r} + {dx_i!r} = {post_spend!r} for "
f"experiment[{i}]"
)
if not math.isfinite(dy_i):
raise ValueError(f"delta_y must be finite; got {delta_ys[i]!r} for experiment[{i}]")
if dy_i == 0:
zero_lift_rows.append(i)
elif (dx_i < 0) != (dy_i < 0):
# Compare signs directly: dx_i * dy_i can underflow to -0.0 for tiny
# magnitudes (e.g. 1e-200 * -1e-200), so a product < 0 check misses
# wrong-signed rows. Both are strictly nonzero here (delta_x validated
# above, delta_y != 0 handled just above), so the comparison is exact.
wrong_sign_rows.append(i)
row: Dict[str, Any] = {"channel": channels[i]}
if dim_cols:
row.update({k: dim_rows[i][k] for k in dim_cols})
row.update({"x": x_i, "delta_x": dx_i, "delta_y": dy_i, "sigma": sig_i})
rows.append(row)
# Two invalid-row classes, one shared disposition. Wrong sign is what
# PyMC-Marketing rejects with NonMonotonicError; zero lift passes its
# monotonicity check but is degenerate for the strictly-positive Gamma lift
# likelihood (an insignificant experiment that cannot calibrate saturation).
if wrong_sign_rows or zero_lift_rows:
parts = []
if wrong_sign_rows:
parts.append(
f"row(s) {wrong_sign_rows} have sign(delta_y) contradicting "
f"sign(delta_x) (PyMC-Marketing rejects these with NonMonotonicError)"
)
if zero_lift_rows:
parts.append(
f"row(s) {zero_lift_rows} have delta_y == 0 (degenerate for "
f"PyMC-Marketing's strictly-positive Gamma lift likelihood, which its "
f"monotonicity check does not catch)"
)
detail = "; ".join(parts)
if on_wrong_sign == "raise":
raise ValueError(
f"{detail}. Pool experiments, re-scope, or exclude the offending "
f"experiment; or pass on_wrong_sign='drop'/'keep' to handle these rows "
f"explicitly."
)
dropped = set(wrong_sign_rows) | set(zero_lift_rows)
if on_wrong_sign == "drop":
if len(dropped) == len(rows):
raise ValueError(f"on_wrong_sign='drop' would remove every row: {detail}.")
warnings.warn(
f"Dropping invalid lift-test row(s): {detail}.", UserWarning, stacklevel=2
)
rows = [row for i, row in enumerate(rows) if i not in dropped]
else: # "keep"
warnings.warn(
f"Keeping invalid lift-test row(s): {detail}. The frame is NOT valid "
f"input for PyMC-Marketing's lift likelihood.",
UserWarning,
stacklevel=2,
)
columns = ["channel", *dim_cols, "x", "delta_x", "delta_y", "sigma"]
return pd.DataFrame(rows, columns=columns)
def _mask_prelude(arr: np.ndarray) -> str:
"""Serialize a boolean mask into snippet statements (exact for any mask).
Ones-based form mirroring Google's configure-model idiom: initialize
all-True, then for each group of columns sharing the same row pattern,
clear the group and set its True rows. Meridian's contract makes the
all-True base the natural one - channels without an experiment use all
periods - and every column has at least one True by the time this runs
(all-False columns are rejected upstream), so groups stay small. Position
lists go through ``.tolist()`` so plain ints are interpolated (numpy 2.x
reprs ``np.int64`` elements otherwise). Long lines for large masks are an
accepted trade-off: this is generated paste-code, not black-formatted
source.
"""
n_rows, n_cols = arr.shape
lines = [f"roi_calibration_period = np.ones(({n_rows}, {n_cols}), dtype=bool)"]
groups: Dict[Tuple[int, ...], List[int]] = {}
for col in range(n_cols):
if arr[:, col].all():
continue
key = tuple(np.flatnonzero(arr[:, col]).tolist())
groups.setdefault(key, []).append(col)
for rows_key, cols in groups.items():
lines.append(f"roi_calibration_period[:, {cols!r}] = False")
lines.append(f"roi_calibration_period[np.ix_({list(rows_key)!r}, {cols!r})] = True")
body = "\n".join(lines)
return f"import numpy as np\n\n{body}\n\n"
@dataclass(frozen=True)
class ExperimentROI:
"""Per-experiment ROI contribution inside a :class:`MeridianROIPrior`."""
roi: float
roi_sd: float
spend: float
weight: float
def to_dict(self) -> Dict[str, float]:
return {
"roi": self.roi,
"roi_sd": self.roi_sd,
"spend": self.spend,
"weight": self.weight,
}
@dataclass(frozen=True)
class MeridianROIPrior:
"""Lognormal prior parameters for Meridian's ``roi_m``/``mroi_m`` calibration.
``mu``/``sigma`` reproduce ``meridian.model.prior_distribution.
lognormal_dist_from_mean_std(roi_mean, roi_sd)`` exactly (closed form; no Meridian
dependency). ``parameter`` records which Meridian prior the caller is informing
(``"roi_m"`` or ``"mroi_m"``). ``per_experiment`` records each pooled experiment's
ROI, widened sd, spend, and spend weight.
"""
roi_mean: float
roi_sd: float
mu: float
sigma: float
parameter: str = "roi_m"
per_experiment: Tuple[ExperimentROI, ...] = field(default=())
def to_dict(self) -> Dict[str, Any]:
return {
"distribution": "LogNormal",
"parameter": self.parameter,
"roi_mean": self.roi_mean,
"roi_sd": self.roi_sd,
"mu": self.mu,
"sigma": self.sigma,
"per_experiment": [e.to_dict() for e in self.per_experiment],
}
def to_code(
self,
*,
channel: Optional[str] = None,
media_channels: Optional[Sequence[str]] = None,
single_channel: bool = False,
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,
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
prior must not silently calibrate search, social, etc. Channel scope is
therefore always explicit:
- ``to_code(channel="tv", media_channels=["search", "tv"])`` emits vector
``mu``/``sigma`` in exactly the ``media_channels`` order (which must match
the Meridian ``InputData`` channel order); the experiment channel carries
this prior and every other channel keeps Meridian's default.
- ``to_code(single_channel=True)`` emits the scalar snippet for
single-channel models, marked as such in the generated code.
- Calling with neither raises ``ValueError``.
The prior's TIME scope is also required: Meridian's default
``roi_calibration_period=None`` applies the prior over all model times, but
the prior was estimated on the experiment window. Three routes:
- ``roi_calibration_period=<boolean numpy array>`` - the
``(n_media_times, n_media_channels)`` mask, typically built by
:func:`meridian_calibration_mask`. The array is serialized into the
snippet as a short ``np.ones`` + per-column-group assignment prelude.
Boolean dtype, or numeric containing only 0/1 (Google's own docs build
the mask with float zeros), is accepted and cast to bool; masked
arrays are rejected (fill or drop the mask explicitly). An all-False
mask is rejected, and so is ANY entirely-False column: Meridian
aggregates each channel's calibration spend through its mask column,
so an all-False column zeroes it - channels without an experiment
must use ALL periods (Google's documented convention; the builder
sets non-experiment channels all-True). Only ``roi_m`` priors can be
time-scoped: Meridian 1.7.0 rejects a non-None
``roi_calibration_period`` unless the media prior type is ``'roi'``,
so ``parameter="mroi_m"`` priors must use ``full_model_window=True``
(applies to the expression-string route too). Two things this method
cannot verify and the caller owns: the ROW count (no time coordinate
is passed here - the builder guarantees consistency when its
``media_times`` matches the model's coordinates; hand-built arrays
are the caller's responsibility) and the column ORDER/identity (the
mask carries no channel labels - the ``media_channels`` given to the
builder and to this method must be the same list in the same order;
only the count is machine-checked).
- ``roi_calibration_period="<expression>"`` - a Python expression string
interpolated verbatim into the snippet (the pre-existing route).
- ``full_model_window=True`` - acknowledges that the experiment window and
the MMM window coincide. Note Meridian's own guidance: the
configure-model guide states the use of ``roi_calibration_period``
"is not generally recommended because calibrating the ROI of a
specific time period does not necessarily improve estimation of the
overall ROI" - prefer ``full_model_window=True`` when the experiment
evidence reasonably transfers to the full window, and reserve the
mask for evidence genuinely specific to a narrower period.
Snippets use the TensorFlow substrate of TensorFlow Probability; JAX-backed
Meridian users should swap the import for
``tensorflow_probability.substrates.jax`` (noted in the generated code).
"""
if roi_calibration_period is None and not full_model_window:
if self.parameter == "roi_m":
remedy = (
"Pass roi_calibration_period=<expression building the boolean "
"(n_media_times, n_media_channels) mask for your experiment "
"window>, or full_model_window=True to acknowledge that the MMM "
"window and the experiment window coincide, or build the array "
"with meridian_calibration_mask(media_times=..., "
"media_channels=..., channel=..., window=...) and pass it here."
)
else:
# Meridian 1.7.0 accepts roi_calibration_period only for 'roi'
# priors, so the mask routes would fail the next validation -
# recommend the one route that works for this parameter.
remedy = (
f"Meridian 1.7.0 accepts roi_calibration_period only when the "
f"media prior type is 'roi', so a {self.parameter!r} prior has "
f"exactly one route: pass full_model_window=True to acknowledge "
f"the full-window interpretation."
)
raise ValueError(
"to_code() needs the prior's time scope: Meridian's default "
"roi_calibration_period=None applies the prior over ALL model times, "
"but this prior was estimated on the EXPERIMENT window, and ROI "
"differs across windows under varying spend and saturation. " + remedy
)
if roi_calibration_period is not None and full_model_window:
raise ValueError(
"pass either roi_calibration_period or full_model_window=True, not both"
)
if roi_calibration_period is not None and self.parameter != "roi_m":
raise ValueError(
f"Meridian 1.7.0 rejects a non-None roi_calibration_period unless "
f"the media prior type is 'roi' "
f"(ModelSpec._validate_roi_calibration_period), so a "
f"{self.parameter!r} prior cannot be time-scoped via this argument "
f"- pass full_model_window=True to acknowledge the full-window "
f"interpretation instead"
)
mask_prelude = ""
mask_arr: Optional[np.ndarray] = None
if roi_calibration_period is None:
calibration_period = "None"
elif isinstance(roi_calibration_period, np.ndarray):
if isinstance(roi_calibration_period, np.ma.MaskedArray):
raise TypeError(
"roi_calibration_period masked arrays are not accepted (np.asarray "
"would silently drop the mask, turning masked cells into "
"calibration values); fill or drop the mask explicitly"
)
arr = np.asarray(roi_calibration_period)
if arr.ndim != 2:
raise ValueError(
f"roi_calibration_period array must be 2-D with shape "
f"(n_media_times, n_media_channels); got shape {arr.shape}"
)
if arr.size == 0:
raise ValueError(
f"roi_calibration_period array must be non-empty; got shape " f"{arr.shape}"
)
if arr.dtype != np.bool_:
if np.issubdtype(arr.dtype, np.number):
if not np.isin(arr, (0, 1)).all():
raise ValueError(
"roi_calibration_period array must be boolean, or numeric "
"containing only 0/1 (Google's example builds it with "
"np.zeros); it contains values other than 0/1"
)
arr = arr.astype(bool)
else:
raise ValueError(
f"roi_calibration_period array must be boolean, or numeric "
f"containing only 0/1 (Google's example builds it with "
f"np.zeros); got dtype {arr.dtype}"
)
if not arr.any():
raise ValueError(
"roi_calibration_period array is all False, which disables ROI "
"calibration at every time and silently discards the experiment "
"prior's time scope; build the mask with "
"meridian_calibration_mask(...) for your experiment window, or "
"pass full_model_window=True"
)
# Meridian aggregates each channel's calibration spend through its
# mask column (input_data._aggregate_spend einsum), so an all-False
# column zeroes that channel's calibration spend - Google's own
# example sets channels without an experiment to ALL periods.
for col in range(arr.shape[1]):
if not arr[:, col].any():
raise ValueError(
f"roi_calibration_period mask column {col} is entirely "
f"False, which zeroes that channel's aggregated calibration "
f"spend in Meridian; channels without an experiment must "
f"use ALL periods (Google's documented convention) - build "
f"the mask with meridian_calibration_mask(...), which sets "
f"non-experiment channels all-True"
)
mask_arr = arr
mask_prelude = _mask_prelude(arr)