Synthetic Regressing Control (SRC)#

When to use SRC – and when not to#

Synthetic Regressing Control (Zhu 2023) is for the case where the ordinary synthetic control fits the treated unit’s pre-treatment path poorly. The canonical synthetic control restricts the donor weights to the unit simplex (non-negative, summing to one), which keeps the counterfactual inside the convex hull of the donors but cannot track a treated unit that sits near or outside that hull – the interpolation-bias / imperfect-pre-fit problem. SRC relaxes that restriction in a disciplined way and is a good choice when:

You have a single treated unit and a donor pool with real heterogeneity, and a simplex synthetic control leaves a visible pre-treatment gap.
You are willing to let the counterfactual extrapolate a little beyond the donor hull, but want the amount of extrapolation regularised rather than unbounded (as an unconstrained regression would give).
You want a deterministic, reproducible estimate. By default SRC’s weights are pinned by the risk criterion, not by an Abadie predictor-weight (V) search, so there is no optimiser-seed dependence. (The paper’s covariate + V-search variant is available as a seeded opt-in; see below.)

Do not use SRC when:

A simplex synthetic control already fits well and you want an interpretable convex weighting. Use canonical SCM (Vanilla Synthetic Control (VanillaSC)) or Two-Step Synthetic Control; SRC’s combined coefficients can be negative and do not sum to one.
You want an explicit bias–variance trade-off between matching and synthetic control chosen by cross-validation. That is Matching and Synthetic Control (MASC).
You want denoising of a low-rank donor matrix before weighting. Use Cluster Synthetic Controls (CLUSTERSC) (robust PCA / PCR) or Synthetic Nearest Neighbors / Causal Matrix Completion (SNN).
You need posterior uncertainty on the weights. Use Bayesian Synthetic Control with a Soft Simplex Constraint (BVS-SS).
There are multiple treated units, spillovers, or a continuous dose – SRC encodes a single binary intervention on one unit.

Notation#

We use the synthetic-control canon. Let \(\mathcal{N} \coloneqq \{1, \ldots, N\}\) index the units, with \(j = 1\) the treated unit and \(\mathcal{N}_0 \coloneqq \mathcal{N} \setminus \{1\} = \{2, \ldots, N\}\) the donor pool of \(N_0 \coloneqq N - 1\) controls. Time is \(t \in \mathcal{T} \coloneqq \{1, \ldots, T\}\); the intervention takes effect after \(T_0\), so the pre-period is \(\mathcal{T}_1 \coloneqq \{t \in \mathcal{T} : t \le T_0\}\) with \(|\mathcal{T}_1| = T_0\) and the post-period is \(\mathcal{T}_2 \coloneqq \{t \in \mathcal{T} : t > T_0\}\) with length \(T_2 \coloneqq |\mathcal{T}_2|\).

Notation bridge. Zhu (2023) indexes the treated unit as one of \(j = 1, \ldots, J + 1\) with \(J\) controls, and writes the pre-period as \(\mathcal{T}_0\); his \(J\) is our \(N_0\), and his \(\mathcal{T}_0\) is our \(\mathcal{T}_1\). He names the method Synthetic Regressing Control (SRC); mlsynth ships it as SRC.

Panel objects. Write the treated unit’s pre-period outcome path as \(\mathbf{y}_1 \coloneqq (y_{1t})_{t \in \mathcal{T}_1} \in \mathbb{R}^{T_0}\) and stack the donors’ pre-period paths columnwise into the donor matrix \(\mathbf{Y}_0 \coloneqq [\mathbf{y}_j]_{j \in \mathcal{N}_0} \in \mathbb{R}^{T_0 \times N_0}\), one column \(\mathbf{y}_j = (y_{jt})_{t \in \mathcal{T}_1} \in \mathbb{R}^{T_0}\) per donor. These \(n \coloneqq T_0\) matching rows are the pre-period outcomes in the default specification; the covariate variant (below) stacks standardized covariate rows on top, giving an \(n \times N_0\) matching matrix with \(n > T_0\), and every formula below carries over with \(T_0\) replaced by that row count \(n\). Centering (demeaning) uses the projector \(\mathbf{Q} \coloneqq \mathbf{I}_{T_0} - T_0^{-1} \mathbf{1}\mathbf{1}^{\top}\), and \(\bar{y}_j \coloneqq T_0^{-1} \sum_{t \in \mathcal{T}_1} y_{jt}\) is donor \(j\)’s pre-period mean, so \(\mathbf{Q}\mathbf{y}_j = \mathbf{y}_j - \bar{y}_j \mathbf{1}\) is the centered column.

Unit regression. SRC first rescales each donor to the treated unit by a separate demeaned univariate least-squares fit – Zhu’s “unit regression”:

\[\widehat{\theta}_j \coloneqq \frac{(\mathbf{y}_j - \bar{y}_j \mathbf{1})^{\top} (\mathbf{y}_1 - \bar{y}_1 \mathbf{1})} {\lVert \mathbf{y}_j - \bar{y}_j \mathbf{1} \rVert_2^2}, \qquad j \in \mathcal{N}_0.\]

The regressed (matched) control for donor \(j\) is the centered, rescaled column \(\widehat{\theta}_j (\mathbf{y}_j - \bar{y}_j \mathbf{1})\); because \(\widehat{\theta}_j\) is unconstrained it can already extrapolate. Collect these columnwise into the regressed-control matrix

\[\widetilde{\mathbf{Y}}_0 \coloneqq \bigl[\, \widehat{\theta}_j (\mathbf{y}_j - \bar{y}_j \mathbf{1}) \,\bigr]_{j \in \mathcal{N}_0} = (\mathbf{Q}\mathbf{Y}_0)\operatorname{diag}(\widehat{\boldsymbol{\theta}}) \in \mathbb{R}^{T_0 \times N_0},\]

whose column \(j\) is donor \(j\)’s matched control (the object the previous version of this page wrote \(\mathbf{E}\)).

Synthesis. The matched controls are combined by box weights \(\mathbf{w} \in \mathcal{C}_{\mathrm{box}} \coloneqq [0, 1]^{N_0} = \{\mathbf{w} \in \mathbb{R}^{N_0} : 0 \le w_j \le 1\}\) – note there is no sum-to-one constraint, which is what separates SRC from a simplex synthetic control (Zhu writes this set \(\mathcal{H}_J\)). With the level \(\widehat{b}_0 \coloneqq \bar{y}_1 - \sum_{j \in \mathcal{N}_0} \widehat{w}_j \widehat{\theta}_j \bar{y}_j\) that restores the mean removed by centering, the SRC counterfactual is

\[\widehat{y}_{1t} = \bar{y}_1 + \sum_{j \in \mathcal{N}_0} w_j\, \widehat{\theta}_j\, (y_{jt} - \bar{y}_j) = \widehat{b}_0 + \sum_{j \in \mathcal{N}_0} w_j \widehat{\theta}_j\, y_{jt}, \qquad t \in \mathcal{T},\]

which estimates the no-intervention outcome \(y_{1t}^N\). The effective coefficient on donor \(j\) is \(\widehat{\theta}_j w_j\), which may be negative. The per-period effect is the scalar \(\tau_t \coloneqq y_{1t} - \widehat{y}_{1t}\), and the ATT is \(\widehat{\tau} \coloneqq T_2^{-1} \sum_{t \in \mathcal{T}_2} \tau_t\).

Weights. The box weights minimise a Mallows \(C_p\) criterion – an unbiased estimator of the synthetic estimator’s squared risk – written in the canon’s loss-plus-penalty form,

\[\widehat{\mathbf{w}} \in \operatorname*{argmin}_{\mathbf{w} \in \mathcal{C}_{\mathrm{box}}} \; \underbrace{\bigl\lVert \mathbf{Q}\mathbf{y}_1 - \widetilde{\mathbf{Y}}_0 \mathbf{w} \bigr\rVert_2^2 }_{\mathcal{L}(\mathbf{w})} + \underbrace{2\,\widehat{\sigma}^2 \lVert \mathbf{w} \rVert_1 }_{\mathcal{P}(\mathbf{w})},\]

where \(\lVert \mathbf{w} \rVert_1 = \sum_{j \in \mathcal{N}_0} w_j\) because \(\mathbf{w} \ge \mathbf{0}\) on the box. The noise level \(\widehat{\sigma}^2\) is the plug-in residual variance of the full (multivariate) demeaned regression of the treated path on the donor matrix,

\[\widehat{\sigma}^2 \coloneqq \frac{1}{\max(T_0 - N_0,\, 1)} \bigl\lVert \mathbf{Q}\mathbf{y}_1 - (\mathbf{Q}\mathbf{Y}_0)\,\widehat{\boldsymbol{\beta}} \bigr\rVert_2^2, \qquad \widehat{\boldsymbol{\beta}} \coloneqq (\mathbf{Q}\mathbf{Y}_0)^{+}\, \mathbf{Q}\mathbf{y}_1,\]

with \((\cdot)^{+}\) the Moore–Penrose pseudoinverse and the degrees of freedom floored at one when \(N_0 \ge T_0\) (a minimum-norm fit; a small ridge then keeps the synthesis QP strictly convex). This \(\widehat{\sigma}^2\) estimates the working-model noise level \(\sigma^2\) and enters the penalty as the per-donor complexity charge \(2\widehat{\sigma}^2\), which shrinks weight toward zero unless a matched control earns its place. That charge – not a sum-to-one constraint – is what identifies \(\widehat{\mathbf{w}}\).

A note on the paper’s display formula. Zhu’s eq. 19 writes \(\widehat{\sigma}^2 = \lVert \mathbf{Q}\mathbf{y}_1 - \mathbf{Q}\mathbf{Y}_0[\operatorname{diag}(\mathbf{Y}_0^{\top} \mathbf{Q}\mathbf{Y}_0)]^{-1}\mathbf{Y}_0^{\top}\mathbf{Q}\mathbf{y}_1 \rVert_2^2\), which reduces to the residual of the equal-weight sum of all \(N_0\) univariate fits, \(\lVert \mathbf{Q}\mathbf{y}_1 - \widetilde{\mathbf{Y}}_0\mathbf{1} \rVert_2^2\). Because that sum stacks \(N_0\) controls that each already track the treated unit, it overshoots by roughly a factor of \(N_0\), so the quantity scales with \(N_0^2\), not with the noise level: in simulations with a known \(\sigma^2\) it over-estimates by two orders of magnitude, and the resulting penalty is large enough to drive every weight to zero (a degenerate flat-at-the-pre-mean counterfactual) in both the paper’s donors-fewer-than-periods regime and out of it. The author’s reference Code_SMC.R does not use eq. 19; it uses the full-regression residual variance given above, which recovers the true \(\sigma^2\) to within a few percent. mlsynth follows the reference implementation (cross-validated to \(< 2\text{e-}13\); see Verification), and reads the printed eq. 19 as an uncorrected typo rather than the estimator the paper’s own non-degenerate Basque results were produced with.

Assumptions#

Single treated unit, balanced panel. One unit is treated after \(T_0\); every unit is observed at every period.

Remark. The setup boundary (mlsynth.utils.src_helpers.setup) enforces both, translating a violation to MlsynthDataError.
At least two pre-treatment periods. The univariate matches and the risk criterion need pre-period variation; the plug-in \(\widehat{\sigma}^2\) uses the full-model residual, so more pre-periods than donors gives it degrees of freedom.

Remark. When donors are not fewer than matching rows (\(N_0 \ge T_0\), or \(N_0 \ge n\) with stacked covariates) the residual degrees of freedom \(\max(T_0 - N_0, 1)\) are floored at one and \(\widehat{\boldsymbol{\beta}}\) is the minimum-norm (pseudoinverse) fit; a small ridge keeps the synthesis QP strictly convex. For a high-dimensional donor pool, screen first (screen="sirs", below) or use Cluster Synthetic Controls (CLUSTERSC).
Regularised extrapolation is acceptable. The combined coefficients can leave the simplex, so the counterfactual can extrapolate; the \(C_p\) penalty bounds how much.

Remark. If extrapolation is unacceptable for the application (you need a convex, interpretable weighting), SRC is the wrong tool – use a simplex SCM.

Screening a wide donor pool#

When the pool is not small relative to the pre-period (\(N_0 \ge T_0\), and the paper recommends screening already at \(N_0 \ge 4T_0/5\)), the unit regression is ill-posed: the centered pre-period system \(\mathbf{Q} \mathbf{y}_1 \approx \mathbf{Q}\mathbf{Y}_0\boldsymbol{\beta}\) has more donor coefficients than periods, so the full-model residual collapses, \(\widehat{\sigma}^2 \to 0\), and the \(C_p\) penalty switches off – the box fit then interpolates the pre-period with no complexity charge. The paper’s Algorithm 2 handles this by screening the pool before Algorithm 1. Set screen="sirs":

import pandas as pd
from mlsynth import SRC

df = pd.read_csv("basedata/basque_data.csv")
df["treat"] = ((df["regionname"] == "Basque Country (Pais Vasco)")
               & (df["year"] >= 1975)).astype(int)
res = SRC({"df": df, "outcome": "gdpcap", "treat": "treat",
           "unitid": "regionname", "time": "year",
           "screen": "sirs", "display_graphs": False}).fit()

Donors are ranked by the sure-independent-ranking-and-screening (SIRS) marginal utility of Zhu, Fan, Li & Zhu (2011) on the pre-period outcomes – a model-free measure of dependence between each donor and the treated path – and the top \(k \coloneqq \min\bigl(\lfloor T_0/\log(T_0/2)\rfloor,\, T_0 - 1\bigr)\) are kept (override with n_screen). The \(T_0 - 1\) cap guarantees the screened pool is smaller than the pre-period, so the fit is well posed and \(\widehat{\sigma}^2 > 0\). On the 1975-treated Basque panel (\(T_0 = 20\), \(N_0 = 16\)) screening keeps \(k = 8\) donors and lifts \(\widehat{\sigma}^2\) by roughly two orders of magnitude, so the penalty is active rather than switched off.

A shape-based alternative: screen="fpca". SIRS ranks donors by marginal dependence; screen="fpca" instead selects an economically coherent peer group by clustering the pre-period trajectories and keeping the donors in the treated unit’s cluster. It reuses Cluster Synthetic Controls (CLUSTERSC)’s own donor-clustering routines – functional PCA (cubic-B-spline smoothing then a PCA truncated at 95% variance) followed by silhouette-chosen k-means (Bayani 2021, compute_fpca_features() and assign_clusters()) – so the two estimators screen donors by the same principle. The cluster size is data-driven (no n_screen), and the selection is deterministic. On the 1975-treated Basque panel it keeps the developed / service regions (Baleares, Cataluña, Madrid) in the Basque Country’s cluster, where SIRS instead ranks in Aragón; the two rules give genuinely different pools, and hence different counterfactuals – a reminder that donor selection is a modelling choice, not an objective one.

Two cautions. Screening is a well-posedness and reproducibility device, not a route to the paper’s Basque weight cells: even a screened pool that contains Cantabria still lands on the same Rioja/Madrid decomposition, because the remaining freedom lives on the non-identified \(\mathbf{V}\) manifold (see below), not in the pool. And the SIRS statistic is scale-free by construction (donor paths are standardised); the paper’s printed Eq. 22 carries the inner term as \(Y_{jt}\), which we read as the cited-form \(Y_{j\ell}\).

Inference and diagnostics#

SRC returns a point estimate. The primary diagnostic is the pre-treatment fit (res.fit_diagnostics.rmse_pre): SRC’s reason for being is a tighter pre-fit than a simplex SCM, so compare the two. The combined coefficients res.weights_vector (which may be negative) and the box weights res.box_weights are both exposed; a donor with a hard zero box weight is genuinely dropped. For placebo inference, run the estimator on each control region in turn and compare the treated gap to the placebo distribution (the fit is deterministic, so placebos are exact).

Example#

import pandas as pd
from mlsynth import SRC

df = pd.read_csv("basedata/basque_data.csv")
df["treat"] = ((df["regionname"] == "Basque Country (Pais Vasco)")
               & (df["year"] >= 1970)).astype(int)

res = SRC({
    "df": df, "outcome": "gdpcap", "treat": "treat",
    "unitid": "regionname", "time": "year",
    "display_graphs": False,
}).fit()

print(f"pre-period RMSE : {res.fit_diagnostics.rmse_pre:.3f}")
print(f"mean post ATT   : {res.att:+.3f}")
for u, c in sorted(res.donor_weights.items(), key=lambda kv: -abs(kv[1]))[:3]:
    print(f"  {u:<28s} {c:+.3f}")

This prints:

pre-period RMSE : 0.048
mean post ATT   : -0.858
  Murcia (Region de)           +0.625
  Madrid (Comunidad De)        +0.370
  Castilla Y Leon              +0.242

The Basque Country tracks its SRC counterfactual to a pre-period RMSE of about 0.048 (thousand-1986-USD per capita), then diverges steadily – the familiar Abadie-Gardeazabal shape of the economic cost of ETA terrorism – recovered here by the box-weighted matched controls rather than a simplex.

Note that this deterministic Algorithm 1 is not the specification behind the paper’s Basque table: it selects Murcia / Madrid / Castilla y León and a mean ATT of about \(-0.86\), whereas the paper (Table 5 / Figure 1) reports Rioja / Madrid / Cantabria and an ATT reaching \(\approx -1.9\). The difference is the predictor-weight layer, covered next.

Reproducing the paper’s covariate specification#

The paper’s Basque result uses Algorithm 3: covariate matching with an Abadie predictor-weight optimisation over the matching rows. Stack the standardized covariate rows on top of the outcome rows to form the \(n \times N_0\) matching matrix, and let \(\mathbf{V} \coloneqq \operatorname{diag}(v_1, \ldots, v_n)\) with \(v_i \ge 0\) be the diagonal predictor weights on those rows. Every quantity above is then computed on the \(\mathbf{V}^{1/2}\)- scaled rows – the treated matching vector becomes \(\mathbf{V}^{1/2} \mathbf{y}_1\) and the donor matrix \(\mathbf{V}^{1/2}\mathbf{Y}_0\) before centering – so \(\mathbf{V}\) sets each matching row’s relative importance in the unit regression and the \(C_p\) synthesis. Choosing \(\mathbf{V}\) is an explicit, seeded opt-in in mlsynth – covariates with per-covariate covariate_windows, the fit_window for the outcome rows, and v_search="de":

res = SRC({
    "df": df, "outcome": "gdpcap", "treat": "treat",
    "unitid": "regionname", "time": "year",
    "covariates": ["school.illit", "school.prim", "school.med", "school.high",
                   "invest", "sec.agriculture", "sec.energy", "sec.industry",
                   "sec.construction", "sec.services.venta",
                   "sec.services.nonventa", "popdens"],
    "covariate_windows": {"invest": (1964, 1969), "sec.industry": (1961, 1969)},  # ...
    "fit_window": (1960, 1969),          # Abadie's time.optimize.ssr
    "v_search": "de", "v_seed": 0,       # seeded global V search
    "display_graphs": False,
}).fit()

This recovers the paper’s donor structure – Rioja dominant, then Madrid – and the Figure 1 ATT magnitude (\(\approx -1.4\) mean, \(-2.4\) by 1997).

A caveat stated plainly. The \(\mathbf{V}\) optimum is not identified: a manifold of \(\mathbf{V}\) fits the pre-period equally well while disagreeing out of sample, so the exact split among the top donors is seed-dependent (a global search reproduces the paper’s average placebo MSPE but not its per-region cells). The search is seeded so a given call is reproducible, and it is opt-in so the identified Algorithm 1 stays the default. Read the v_search weights as “the paper’s specification, one draw from its non-identified \(\mathbf{V}\) manifold”, not as identified quantities. When you want a single reproducible weighting, use the default.

Verification#

The SRC weight computation is cross-validated against the author’s reference R implementation (Code_SMC.R): on the identical Basque matching matrix, the per-donor \(\widehat{\theta}_j\), the box weights, the combined coefficients, the intercept \(\widehat{b}_0\) (bias) and \(\widehat{\sigma}^2\) all match the reference SMCV (whose synthesis QP is solved by solve.QP) to < 2e-13. The active-set box QP in mlsynth.utils.src_helpers.solver reproduces solve.QP to machine precision while pinning the box bounds exactly. With the paper’s covariate specification and v_search="de" the estimator additionally reproduces the Table 5 / Figure 1 donor structure and ATT magnitude (subject to the \(\mathbf{V}\) non-identification noted above). See the replication page SRC — Synthetic Regressing Control (Zhu 2023) and the durable case benchmarks/cases/src_basque.py. The solver, weight computation, the covariate / V-search paths, setup and result contract are unit-tested (mlsynth/tests/test_src.py, full coverage).

Core API#

Synthetic Regressing Control (SRC) estimator.

A thin, NumPy-first orchestration over mlsynth.utils.src_helpers. SRC (Zhu 2023) addresses the imperfect-pre-fit problem of the simplex synthetic control in two deterministic steps. First, unit matching: each donor is rescaled to the treated unit by a univariate OLS, theta_j, giving a matched control theta_j x_j that can already extrapolate (theta_j is unconstrained). Second, synthesis: the matched controls are combined with box-[0, 1] weights w chosen to minimise a Mallows/Cp unbiased-risk criterion; the combined donor coefficient is theta_j w_j. The Cp penalty on sum(w) – not an Abadie predictor-weight search – identifies the weights, so the estimator is reproducible.

The synthesis QP is solved exactly by an active-set box solver (solve_box_qp()) that matches the reference R solve.QP to machine precision. Optional covariates add predictor rows (Algorithm 3) at equal predictor weight.

References

Zhu, Rong J. B. (2023). Synthetic Regressing Control. arXiv:2306.02584. https://arxiv.org/abs/2306.02584

class mlsynth.estimators.src.SRC(config: SRCConfig | dict)#

Bases: object

Synthetic Regressing Control estimator.

Parameters:: config (SRCConfig or dict) – Validated configuration. Beyond the common fields (df, outcome, treat, unitid, time, display_graphs, save, colours), SRC reads ridge (the Cp box-QP stabiliser) and covariates (optional predictor columns for Algorithm 3).

fit() → SRCResults#

Run SRC end to end and return SRCResults.

Returns:

SRCResults – Container with the ATT, the counterfactual and gap paths, the combined donor coefficients (theta * w), the box weights, the pre-period RMSE, and the standardized sub-models.

Raises:

MlsynthDataError – If the panel violates SRC’s requirements (single treated unit, balanced panel, at least two pre-treatment periods, a donor pool).
MlsynthEstimationError – If the weight computation fails at runtime.
MlsynthPlottingError – If plotting raises when display_graphs=True.

Configuration#

class mlsynth.config_models.SRCConfig(*, 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>, ridge: ~typing.Annotated[float, ~annotated_types.Gt(gt=0)] = 0.001, covariates: ~typing.List[str] | None = None, covariate_windows: ~typing.Dict[str, ~typing.Any] | None = None, fit_window: ~typing.Tuple[~typing.Any, ~typing.Any] | None = None, screen: ~typing.Literal['none', 'sirs', 'fpca'] = 'none', n_screen: ~typing.Annotated[int | None, ~annotated_types.Ge(ge=1)] = None, v_search: ~typing.Literal['none', 'de'] = 'none', v_seed: int = 0, v_maxiter: ~typing.Annotated[int, ~annotated_types.Ge(ge=1)] = 60, v_popsize: ~typing.Annotated[int, ~annotated_types.Ge(ge=2)] = 12)#

Configuration for the Synthetic Regressing Control estimator (Zhu 2023).

SRC matches each donor to the treated unit by a univariate regression, then synthesises the matched controls with box-[0, 1] weights chosen by a Mallows/Cp unbiased-risk criterion. The combined donor coefficients may be negative (controlled extrapolation), and the Cp penalty – not a predictor (V) search – identifies the weights, so the default estimator is deterministic.

The covariate + predictor-weight (V) variant of the paper’s Basque application (Algorithm 3 / Table 5) is available as an explicit, seeded opt-in (covariates=[...] with v_search="de"); see v_search.