CMBSTS — Multivariate Bayesian Structural Time Series (Menchetti & Bojinov 2022)

Estimator:

CMBSTS — Causal Multivariate Bayesian Structural Time Series — mlsynth.CMBSTS

Source:

Menchetti, Fiammetta, and Iavor Bojinov (2022), “Estimating the Effectiveness of Permanent Price Reductions for Competing Products Using Multivariate Bayesian Structural Time Series Models,” Annals of Applied Statistics 16(1): 414–435 [MenchettiBojinov2022].

Reference implementation:

the authors’ CausalMBSTS R package (Bojinov & Menchetti 2020).

Replication type:

Cross-validation — the mlsynth port is checked cell-by-cell against the R package.

Status:

Verified — the port reproduces CausalMBSTS within Monte-Carlo error on the package’s vignette and on the Florence supermarket study.

Validation strategy

CMBSTS is a faithful NumPy/SciPy port of the CausalMBSTS R package: a multivariate structural state space (trend, optional seasonal and cycle, and a spike-and-slab regression), a Gibbs sampler with a Durbin–Koopman simulation smoother for the latent states, and a posterior-predictive counterfactual. The port is validated by running the genuine R package and comparing its output to mlsynth’s on identical inputs.

Seed-faithful agreement is not achievable: the R and Python samplers draw from independent random-number streams, so individual draws differ. The honest target is therefore agreement of the posterior summaries — the per-series temporal- average effect and its credible bounds — within Monte-Carlo error.

The vignette (durable benchmark)

The package vignette generates a bivariate weekly series, \(y_{1t} = 3\sin(2t) + e\), \(y_{2t} = 2\cos(2t) + e\), with a fictional \(+2\) intervention, and fits a trend-plus-cycle model (cycle period 75). The data ship as basedata/cmbsts_vignette.csv (generated by the R DGP with set.seed(1), so the noise is identical to the reference run). With the estimator’s default Inverse-Wishart prior (\(0.01 \cdot \mathrm{diag}(\widehat{\mathrm{var}}(\mathbf{y}_{\text{pre}}))\)) and 2000 Gibbs iterations:

Quantity	R CausalMBSTS	mlsynth CMBSTS
series 1 (treated) ATT	2.4469	2.457
series 2 (group) ATT	2.2756	2.287
series 1 95% CI lower	0.352	0.323
series 1 95% CI upper	4.567	4.642

The posterior-mean effects agree to about \(0.01\) and the credible bounds to a few hundredths, and the injected \(+2\) sits inside both bands. This is captured durably in benchmarks/cases/cmbsts_vignette.py.

The supermarket study

The paper’s empirical application is the Florence supermarket chain’s permanent discount on store-brand cookies, modelled pair-by-pair (a store brand and its direct competitor) with a trend-plus-weekly-seasonal model and a regression block of calendar dummies, a frozen store price, the competitor price, and ten wine-control series. Table 3 reports the temporal-average effect per pair at the one-month horizon, finding significant positive store-brand effects on pairs 4, 7 and 10 and no significant competitor effects. The authors’ Supplement B data ship in basedata/cmbsts_supermarket/, and benchmarks/cases/cmbsts_supermarket.py reproduces this at the one-month horizon.

Two things hold. The mlsynth store-brand effects cross-validate against the R CausalMBSTS package on identical controls and prior — pair 4 48.9 vs 47.4, pair 7 80.0 vs 78.1, pair 10 12.4 vs 12.3 — within Monte-Carlo error. And the substantive Table 3 result reproduces: large positive store effects with no significant competitor effect. The controls here are screened by dynamic time warping (the package’s method) through the optional fastdtw package rather than the authors’ MarketMatching, so the control set differs; the pair-10 effect is strictly significant in both implementations, while pairs 4 and 7 sit on the zero boundary where Monte-Carlo noise alone moves the lower credible bound across (the package and the port land on opposite sides). The published strict significance on all three pairs uses the MarketMatching controls at niter = 2200.

A reference-implementation difference worth noting: CausalMBSTS adds a single observation-noise draw per posterior draw and recycles it across the forecast horizon, so that noise does not average out of the temporal-average effect and widens its credible interval. The mlsynth port reproduces this behaviour to match the reference; a per-period-independent draw would give a narrower, arguably more defensible interval.