Distribution-Regression Synthetic Control (DRSC)#

When to Use This Estimator#

Most synthetic control methods return one number per period: the treated unit’s outcome minus its synthetic counterpart. That is the right summary when the policy shifts everyone in the same direction. It is the wrong summary when the policy reshapes the outcome distribution – a minimum wage that lifts the bottom tail and leaves the median alone, a transfer that compresses inequality without moving the mean.

Distributional Synthetic Control (DSC) already answers part of that: it matches whole distributions instead of means, and returns quantile treatment effects. But it works on the unconditional distribution, the one you get by pooling everybody. If the outcome depends on individual characteristics – and wages depend heavily on education and experience – then pooling can hide the effect entirely. An increase that raises wages for young low-education workers and nobody else barely moves the overall wage distribution, because that group is a small slice of it.

DRSC conditions. It estimates the counterfactual distribution of the outcome given covariates, so you can ask what the policy did to a 25-year-old with a high-school education specifically, not to the workforce on average.

Reach for it when all of these hold: you have micro-data (many individuals per unit-period, not one aggregate number), one treated unit, and a reason to think the effect varies across people in a way you can describe with covariates.

Notation#

There are \(J + 1\) groups observed over \(T\) periods. Group \(i = 1\) is treated from period \(T_0 + 1\) onward; groups \(i = 2, \ldots, J+1\) are the donor pool. Within each group-period cell \((i, t)\) we observe \(n_{it}\) individuals, each with an outcome \(Y\) and a covariate vector \(x \in \mathbb{R}^{p}\).

The conditional distribution is modelled semiparametrically as

\[F_{it}(y \mid x) \;=\; \Lambda\bigl(x^{\top}\theta_{it}(y)\bigr),\]

where \(\Lambda\) is a known link (the normal CDF by default) and \(\theta_{it} : \mathcal{Y} \to \mathbb{R}^{p}\) is an unknown parameter function. This is the distribution regression of Foresi and Peracchi (1995). For a fixed \(y\) it is just a binary regression of \(\mathbf{1}\{Y \le y\}\) on \(x\), so estimation reduces to running one probit per outcome grid point per cell.

The object of interest is the pointwise difference between the observed and counterfactual conditional distributions,

\[\Delta_t(y \mid x) \;=\; F_{1,t}(y \mid x) - F^{0}_{1,t}(y \mid x),\]

and its integrated square over a region \(\mathcal{Y}_0\),

\[f_t(x, \mathcal{Y}_0) \;=\; \int_{\mathcal{Y}_0} \Delta_t(y \mid x)^2 \, dy .\]

Setting \(\mathcal{Y}_0\) to the whole support gives the total effect; restricting it to a policy-relevant window – the corridor between the old and new minimum wage – concentrates the test where the policy actually operates.

Assumptions#

Assumption 1 (Parallel trends in parameters). There exist weights \(w_2, \ldots, w_{J+1}\) summing to one such that, for all \(y \in \mathcal{Y}_0\) and all post-treatment \(t\),

\[\theta^{0}_{1,t}(y) - \theta_{1,T_0}(y) \;=\; \sum_{i=2}^{J+1} w_i \bigl[\theta_{i,t}(y) - \theta_{i,T_0}(y)\bigr].\]

Remark. This is the familiar parallel-trends idea moved into the parameter space of the model, not imposed on the outcome or the CDF. That choice does real work. The parameter function is linear in the common factors whereas the CDF is not, so the weights come out constant across the whole distribution, and the counterfactual is guaranteed to still be a valid distribution in the model class. Economically it says the treated unit’s exposure to aggregate shocks – skill-biased technical change, the business cycle – can be written as a weighted blend of the donors’ exposures.

Assumption 2 (Sampling). Within each cell the individuals are i.i.d., and cells are independent, with \(n_{it}/n \to r_{it} \in (0, \infty)\).

Remark. The asymptotics run in \(n\), the number of individuals, with the donor count and the number of periods held fixed. This is unlike every other estimator in mlsynth, where precision comes from a long panel. Here one pre-period is enough for consistency, and extra pre-periods buy precision in the weights, not in the distribution regressions. The practical consequence: worry about thin cells, not short panels.

Assumption 3 (Pre-treatment balance). The treated unit’s pre-treatment parameter function lies in the affine span of the donors’.

Remark. The functional analogue of the usual requirement that the treated unit sit inside the donor hull. It is only affine, not convex – see the note on negative weights below – and it is testable, which is what the pre-trend tests in the next section do.

Estimation#

The weights minimise the pre-treatment discrepancy, averaged over pre-periods and integrated over the outcome grid, subject only to adding up to one:

\[\widehat w \;=\; \operatorname*{arg\,min}_{\mathbf{1}'w = 1} \frac{1}{T_0 m} \sum_{t=1}^{T_0} \sum_{l=1}^{m} \bigl\| \widehat\theta_{1t}(y_l) - \sum_{i=2}^{J+1} w_i \widehat\theta_{it}(y_l) \bigr\|^{2} .\]

Because the objective is quadratic and the single constraint is linear, this has a closed form – no solver, no iteration:

\[\widehat w \;=\; \widehat G^{-1}\widehat c \;-\; \widehat G^{-1}\mathbf{1}\, \frac{\mathbf{1}'\widehat G^{-1}\widehat c - 1} {\mathbf{1}'\widehat G^{-1}\mathbf{1}} .\]

The grid \(\{y_l\}_{l=1}^{m}\) is a set of quantiles of the pooled outcome between grid_lo and grid_hi. One detail there matters before you change it: n_grid is the number of points requested, and ties are collapsed afterwards, so the active grid is generally shorter. Wages reported to the cent produce a lot of ties – on the New Jersey panel the paper’s 38 requested points leave 32 active. Passing the active count instead of the requested one is a silent mistake, not a loud one: asking for 32 leaves 28 active and moves Florida’s weight from 0.515 to 0.496.

Two consequences follow. Both look like faults, and neither is.

Negative weights are normal. Non-negativity is deliberately dropped (Doudchenko and Imbens), which lets the synthetic unit sit outside the donor hull and is what makes a good pre-treatment fit attainable when the treated unit is unusual. In the New Jersey application 20 of 42 donors receive negative weight.

The Gram matrix is ill-conditioned, and that requires float64. Donors’ parameter functions are similar to each other – which is precisely why a weighted combination can track the treated unit – so \(\widehat G\) is nearly singular. On the New Jersey panel its condition number is \(9.6 \times 10^{4}\), and the solve amplifies relative input error by roughly \(10^{6}\). Passing float32 data (relative error \(\sim 6 \times 10^{-8}\)) moves Florida’s weight from 0.515 to 0.505, and the largest donor weight by about 0.04; random perturbation of the same magnitude, which does not partially cancel the way rounding does, moves weights by 0.12 to 0.17. The estimand is far more stable than the weights, so the reported effects barely move, but the weights themselves are what people quote. Keep the inputs in float64. Setting ridge to a small positive value stabilises the solve at the cost of shrinking the weights toward \(1/J\).

Inference and diagnostics#

The null of no effect, \(H_0 : \Delta_t(\cdot \mid x) = 0\) on \(\mathcal{Y}_0\), is tested with the supremum statistic

\[T_n(x, t, \mathcal{Y}_0) \;=\; \sqrt{n} \, \sup_{y_l \in \mathcal{Y}_0} \bigl| \widehat\Delta_t(y_l \mid x) \bigr| ,\]

whose limit under the null is the supremum of a mean-zero Gaussian process. Critical values come from simulating that process using a plug-in estimate of its covariance kernel. The kernel carries two sources of uncertainty at the same \(\sqrt{n}\) rate: the distribution-regression error in the post period, and the weight-estimation error inherited from the pre-periods.

When the test rejects, a one-sided lower confidence bound for \(f_t(x)\) is reported. It is only meaningful conditional on rejection: under the null the delta-method derivative vanishes, the variance degenerates, and the bound collapses to zero.

With two or more pre-periods, parallel trends is testable directly by treating an earlier period as a pseudo-post period and asking whether the same machinery finds an effect where none can exist.

One reproducibility caveat before comparing numbers across machines. The simulated critical values require a matrix square root of the estimated kernel, obtained from an eigendecomposition. Eigenvector signs are implementation-defined, so a different LAPACK build produces a different factor and hence a different realisation from identical Gaussian draws, even with the seed fixed. The effect is a few tenths of a percent on the critical value. The test statistic itself involves no eigendecomposition and is exact.

Example#

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

# micro-panel: one row per (state, policy year, worker)
df = pd.read_parquet("benchmarks/reference/wied_nj_minwage/"
                     "nj_estimation_sample.parquet")
df["treat"] = ((df.STATEFIP == 34) & (df.t == 4)).astype(int)
df["x_std_sq"] = df["x_std"] ** 2

res = DRSC({
    "df": df, "outcome": "logwage", "treat": "treat",
    "unitid": "STATEFIP", "time": "t",
    "covariates": ["e_std", "x_std", "x_std_sq"],
    "evaluation_points": {
        "low_skill_young": {"e_std": -0.98, "x_std": -1.17,
                            "x_std_sq": 1.37},
        "high_skill":      {"e_std":  1.64, "x_std":  1.44,
                            "x_std_sq": 2.08},
    },
    "focus_region": (np.log(4.25), np.log(5.10)),   # the MW corridor
    "display_graphs": False,
}).fit()

e = res.conditional_effects["low_skill_young"]
e.f_hat, e.p_value, e.p_value_focused     # 1.71e-3, 0.054, 0.012
res.gram_condition_number                  # 9.6e4 -- see the float64 note
res.n_negative_weights                     # 20 of 42, expected

Verification#

Reproduced against Wied (2026), Tables 1 and 2, on the author’s own CPS estimation sample: all five donor weights and all twenty cells of Table 2 match to the paper’s printed precision, along with the active grid size (\(m = 32\)), the Gram condition number and the negative-weight count. See DRSC — New Jersey’s 1992 minimum wage (Wied 2026) and benchmarks/reference/wied_nj_minwage/.

Not to be confused with#

Distributionally Robust Synthetic Control (DROSC) is Distributionally Robust synthetic control – robustness of a conventional estimate to distributional shift, a different question entirely. Distributional Synthetic Control (DSC) is the unconditional distributional estimator DRSC extends.

Core API#

class mlsynth.DRSC(config: DRSCConfig | dict)#

Distribution-regression synthetic control.

Parameters:: config (DRSCConfig or dict) – See mlsynth.utils.drsc_helpers.config.DRSCConfig.

Examples

>>> from mlsynth import DRSC
>>> res = DRSC({
...     "df": micro_panel, "outcome": "logwage", "treat": "treat",
...     "unitid": "state", "time": "policy_year",
...     "covariates": ["educ_std", "exper_std", "exper_std_sq"],
...     "evaluation_points": {
...         "low_skill_young": {"educ_std": -0.98, "exper_std": -1.17,
...                             "exper_std_sq": 1.37}},
... }).fit()
>>> res.conditional_effects["low_skill_young"].p_value

fit() → DRSCResults#: Estimate the conditional distributional treatment effect.

class mlsynth.utils.drsc_helpers.config.DRSCConfig(*, 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>, covariates: ~typing.List[str], evaluation_points: ~typing.Dict[str, ~typing.Dict[str, float]], link: ~typing.Literal['probit', 'logit'] = 'probit', n_grid: ~typing.Annotated[int, ~annotated_types.Gt(gt=1)] = 38, grid_lo: ~typing.Annotated[float, ~annotated_types.Gt(gt=0.0), ~annotated_types.Lt(lt=1.0)] = 0.1, grid_hi: ~typing.Annotated[float, ~annotated_types.Gt(gt=0.0), ~annotated_types.Lt(lt=1.0)] = 0.9, focus_region: ~typing.Tuple[float, float] | None = None, ridge: ~typing.Annotated[float, ~annotated_types.Ge(ge=0)] = 0.0, alpha: ~typing.Annotated[float, ~annotated_types.Gt(gt=0.0), ~annotated_types.Lt(lt=1.0)] = 0.1, n_gp_draws: ~typing.Annotated[int, ~annotated_types.Gt(gt=0)] = 10000, seed: int = 1992)#

Typed configuration for mlsynth.DRSC.

model_config: ClassVar[ConfigDict] = {'arbitrary_types_allowed': True, 'extra': 'forbid'}#: Configuration for the model, should be a dictionary conforming to [ConfigDict][pydantic.config.ConfigDict].

class mlsynth.utils.drsc_helpers.structures.DRSCResults(*, 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, conditional_effects: ~typing.Dict[str, ~mlsynth.utils.drsc_helpers.structures.ConditionalEffect] = <factory>, outcome_grid: ~numpy.ndarray | None = None, gram_condition_number: float | None = None, n_negative_weights: int | None = None, metadata: ~typing.Dict[str, ~typing.Any] = <factory>)#

Container returned by mlsynth.DRSC.fit().

An EffectResult: the headline scalar is lifted into the standardized sub-models so the flat accessors resolve through the base contract, while the conditional objects – which have no analogue in a scalar-ATT estimator – stay on conditional_effects.

effects.att is the integrated squared effect at the first evaluation point. That is a summary of convenience: this estimator’s output is a function of the covariate value, and reading only att discards what it was built to show.

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].

class mlsynth.utils.drsc_helpers.structures.ConditionalEffect(label: str, x: ndarray, delta: ndarray, f_hat: float, std_error: float, lower_bound: float, t_stat: float, crit_value: float, p_value: float, t_stat_focused: float = nan, crit_value_focused: float = nan, p_value_focused: float = nan, lower_bound_focused: float = nan)#

The estimate at one evaluation point, for one post-treatment period.

label#

The caller’s name for the evaluation point.

Type:: str

x#

The design vector, intercept first.

Type:: np.ndarray

delta#

Pointwise CDF difference over the outcome grid (eq. 3).

Type:: np.ndarray

f_hat#

Integrated squared effect (eq. 2), the scalar summary.

Type:: float

std_error#

Standard error of f_hat.

Type:: float

lower_bound#

One-sided lower confidence bound (eq. 10); 0 when uninformative.

Type:: float

t_stat#

Supremum statistic on the full grid (eq. 14).

Type:: float

crit_value#

Simulated critical value at the configured level.

Type:: float

p_value#

Full-support p-value.

Type:: float

t_stat_focused, crit_value_focused, p_value_focused, lower_bound_focused: The same quantities restricted to focus_region; NaN when no focus region was given.

Distribution-Regression Synthetic Control (DRSC)

Contents