Gaussian-Process Interrupted Time Series (GPITS)#

Overview#

Some interventions reach everybody at once. A Supreme Court ruling binds every state on the day it is handed down; a national lockdown applies to the whole country; a platform-wide product change ships to all users. Difference-in-differences and synthetic control both work by comparing the treated to somebody who stayed untreated, and in these settings nobody did. The comparison has to come from the same unit at a different time.

That design is the interrupted time series: learn how the outcome behaved before the intervention, project that forward, and read the effect off the gap between what happened and the projection. GPITS is Cho’s (2026) version of it. Instead of committing to a straight line or a polynomial, it puts a Gaussian process over the possible trends and conditions on the pre-treatment record. What survives is the set of trajectories consistent with what was observed, and the counterfactual is their posterior mean.

The reason to carry the whole set, and not one fitted curve, is uncertainty. Many trends fit a pre-period equally well and disagree about the future, and they disagree more the further out you look. A segmented regression reports the uncertainty of its own fitted line, which does not grow with the horizon, so its interval is the same width one month out and two years out. A Gaussian process reports how much the trends it still admits disagree at the point you are asking about, and that grows with distance from the data.

When to use this estimator#

Reach for GPITS when there is no untreated comparison unit at all and the pre-treatment series has structure a kernel can name:

  • Universal treatments. A nationwide statute, a court decision, a pandemic — the paper’s application is the 2008 ruling in District of Columbia v. Heller, which bound every U.S. jurisdiction but in practice struck down only D.C.’s handgun ban.

  • Seasonal series with a known cycle. Monthly or daily data where the periodic kernel can carry the cycle and the linear component the trend.

  • Short post-treatment windows. The identifying assumption weakens as the horizon lengthens, so the design is most credible over a handful of periods.

If you have credible untreated donors, use them: Cluster Synthetic Controls (CLUSTERSC), Panel Data Approach (PDA) and Synthetic Business Cycle (SBC) will be more precise. If every unit is treated but the series has recurring local structure without a clean period, Synthetic Historical Control (SHC) matches overlapping historical blocks of the same series instead of imposing a kernel, and infers by conformal permutation. The two divide as follows: reach for Synthetic Historical Control (SHC) when you cannot name the cycle, and for GPITS when you can, or when you want an interval that widens with the forecast horizon.

Notation#

Let \(Y_{it}\) be the outcome for unit \(i\) at time \(t\), and \(D_{it} = \mathbf{1}(t \ge t_0)\) the treatment indicator, which turns on at the same \(t_0\) for every unit. Write \(Y_{it}(1)\) and \(Y_{it}(0)\) for the potential outcomes. Let \(X_{it}\) collect the observed inputs — the calendar time itself, plus any covariates — and \(U_{it}\) the unobservables.

The untreated outcome is governed by

\[Y_{it}(0) = h_i(X_{it}, U_{it}) + \varepsilon_{it}, \qquad \mathbb{E}[\varepsilon_{it} \mid X_{it}, U_{it}] = 0 ,\]

and the object the design targets is the best prediction available from observables,

\[g_i(X_{it}) := \mathbb{E}[Y_{it}(0) \mid X_{it}] .\]

The estimand is the unit-period effect against that conditional mean,

\[\tau_{it^*} = Y_{it^*}(1) - g_i(X_{it^*}), \qquad t^* \ge t_0 ,\]

with the average and cumulative versions over the post-treatment window following from it.

Assumptions#

  1. Consistency and no anticipation. \(Y_{it} = D_{it} Y_{it}(1) + (1 - D_{it}) Y_{it}(0)\), and before \(t_0\) the treated and untreated potential outcomes coincide.

    Remark. The second half is what makes the pre-period usable as training data. It rules out anticipatory drift, not stable expectations: a belief held at a constant level throughout the pre-period is absorbed into the learned relationship. Behaviour that shifts in advance of the intervention is the failure mode, and it is what the placebo checks below are built to detect.

  2. Mean sufficiency. For every unit and period, \(\mathbb{E}[Y_{it}(0) \mid X_{it}, U_{it}] = \mathbb{E}[Y_{it}(0) \mid X_{it}]\).

    Remark. This is the assumption that replaces the donor pool, and it is a restriction on first moments only, weaker than conditional independence. It says that once you know the calendar date and the covariates, knowing the unobservables would not change your expectation of the untreated outcome. It fails when something other than the intervention starts moving the series after \(t_0\) — a concurrent policy, a structural break, a shock that took values in the post-period it never took before. Nothing in the pre-period can confirm it, which is why shorter post-treatment windows are more credible and why domain knowledge about concurrent events is doing real work.

  3. A kernel that can represent the counterfactual. The trend must lie in, or close to, the space of functions the kernel treats as plausible, at a complexity the pre-period can pin down.

    Remark. This is where the choice between kernel="gaussian" and kernel="gaussian_periodic_linear" bites, and it is not a cosmetic setting. The Gaussian kernel is stationary: far from the training data its posterior reverts to the prior mean and its band flattens at a ceiling, so it cannot carry a trend forward. The combined kernel adds a periodic component for the cycle and a linear component for the trend, and it is the working form for any series that is going somewhere. Cho reports the practical size of this: on the Heller series a unit-variance trend costs a complexity budget of 0.96 under the combined kernel against 3.52 under the Gaussian alone, and the reported interval is short by that factor when the budget exceeds one.

Inference and diagnostics#

The counterfactual band is the Gaussian-process posterior variance. Cho’s Proposition 2 gives it a frequentist reading: it bounds the worst-case disagreement, at the period you are asking about, among the functions in the kernel’s unit ball that agree with everything observed before the intervention. The interval is therefore calibrated to the worst case the model class admits, not to the typical case.

Two consequences follow, and both matter in practice. The intervals cover at or above their nominal rate, which is what you want from an extrapolation. And they are wide: in replication of the paper’s own simulation the coverage sits at 1.000 in most cells, with intervals a median of 2.3 times wider than a segmented regression’s. Coverage is bought with power, so an effect that is small relative to the pre-period noise will sit inside the band. The Heller effect survives because it is roughly 20 times D.C.’s pre-period standard deviation.

Set placebo_periods to run the temporal placebo check. It refits on everything before each of the last few pre-treatment periods and predicts one step ahead, where the true effect is zero. The check is one-sided by construction: a confounder inside the training window is absorbed into the fit and leaves the placebo clean even when the same disturbance breaks mean sufficiency after \(t_0\). A clean placebo reports the absence of detected instability, not the presence of identification.

result.inference.ci_lower / ci_upper give the ATT interval, and result.cumulative_ci the interval on each running cumulative total. Both use the full post-period posterior covariance, not its diagonal, because successive counterfactual errors covary and summing variances alone would understate a running total.

Example#

import pandas as pd
from mlsynth import GPITS

# Monthly handgun background checks in D.C., 2002-07 to 2008-10;
# the Heller decision lands 2008-07.
df = pd.read_csv("basedata/dc_handgun_heller.csv", parse_dates=["date"])

res = GPITS({
    "df": df,
    "outcome": "handgun_rate",
    "treat": "treated",
    "unitid": "unit",
    "time": "date",
    "covariates": ["month"],
    "categorical_covariates": ["month"],
    "kernel": "gaussian_periodic_linear",
    "period": 12,
    "placebo_periods": 4,
}).fit()

print(res.effects.att)                  # mean monthly effect
print(res.cumulative_effect[-1])         # 15.13 checks per 100k over 4 months
print(res.cumulative_ci[-1])             # (12.97, 17.30)
print(res.placebo.all_cover)             # True

The plot helpers return figures and leave showing and saving to you. plot_gpits draws the fit panel, which is what display_graphs shows:

from mlsynth.utils.gpits_helpers import plot_gpits, plot_gpits_panels

fig = plot_gpits(res)
fig.savefig("heller_fit.png", dpi=150)

The fit panel draws the pre-period band in grey and the post-period band in the counterfactual colour, because the two are different quantities: before the intervention the band is a fit’s uncertainty, after it the band is an extrapolation’s.

plot_gpits_panels returns the four panels of the paper’s own plotting code (plot.gp_its in the replication repository), keyed by its names:

figs = plot_gpits_panels(res)
figs["fit"]          # observed points, fitted trend, counterfactual
figs["pointwise"]    # per-period effects, placebo window shaded separately
figs["cumulative"]   # running total with its interval
figs["average"]      # running average, ending at the ATT

The pointwise panel is the one that carries the paper’s argument: the placebo periods sit flat near zero immediately before an intervention the observed series jumps away from. Time on the effect panels is measured from the intervention, so period 0 is the first treated period.

Verification#

GPITS reproduces the paper’s empirical result exactly and is cross-validated against the author’s own R implementation. See GPITS — Cho (2026), the Heller decision for the numbers, and benchmarks/cases/gpits.py for the durable case.

Core API#

class mlsynth.GPITS(config: GPITSConfig | dict)#

Bases: object

Gaussian-process interrupted time series (Cho 2026).

Estimates a treated unit’s untreated counterfactual from its own pre-treatment history, with no donor pool, and reports a band that widens as the forecast extends past the data. Use it when the treatment reached everyone at once so no comparison unit survives.

Parameters:

config (GPITSConfig or dict) – Configuration object. See mlsynth.utils.gpits_helpers.config.GPITSConfig. Key fields: kernel and period (the covariance structure), covariates / categorical_covariates (extra design columns, typically calendar indicators), and placebo_periods (the Section 3.3 diagnostic).

Returns:

GPITSResults – Counterfactual and pointwise band over every period, the ATT with its interval, the running cumulative effect with intervals from the full post-period covariance, fit diagnostics, and the placebo checks.

Notes

Any untreated units in the frame are ignored: the counterfactual is built from the treated unit’s own series alone. The hyperparameters are chosen in two stages, as the paper specifies – the length-scale by a rule that reads only the covariates, then the noise variance by marginal likelihood with the length-scale fixed. Pin either through length_scale or noise_variance and the result records that you did.

Examples

>>> from mlsynth import GPITS
>>> res = GPITS({"df": panel, "outcome": "handgun_rate", "treat": "treated",
...              "unitid": "unit", "time": "date",
...              "covariates": ["month"],
...              "categorical_covariates": ["month"],
...              "period": 12}).fit()
>>> res.effects.att
fit() GPITSResults#

Run the GPITS pipeline and return structured results.

class mlsynth.utils.gpits_helpers.config.GPITSConfig(*, df: ~pandas.DataFrame, outcome: str, treat: str, unitid: str, time: str, display_graphs: bool = True, save: bool | str = False, counterfactual_color: ~typing.List[str] = <factory>, treated_color: str = 'black', plot: ~mlsynth.config_models.PlotConfig = <factory>, kernel: str = 'gaussian', period: float | None = None, covariates: ~typing.List[str] | None = None, categorical_covariates: ~typing.List[str] | None = None, length_scale: float | None = None, noise_variance: float | None = None, interval_type: str = 'prediction', alpha: float = 0.05, placebo_periods: int | None = None)#

Settings for Gaussian-process interrupted time series (Cho 2026).

alpha: float#
categorical_covariates: List[str] | None#
check_gpits_params() GPITSConfig#
covariates: List[str] | None#
interval_type: str#
kernel: str#
length_scale: float | None#
model_config: ClassVar[ConfigDict] = {'arbitrary_types_allowed': True, 'extra': 'forbid'}#

Configuration for the model, should be a dictionary conforming to [ConfigDict][pydantic.config.ConfigDict].

noise_variance: float | None#
period: float | None#
placebo_periods: int | None#
class mlsynth.utils.gpits_helpers.structures.GPITSResults(*, effects: ~mlsynth.config_models.EffectsResults | None = None, fit_diagnostics: ~mlsynth.config_models.FitDiagnosticsResults | None = None, time_series: ~mlsynth.config_models.TimeSeriesResults | None = None, weights: ~mlsynth.config_models.WeightsResults | None = None, inference: ~mlsynth.config_models.InferenceResults | None = None, method_details: ~mlsynth.config_models.MethodDetailsResults | None = None, sub_method_results: ~typing.Dict[str, ~typing.Any] | None = None, additional_outputs: ~typing.Dict[str, ~typing.Any] | None = None, raw_results: ~typing.Dict[str, ~typing.Any] | None = None, execution_summary: ~typing.Dict[str, ~typing.Any] | None = None, plot_config: ~mlsynth.config_models.PlotConfig | None = None, inputs: ~mlsynth.utils.gpits_helpers.structures.GPITSInputs, design: ~mlsynth.utils.gpits_helpers.structures.GPITSDesign, att_value: float, att_interval: ~typing.Tuple[float, float], observed: ~numpy.ndarray, cf_series: ~numpy.ndarray, gap_series: ~numpy.ndarray, counterfactual_lower: ~numpy.ndarray, counterfactual_upper: ~numpy.ndarray, cumulative_effect: ~numpy.ndarray, cumulative_ci: ~typing.List[~typing.Tuple[float, float]], time_labels: ~numpy.ndarray, fit_diagnostics_detail: ~typing.Dict[str, ~typing.Any], placebo: ~mlsynth.utils.gpits_helpers.structures.GPITSPlacebo | None = None, metadata: ~typing.Dict[str, ~typing.Any] = <factory>)#

Public container returned by mlsynth.GPITS.fit().

Parameters:
  • inputs (GPITSInputs) – Preprocessed series and design.

  • design (GPITSDesign) – Fitted hyperparameters, posterior mean, and post-period covariance.

  • att_value (float) – Mean post-treatment gap.

  • att_interval (tuple of float) – Interval for the ATT, from the post-period covariance block. The base class exposes it as att_ci.

  • observed, cf_series, gap_series (np.ndarray) – Length-T series. The base class exposes the latter two as counterfactual and gap.

  • counterfactual_lower, counterfactual_upper (np.ndarray) – Pointwise band on the counterfactual.

  • cumulative_effect (np.ndarray) – Running sum of the post-period gaps.

  • cumulative_ci (list of (float, float)) – Interval for each cumulative total, using the full covariance block up to that period.

  • time_labels (np.ndarray) – Period labels, length T.

  • fit_diagnostics_detail (dict) – Pre/post RMSE and pre-period R-squared.

  • placebo (GPITSPlacebo or None) – Placebo checks, when requested.

  • metadata (dict) – Free-form diagnostics.

att_interval: Tuple[float, float]#
att_value: float#
cf_series: np.ndarray#
counterfactual_lower: np.ndarray#
counterfactual_upper: np.ndarray#
cumulative_ci: List[Tuple[float, float]]#
cumulative_effect: np.ndarray#
design: GPITSDesign#
fit_diagnostics_detail: Dict[str, Any]#
gap_series: np.ndarray#
inputs: GPITSInputs#
metadata: Dict[str, Any]#
model_config: ClassVar[ConfigDict] = {'arbitrary_types_allowed': True, 'extra': 'forbid', 'frozen': True, 'json_encoders': {<class 'numpy.ndarray'>: <function BaseEstimatorResults.Config.<lambda>>}}#

Configuration for the model, should be a dictionary conforming to [ConfigDict][pydantic.config.ConfigDict].

observed: np.ndarray#
placebo: GPITSPlacebo | None#
time_labels: np.ndarray#

References#

Cho, S. (2026). “Let Time Tell: Identification and Gaussian Process Estimation for Interrupted Time Series.” arXiv:2608.20610.

Rasmussen, C. E., & Williams, C. K. I. (2006). Gaussian Processes for Machine Learning. MIT Press.