Confidence Intervals under Staggered Adoption (CAST)

Confidence Intervals under Staggered Adoption (CAST)#

Overview#

Policies rarely arrive everywhere at once. States expand Medicaid in different years, countries join a trade bloc at different times, stores roll out a promotion region by region. This is staggered adoption: every unit is eventually treated (or never), but each on its own schedule.

Most panel estimators answer the pooled question – what was the average effect across the treated? CAST (Xia, Yan and Wainwright 2025) answers a finer one. It returns a confidence interval for every treated unit in every period: not just “expansion reduced uninsurance by 2.6 percent on average in 2019”, but “Louisiana’s reduction in 2019 was 3.1 points, with an interval that excludes zero”. Aggregates are still available, and because they are built from the same entrywise variances they come with standard errors of their own.

When to use this estimator#

  • Your treatment is staggered and you care about which units responded, not only the average. Reporting a per-unit effect with an interval is the natural output when the units are themselves of interest (states, countries, large accounts).

  • You want a population-weighted total rather than an unweighted mean – “how many people gained coverage”, not “the average percentage-point change” – with a standard error attached.

  • Your panel is reasonably large. The intervals are asymptotic; see the coverage remark under Assumption 3.

Reach for something else when you only need a pooled ATT (Synthetic Difference-in-Differences (SDID), Sequential Synthetic Difference-in-Differences (Sequential SDiD) and Partially Pooled SCM (PPSCM) are purpose-built for that and report event-study aggregates), or when there is no never-treated unit to anchor the comparison.

Notation#

There are \(N\) units indexed \(i\) observed over \(T\) periods indexed \(t\). Unit \(i\) adopts treatment at time \(t_i\), and treatment is absorbing: once treated, always treated. A unit with \(t_i > T\) is never treated.

Write \(Y_{i,t}\) for the observed outcome and \(M^\star_{i,t}\) for the mean outcome that would have been observed without treatment. The individual treatment effect is

\[\tau_{i,t} = Y_{i,t} - M^\star_{i,t}, \qquad t \ge t_i .\]

Because \(Y_{i,t}\) is observed, estimating \(\tau_{i,t}\) is the same problem as estimating the counterfactual mean \(M^\star_{i,t}\), and a confidence interval for one is an interval for the other. More general targets are the weighted effects

\[\tau_{w,t} = \sum_{i=1}^{N} w_i \, \tau_{i,t} \, \mathbf{1}(t_i \le t),\]

which recover the average effect on the treated when \(w\) is uniform over the treated units, and a population total when \(w\) counts people.

Mathematical formulation#

The staircase#

Sort the units by how long they are treated. The observed cells then form a descending staircase: the never-treated units contribute a full row, the earliest adopters the shortest. Everything above the staircase is counterfactual and must be imputed.

The paper’s device is to notice that such a staircase decomposes into overlapping four-block problems. Each has the shape

\[\begin{split}Y = \begin{pmatrix} Y_a & Y_b \\ Y_c & ? \end{pmatrix},\end{split}\]

where \(Y_a\) (\(N_1 \times T_1\)) is fully observed, \(Y_b\) and \(Y_c\) extend it right and down, and the bottom-right block is the unknown. Every unobserved cell in the original panel falls inside one such sub-problem, so solving the four-block case solves the general one.

Estimation#

Assume the counterfactual means have low rank \(r\), that is \(M^\star = U^\star \Sigma^\star (V^\star)^\top\) with \(r\) factors. The four-block routine estimates the two subspaces from the observed blocks and transports one onto the other:

  1. the column subspace \(\hat U\) from a rank-\(r\) SVD of the left blocks \([Y_a; Y_c]\);

  2. the row subspace from a rank-\(r\) SVD of the upper blocks \([Y_a \; Y_b]\), which also gives a denoised \(\hat M_b\);

  3. the imputation, by regressing the denoised block onto the shared subspace,

\[\hat M_d = \hat U_2 (\hat U_1^\top \hat U_1)^{-1} \hat U_1^\top \hat M_b ,\]

where \(\hat U_1\) and \(\hat U_2\) are the top and bottom row blocks of \(\hat U\). If the mean matrix is exactly rank \(r\) and the observations are noiseless, this recovers the missing block exactly.

Inference#

The same subspaces give the uncertainty. Residuals are formed on the observed entries, \(\hat E = (Y - \hat M) \odot A\) with \(A\) the observation indicator, and the entrywise variance is

\[\hat\gamma_{i,t} = \sum_{k \le N_1} \hat E_{k,t}^2 \Bigl[\hat U_{i,\cdot} (\hat U_1^\top \hat U_1)^{-1} \hat U_{k,\cdot}^\top\Bigr]^2 + \sum_{s \le T_1} \hat E_{i,s}^2 \Bigl[\hat V_{t,\cdot} (\hat V_1^\top \hat V_1)^{-1} \hat V_{s,\cdot}^\top\Bigr]^2 ,\]

giving the interval \(\hat M_{i,t} \pm \Phi^{-1}(1 - \alpha/2)\,\hat\gamma_{i,t}^{1/2}\). The two terms are the noise propagated through the column and row subspace projections respectively.

Assumptions#

Assumption 1 (Absorbing staggered treatment). Each unit adopts at some \(t_i\) and stays treated; at least one unit is never treated.

Remark. The never-treated units anchor the staircase – they are the only rows observed for the full horizon, so they identify the row subspace over the post-treatment periods. mlsynth raises MlsynthDataError when none is present.

Assumption 2 (Low-rank counterfactual means). \(M^\star\) has rank \(r\), small relative to \(N\) and \(T\).

Remark. The rank is the method’s one substantive tuning choice. Too small under-fits the factor structure and leaves signal in the residual; too large lets the factors absorb the treatment effect itself. Set it with rank or choose rank_method: "broken_stick" (the default, matching the authors’ reference implementation), "svht" (Gavish-Donoho optimal hard threshold) or "ratio" (largest spectral gap). Inspecting the singular values yourself is worthwhile – the paper selects by eye from a scree plot.

Assumption 3 (Independent sub-Gaussian noise). The errors \(\varepsilon_{i,t} = Y_{i,t} - M^\star_{i,t}\) are independent, mean-zero and sub-Gaussian across both units and time.

Remark. Independence across time is the demanding half. Macro panels such as state-year outcomes are typically serially correlated, which this rules out; where that is a concern the intervals should be read as optimistic. Coverage is also asymptotic: in our simulations on a rank-2 design, coverage of \(M^\star\) was about 0.82 at \(N=30, T=24\), rising to 0.91 at \(N=60, T=40\) and 0.94 at \(N=240, T=160\) against a nominal 0.95. Small panels undercover.

Assumption 4 (Incoherence). The singular vectors of \(M^\star\) are delocalised, so no single unit or period dominates the factor structure.

Remark. Standard in matrix completion. It fails when one unit is effectively its own factor – an outlier so idiosyncratic that no combination of the others describes it.

Inference and diagnostics#

The entrywise panel. effects_panel, std_errors, ci_lower_panel and ci_upper_panel are \(N \times T\) arrays, finite exactly on the treated cells and nan elsewhere. CASTResults.effect_table flattens them into tidy (unit, period, effect, se, ci_lower, ci_upper) rows.

Per-period aggregates. period_effects maps each post-treatment period to (effect, standard error) for the weighted average of the units treated by then. Supply weights for a population-weighted version, and set renormalize_weights=False to report a total rather than a mean.

The significance map. significance counts, per period, how many treated units have a significantly positive, significantly negative, or null effect at alpha; CASTResults.significant_units() returns the units themselves.

An important caveat on what is being tested. The interval covers the counterfactual mean \(M^\star_{i,t}\), so the effect \(\tau_{i,t} = Y_{i,t} - M^\star_{i,t}\) still contains that cell’s own idiosyncratic noise \(\varepsilon_{i,t}\), which the variance deliberately excludes. A single unit-period can therefore be flagged as significant on the strength of its own noise draw. Per-unit significance is best read as suggestive, with the per-period aggregates – where the noise averages out – carrying the evidential weight.

Example#

import numpy as np, pandas as pd
from mlsynth import CAST

# staggered panel: two cohorts adopt at t=16 and t=20, ten units never treated
rng = np.random.default_rng(0)
N, T, r = 30, 24, 2
U, V = rng.normal(size=(N, r)), rng.normal(size=(T, r))
Y = U @ V.T + 0.2 * rng.normal(size=(N, T))
adopt = np.full(N, T); adopt[10:20] = 16; adopt[20:] = 20
for i, a in enumerate(adopt):
    if a < T:
        Y[i, a:] += 3.0                      # planted effect

df = pd.DataFrame(
    [dict(u=f"u{i:02d}", t=t, y=Y[i, t], d=int(t >= adopt[i]))
     for i in range(N) for t in range(T)]
)

res = CAST({"df": df, "outcome": "y", "treat": "d", "unitid": "u",
            "time": "t", "rank": 2}).fit()

print(res.att, res.att_ci)                   # treated average, ~3.0
print(res.period_effects[16])                # (effect, se) in the first post period
print(res.significance[16])                  # positive / negative / null counts
print(res.effect_table[:3])                  # entrywise rows

For a population-weighted total instead of a mean, pass per-unit weights and turn off renormalisation:

res = CAST({..., "weights": population.tolist(),
            "renormalize_weights": False}).fit()

Verification#

CAST is cross-validated against the authors’ own CAST-panel package and reproduces their Affordable Care Act analysis. The durable case lives in benchmarks/cases/cast_aca.py and the full replication story, including a standard-error discrepancy we found in the reference package, is on CAST — Inference under Staggered Adoption (Xia, Yan & Wainwright 2025).

Core API#

class mlsynth.CAST(config: CASTConfig | dict)#

Confidence intervals under staggered adoption (Xia et al. 2025).

Parameters:

config (CASTConfig or dict) – Configuration. See mlsynth.config_models.CASTConfig.

Returns:

CASTResults – An EffectResult whose flat accessors summarize the treated average, plus the entrywise effect / standard-error / interval panels.

fit() → CASTResults#

Run the CAST pipeline end to end.

Configuration#

class mlsynth.config_models.CASTConfig(*, 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>, rank: ~typing.Annotated[int | None, ~annotated_types.Ge(ge=1)] = None, rank_method: ~typing.Literal['broken_stick', 'svht', 'ratio'] = 'broken_stick', max_rank: ~typing.Annotated[int, ~annotated_types.Ge(ge=1)] = 10, alpha: ~typing.Annotated[float, ~annotated_types.Gt(gt=0.0), ~annotated_types.Lt(lt=1.0)] = 0.05, weights: list | None = None, renormalize_weights: bool = True)#

Configuration for CAST – confidence intervals for treatment effects under staggered adoption (Xia, Yan & Wainwright 2025, Inference under Staggered Adoption: Case Study of the Affordable Care Act).

CAST estimates the counterfactual outcome matrix by a spectral four-block routine and returns a confidence interval for every treated unit-period, not only for the pooled effect. Aggregates (the average effect on the treated, or any weighted version of it) are reported with standard errors derived from the same entrywise variances.

Parameters:

  • rank (int or None) – Rank r of the underlying mean matrix. None (default) selects it from the spectrum of the observed control block via rank_method. The rank is the method’s one substantive tuning choice: too small under-fits the factor structure, too large absorbs the treatment effect.

  • rank_method ({“broken_stick”, “svht”, “ratio”}) – Selector used when rank is None. "broken_stick" matches the authors’ reference package; "svht" is the Gavish-Donoho optimal hard threshold; "ratio" takes the largest singular-value gap.

  • max_rank (int) – Upper bound for the automatic search.

  • alpha (float) – Two-sided level for the entrywise and aggregate confidence intervals.

  • weights (list or None) – Optional weights for the aggregate effect tau_w,t = sum_i w_i tau_i,t, in the order the units appear in the panel (sorted unit labels). Either a length-N vector (one weight per unit, constant over time) or an N x T matrix for a weight that varies by period – the paper’s general bilinear functional, which is what a population-weighted total needs when the population itself moves over time. None gives the equally weighted average over treated units.

  • renormalize_weights (bool) – Divide the supplied weights by their total over the units treated at each period, so the aggregate is a weighted mean. Set False to report a population total (the paper’s “population effect”).

Results#

class mlsynth.utils.cast_helpers.structures.CASTResults(*, 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.cast_helpers.structures.CASTInputs, rank: int, effects_panel: ~numpy.ndarray, std_errors: ~numpy.ndarray, ci_lower_panel: ~numpy.ndarray, ci_upper_panel: ~numpy.ndarray, counterfactual_panel: ~numpy.ndarray, period_effects: ~typing.Dict[~typing.Any, ~typing.Any], significance: ~typing.Dict[~typing.Any, ~typing.Any], metadata: ~typing.Dict[str, ~typing.Any] = <factory>)#

Top-level container returned by mlsynth.CAST.fit().

An EffectResult: the standardized sub-models describe the average treated unit (so att / att_ci / counterfactual / gap resolve through the base contract), while the entrywise panel below carries the per-unit-per-period inference that is the method’s reason for existing.

Parameters:

  • inputs (CASTInputs) – Preprocessed panel.

  • rank (int) – Rank used for the factor structure.

  • effects_panel (np.ndarray) – (N, T) individual treatment effects; nan at untreated cells.

  • std_errors (np.ndarray) – (N, T) entrywise standard errors; nan at untreated cells.

  • ci_lower_panel, ci_upper_panel (np.ndarray) – (N, T) entrywise confidence bounds; nan at untreated cells.

  • counterfactual_panel (np.ndarray) – (N, T) estimated counterfactual outcomes (observed values retained at untreated cells).

  • period_effects (dict) – {period label: (effect, standard error)} – the aggregate treated effect per period, with the weighting requested in the config.

  • significance (dict) – {period label: {"positive": int, "negative": int, "null": int, "n_treated": int}} – the per-period sign breakdown at alpha.

  • metadata (dict) – Free-form pipeline diagnostics (rank selection, weighting, …).

property effect_table: List[Dict[str, Any]]#

Tidy (unit, period, effect, se, ci_lower, ci_upper) rows for the treated cells – the entrywise interval panel, long-format.

significant_units(period: Any, level: float = 0.05) → Dict[Any, float]#

{unit: effect} for units whose effect at period excludes zero.

Uses the stored entrywise intervals when level matches the fitted alpha; otherwise rescales the interval to level.

Helper Modules#

Numerical kernels for CAST (Xia, Yan & Wainwright 2025, arXiv:2412.09482).

A staggered-adoption panel, after sorting units by how long they are treated, takes a staircase shape: the observed (control) cells form a descending staircase and the counterfactual cells sit above it. The paper’s insight is that such a staircase decomposes into a collection of four-block problems

\[\begin{split}Y = \begin{pmatrix} Y_a & Y_b \\ Y_c & ? \end{pmatrix},\end{split}\]

each of which is solved by a spectral routine (four_block_est(), Algorithm 1) that estimates the column subspace from the left blocks, denoises the upper blocks, and regresses one onto the other to impute ?. The same subspaces yield a closed-form entrywise variance (entrywise_variance(), equations 6a-6c), which is what makes per-cell confidence intervals possible.

mlsynth.utils.cast_helpers.four_block.entrywise_variance(Y: ndarray, A: ndarray, r: int) → ndarray#

Entrywise variance of the imputed cells (paper equations 6a and 6c).

The residuals are formed on the observed entries of the sorted staircase (eq. 6a) and combined with the estimated subspaces to give, for each unobserved (i, t), the variance of the imputation error (eq. 6c).

Returns an (N, T) array of variances in the original row order; entries at observed cells are zero.

Notes

The mask is applied in the sorted frame, as eq. (6a) specifies. The authors’ CAST-panel package (1.0.2) masks sorted residuals with the unsorted indicator, which understates the standard errors whenever the sort permutation is non-trivial; see docs/replications/cast.rst.

mlsynth.utils.cast_helpers.four_block.four_block_est(Y: ndarray, N1: int, T1: int, r: int) → ndarray#

Algorithm 1 (FourBlockEst): impute the unobserved bottom-right block.

Parameters:

  • Y (np.ndarray) – (N, T) panel laid out as [[Ya, Yb], [Yc, ?]]; the bottom-right block is unobserved and its content is ignored.

  • N1, T1 (int) – Row/column size of the fully observed top-left block Ya.

  • r (int) – Rank of the underlying mean matrix.

Returns:

np.ndarray – (N - N1, T - T1) estimate of the missing block.

Notes

With a rank-r mean matrix and noiseless observations the output equals the missing block exactly (paper S3.1.1).

mlsynth.utils.cast_helpers.four_block.staircase_blocks(A: ndarray)#

Sort a staggered design into staircase form (paper Appendix A).

Parameters:

A (np.ndarray) – (N, T) indicator, 1 where the outcome is observed under control and 0 where it is an unobserved counterfactual.

Returns:

(order, inv, N_all, T_all, k) – order sorts units by ascending treatment duration (so the longest-observed units come first); inv restores the original row order. T_all holds the distinct control-run lengths in ascending order, N_all[i] the cumulative unit count down to the i-th longest, and k the number of distinct blocks.

CAST estimation pipeline (Xia, Yan & Wainwright 2025).

Selects the factor rank, imputes every counterfactual cell through the staircase of four-block sub-problems, attaches the entrywise variances, and aggregates them into per-period effects with standard errors.

mlsynth.utils.cast_helpers.pipeline.aggregate_effects(effects: ndarray, var: ndarray, treated: ndarray, time_labels: ndarray, weights: ndarray | None, renormalize: bool) → Dict[Any, Any]#

Per-period aggregate effect and its standard error.

The aggregate is the bilinear functional tau_w,t = sum_i w_i tau_i,t of paper eq. (4); its variance follows from the entrywise variances, treating the per-unit imputation errors within a period as independent. weights may be a length-N vector or an N x T matrix (a period-varying weight).

mlsynth.utils.cast_helpers.pipeline.run_cast(inputs, config) → CASTFit#

Run CAST end to end on prepared inputs.

mlsynth.utils.cast_helpers.pipeline.select_rank(Y: ndarray, A: ndarray, method: str = 'broken_stick', max_rank: int = 10) → int#

Choose the factor rank from the spectrum of the observed control block.

broken_stick reproduces the reference package’s default; svht is the Gavish-Donoho optimal hard threshold for a square-ish noisy matrix; ratio takes the largest consecutive singular-value gap.