Synthetic Difference-in-Differences

46.1 When synthetic control can't fit (40 min)

Synthetic control's weights are non-negative and sum to 1, so synthetic R07 always lies inside the range of the donors. If R07 is the biggest region (or the fastest-growing, or the most seasonal), no mix of smaller regions can match it. Toronto, a flagship store, a market leader: the units companies most often run campaigns in are the hardest to match.

Two fixes:

46.2 The biggest region (50 min)

regions(treated_size=250) is the same world, except that R07 sells $250k a week more than the median region. Create day46/big_region.py:

# day46/big_region.py   (python -m day46.big_region)
"""
Synthetic control's weakness: weights are >= 0 and sum to 1, so the synthetic region can never be
BIGGER than the biggest donor. Same world, but now the campaign ran in a region that sells 250 $k a week
more than the median region. Adding an intercept (matching changes, not levels) fixes it.
"""
import matplotlib.pyplot as plt
import polars as pl

from causalkit.regions import public, regions, truth
from causalkit.synth import synthetic_control

full = regions(treated_size=250)
df = public(full)
t = truth(full)
plain = synthetic_control(df, "R07", 80)
demeaned = synthetic_control(df, "R07", 80, intercept=True)
print(f"R07 sells {plain['actual'][:79].mean():.0f} $k a week; the biggest donor sells "
      f"{df.filter(pl.col('region') != 'R07').group_by('region').agg(pl.col('sales').mean())['sales'].max():.0f}")
for name, sc in [("synthetic control", plain), ("synthetic control + intercept", demeaned)]:
    print(f"{name:<31} effect {sc['effect']:6.1f}   pre RMSPE {sc['rmspe_pre']:5.1f}   weights {sc['weights']}")
print(f"Truth: {t['avg_post_lift']:.1f}")

fig, ax = plt.subplots(figsize=(9, 4))
ax.plot(plain["weeks"], plain["actual"], color="#dc2626", label="R07")
ax.plot(plain["weeks"], plain["synthetic"], color="#64748b", ls="--", label="synthetic control")
ax.plot(demeaned["weeks"], demeaned["synthetic"], color="#2563eb", ls="--", label="with intercept")
ax.axvline(80, color="grey", ls="--", lw=0.8)
ax.legend()
ax.set_xlabel("week")
ax.set_ylabel("sales ($k per week)")
fig.tight_layout()
fig.savefig("figures/day46_big_region.png", dpi=120)
print("Saved figures/day46_big_region.png")

python -m day46.big_region

R07 sells 764 $k a week; the biggest donor sells 604
synthetic control               effect  214.8   pre RMSPE 164.9   weights {'R25': 0.987}
synthetic control + intercept   effect   17.3   pre RMSPE   3.9   weights {'R14': 0.211, 'R18': 0.21, 'R13': 0.209, 'R05': 0.171, 'R22': 0.169, 'R16': 0.015, 'R23': 0.012}
Truth: 14.3

Plain synthetic control puts all its weight on the biggest donor and still misses by $165k a week before the campaign. Its "effect" of $215k is just the size gap. The warning is right there in the pre-period RMSPE: never report a synthetic control without its pre-period fit.

With an intercept, the weights and the estimate are exactly Day 45's: the level difference is absorbed, and the donors are chosen on how R07 moves.

46.3 Synthetic DiD (65 min)

Concept. Synthetic DiD (SDID) combines three ideas:

1. Unit weights, like synthetic control, but with an intercept (levels don't need to match) and a ridge penalty that spreads weight over more donors (less dependence on any one).
2. Time weights: choose pre-period weeks whose donor sales best predict the donors' post-period average. Weeks that look like the post-period count more.
3. A weighted DiD: R07's change minus the weighted donors' change, each from the weighted pre-period to the post-period.

It's DiD made robust to non-parallel trends, or synthetic control made robust to level differences, depending on where you start. Type sdid into synth.py. Create day46/sdid.py:

# day46/sdid.py   (python -m day46.sdid)
"""
Synthetic difference-in-differences (Arkhangelsky, Athey, Hirshberg, Imbens & Wager, 2021):
  unit weights: like synthetic control, but with an intercept and a penalty that spreads weight;
  time weights: the pre-period weeks that best resemble the post-period (for the donors);
  then a weighted DiD. It matches trends, not levels, and borrows strength from more donors.
"""
from causalkit.regions import public, regions, truth
from causalkit.synth import did, sdid, synthetic_control

for title, size in [("R07 a typical-size region", 0), ("R07 a very large region", 250)]:
    full = regions(treated_size=size)
    df = public(full)
    s = sdid(df, "R07", 80)
    print(f"{title}:   truth {truth(full)['avg_post_lift']:.1f}")
    print(f"   DiD {did(df, 'R07', 80)['effect']:6.1f}   SC {synthetic_control(df, 'R07', 80)['effect']:6.1f}   "
          f"SDID {s['effect']:6.1f}")
print(f"\nSDID unit weights ({len(s['unit_weights'])} donors):", dict(list(s["unit_weights"].items())[:6]), "...")
print("SDID time weights:", s["time_weights"])

python -m day46.sdid

R07 a typical-size region:   truth 14.3
   DiD   32.1   SC   17.3   SDID   20.9
R07 a very large region:   truth 14.3
   DiD   32.1   SC  214.8   SDID   20.9

SDID unit weights (22 donors): {'R14': 0.094, 'R13': 0.089, 'R18': 0.089, 'R22': 0.084, 'R16': 0.062, 'R05': 0.061} ...
SDID time weights: {49: 0.172, 51: 0.01, 52: 0.056, 55: 0.106, 64: 0.014, 66: 0.059, 74: 0.155, 78: 0.11, 79: 0.319}

SDID is unaffected by R07's size. It spreads weight over 22 donors (synthetic control used 6), and it puts most time weight on recent weeks (78, 79) and on weeks 49–55, a year before the post-period, so the same season.

On this dataset SDID ($20.9k) is further from the truth than synthetic control ($17.3k). One dataset can't tell you which estimator is better: that's the mini-project.

46.4 What do the time weights do? (30 min)

Try it. Rerun SDID over the 20 datasets of the mini-project (regions(seed=100 + s)) with and without time weights (without: compare the post-period with the plain pre-period average). Create day46/no_time_weights.py. You need _split and simplex_ls from synth.py.

# day46/no_time_weights.py   (python -m day46.no_time_weights)
import numpy as np

from causalkit.regions import public, regions, truth
from causalkit.synth import _split, sdid, simplex_ls

errs = {"SDID": [], "SDID without time weights": []}
for s in range(20):
    full = regions(seed=100 + s)
    df, tr = public(full), truth(full)["avg_post_lift"]
    errs["SDID"].append(sdid(df, "R07", 80)["effect"] - tr)
    weeks, donors, y1, Y0, pre = _split(df, "R07", 80)
    zeta = (~pre).sum() ** 0.25 * np.diff(Y0[pre], axis=0).std(ddof=1)
    omega, _ = simplex_ls(Y0[pre], y1[pre], intercept=True, ridge=zeta**2 * pre.sum())
    flat = (y1[~pre].mean() - y1[pre].mean()) - (Y0[~pre].mean(0) - Y0[pre].mean(0)) @ omega
    errs["SDID without time weights"].append(flat - tr)
for k, e in errs.items():
    print(f"{k:<28} average error {np.mean(e):+.1f}   SD {np.std(e):.1f}")

SDID                         average error +2.2   SD 4.7
SDID without time weights    average error +4.3   SD 3.8

46.5 Day 46 mini-project: four estimators, two worlds (65 min)

Create day46/monte_carlo.py:

# day46/monte_carlo.py   (python -m day46.monte_carlo)
"""
Mini-project: four estimators, 20 simulated datasets, two worlds. Average error and its spread.
"""
import numpy as np
import polars as pl

from causalkit.regions import public, regions, truth
from causalkit.synth import did, sdid, synthetic_control

estimators = {
    "DiD": lambda d: did(d, "R07", 80)["effect"],
    "synthetic control": lambda d: synthetic_control(d, "R07", 80)["effect"],
    "SC + intercept": lambda d: synthetic_control(d, "R07", 80, intercept=True)["effect"],
    "SDID": lambda d: sdid(d, "R07", 80)["effect"],
}
rows = []
for world, size in [("typical-size region", 0), ("very large region", 250)]:
    errors = {k: [] for k in estimators}
    for s in range(20):
        full = regions(seed=100 + s, treated_size=size)
        df, tr = public(full), truth(full)["avg_post_lift"]
        for k, f in estimators.items():
            errors[k].append(f(df) - tr)
    for k, e in errors.items():
        rows.append({"world": world, "estimator": k, "average error": round(np.mean(e), 1),
                     "SD of error": round(np.std(e), 1)})
with pl.Config(tbl_rows=10, tbl_width_chars=100):
    print(pl.DataFrame(rows))

python -m day46.monte_carlo

│ world               ┆ estimator         ┆ average error ┆ SD of error │
│ typical-size region ┆ DiD               ┆ 7.9           ┆ 6.1         │
│ typical-size region ┆ synthetic control ┆ 1.0           ┆ 3.8         │
│ typical-size region ┆ SC + intercept    ┆ 2.5           ┆ 3.0         │
│ typical-size region ┆ SDID              ┆ 2.2           ┆ 4.7         │
│ very large region   ┆ DiD               ┆ 7.9           ┆ 6.1         │
│ very large region   ┆ synthetic control ┆ 143.3         ┆ 28.8        │
│ very large region   ┆ SC + intercept    ┆ 2.5           ┆ 3.0         │
│ very large region   ┆ SDID              ┆ 2.2           ┆ 4.7         │

• DiD is biased in both worlds: parallel trends fail.
• Plain synthetic control is the best estimator when R07 is a typical region, and useless when it's the biggest.
• SC with an intercept and SDID are unaffected by size, with a small bias (+$2–2.5k on a $14k effect) and similar spread.

No estimator wins everywhere. The rule that does: check the pre-period fit, then use the estimator whose assumptions the fit supports. If plain synthetic control fits well, it's hard to beat; if it can't fit, use an intercept or SDID, and report both.

46.6 Save your work with Git (20 min)

git add .
git commit -m "Day 46: big-region failure, SC with intercept, synthetic DiD, Monte Carlo"
git push

Day 46 checklist

• I can explain why synthetic control fails when the treated unit is outside the donors' range
• I can fix it with an intercept, and explain what that changes
• I can explain SDID's unit weights, time weights and weighted DiD
• I choose among DiD, SC and SDID using the pre-period fit, and report more than one