MCMC Inference for CausalImpact

Summary

Posterior inference in the BSTS model uses a Gibbs sampler that alternates between a data-augmentation step (sampling states $\alpha$ given parameters $\theta$, using the Kalman filter and fast mean smoother) and a parameter-simulation step (sampling $\theta$ given states, with the spike-and-slab Gibbs draw for variable selection). The algorithm is linear in the number of time points and runs in < 30 seconds for typical datasets.

Overview

The posterior $p(\theta ,\alpha \mid \mathbf{y})$ is not available in closed form due to the spike-and-slab prior. The Gibbs sampler alternates two steps:

Gibbs Sampler Steps

Step 1 — Data Augmentation (State Simulation)

Sample the full state sequence $\alpha =({\alpha }_1,\dots ,{\alpha }_m)$ given parameters $\theta$ and data ${\mathbf{y}}_{1:n}$:

$$(\alpha \mid {\mathbf{y}}_{1:n},\theta )$$

Algorithm: Uses the simulation smoother of Durbin & Koopman (2002), which improves on the earlier forward-filtering, backward-sampling algorithms (Carter & Kohn 1994, Frühwirth-Schnatter 1994).

Key property: Because $p({\mathbf{y}}_{1:n}\mid \alpha ,\theta )$ is jointly multivariate Gaussian, the variance of $p(\alpha \mid {\mathbf{y}}_{1:n},\theta )$ does not depend on ${\mathbf{y}}_{1:n}$. The sampler:

1. Generates $({\tilde{y}}_{1:n}^{\ast },{\alpha }^{\ast })\sim p({\mathbf{y}}_{1:n},\alpha \mid \theta )$
2. Subtracts $E({\alpha }^{\ast }\mid {\tilde{y}}_{1:n}^{\ast },\theta )$ to get zero-mean noise
3. Adds $E(\alpha \mid {\mathbf{y}}_{1:n},\theta )$ (via Kalman filter) to restore the correct mean

Computational complexity: Linear in $m$ (total time points, pre + post), quadratic in $d$ (state dimension). For $m=500$, $J=10$ covariates, 10,000 iterations: < 30 seconds.

Step 2 — Parameter Simulation

Sample $\theta =({\sigma }_{\mu }^2,{\sigma }_{\delta }^2,\dots ,\varrho ,\beta ,{\sigma }_{\epsilon }^2)$ given states $\alpha$ and data.

For variance parameters (${\sigma }_{\mu }^2$, ${\sigma }_{\delta }^2$, etc.): Because error terms ${\eta }_t={\alpha }_{t+1}-T_t{\alpha }_t$ are available given $\alpha$, the posterior is Gamma by conjugacy (from the inverse-Gamma prior in Eq. 2.7).

For static regression coefficients ($\varrho$, ${\beta }_{\varrho }$, ${\sigma }_{\epsilon }^2$): Gibbs sampling from the spike-and-slab posterior (see Spike-and-Slab Prior for Covariate Selection). Each ${\varrho }_j$ is drawn independently given ${\varrho }_{-j}$, then ${\beta }_{\varrho }\mid \varrho ,{\sigma }_{\epsilon }^2$ and $1/{\sigma }_{\epsilon }^2\mid \varrho ,{\beta }_{\varrho }$ are drawn using conjugate formulae.

Posterior Predictive Simulation

After fitting the model on pre-intervention data ${\mathbf{y}}_{1:n}$, the key quantity is the posterior predictive distribution over counterfactuals:

$$p({\tilde{\mathbf{y}}}_{n+1:m}\mid {\mathbf{y}}_{1:n},{\mathbf{x}}_{1:m})\text{\hspace{1em}(2.14)}$$

This is the distribution of what would have happened had no intervention occurred. It:

• Is conditioned only on pre-intervention outcomes and all control series (not on parameter estimates)
• Integrates out all $\beta$ and ${\sigma }_{\epsilon }^2$ — no commitment to any particular set of covariates
• Is a joint distribution over all post-intervention time points (not a collection of marginals) — preserves serial correlation

Sampling: Each Gibbs iteration draws a complete counterfactual trajectory ${\tilde{\mathbf{y}}}_{n+1:m}^{(\tau )}$ using the Kalman filter run forward through the post-intervention period.

Why Integrating Out Parameters Matters

• Integrating out $\beta$ and ${\sigma }_{\epsilon }^2$ means no arbitrary covariate selection — the posterior predictive averages over all candidate subsets weighted by posterior probability
• Integrating out ${\sigma }_{\epsilon }^2$ means no commitment to point estimates of noise — full propagation of uncertainty
• The result is wider but properly calibrated uncertainty intervals

Connections

• Uses Local Linear Trend and Seasonality Kalman filter/smoother
• Uses Spike-and-Slab Prior for Covariate Selection Gibbs update for $\varrho$
• Output feeds Counterfactual Impact Estimation
• Same MCMC paradigm as MCMC Basics (Gibbs sampling)

See Also

• Bayesian Structural Time-Series Model — model being inferred
• Counterfactual Impact Estimation — how posterior predictive draws are used
• The Kalman Filter — the simulation smoother builds on the Kalman recursion