合成コントロール法

🎓 レベル：発展　|　重要度：B（標準）

📎 前提：差分の差分と並行トレンド　|　潜在結果モデル

要点（BLUF）

合成コントロール法（SCM） は、処置ユニットが1つしかない（ある州・国に施策が入った等）状況で、複数の対照ユニット（ドナー）の重み付き平均で「もし介入が無ければ」を作る。
重みは 事前期間の軌跡を最もよく再現するよう、 $w_j\ge 0,\ \sum_j w_j=1$ （凸結合）で最適化する。事後の実測−合成のギャップが効果。
識別は 干渉なし＋処置ユニットが対照の凸包に入り長い事前期間でよく当てはまること。検証はプラセボ検定（対照を偽の処置にして当てはめ過剰でないか）。

概念：対照群が「ちょうどよく」存在しないとき

差分の差分と並行トレンドは処置群と対照群が複数ユニットあり、両者が平行に動くことを前提にした。だが現実には「カリフォルニア州だけがタバコ税を上げた」のように、処置を受けたのが1ユニットだけで、しかも単独の対照州ではうまく似ていないことが多い。

SCMの発想は、複数の対照ユニットを混ぜて処置ユニットそっくりの「合成対照」を作ること。タバコ消費なら、カリフォルニアに似た消費トレンドを、ユタ・ネバダ・コロラド…の加重平均で再現する。

flowchart LR
    F["共通の潜在ファクター F_t（景気・嗜好トレンド）"] --> D["ドナー群の結果"]
    F --> Tr["処置ユニットの結果"]
    I["介入（処置）"] --> Tr

ドナーと処置ユニットは共通のファクターに晒される。合成対照（ドナーの凸結合）が処置ユニットのファクター曝露を再現できれば、介入だけが入った処置ユニットとの事後の差が効果を映す。

識別の仮定

潜在結果がファクター構造に従うとする（ $\mu_j$ はユニット $j$ のローディング、 $F_t$ は共通ファクター）。

Y_{jt}(0) = \boldsymbol{\mu}_j^\top \boldsymbol{F}_t + \varepsilon_{jt}

処置ユニット（ $j=1$ ）には事後 $t> T_0$ で効果 $\tau_t$ が乗る： $Y_{1t} = Y_{1t}(0) + \tau_t$ 。重み $\boldsymbol w=(w_2,\dots,w_{J+1})$ を事前期間の当てはめで決める。

\boldsymbol{w}^\* = \arg\min_{\boldsymbol w}\ \sum_{t\le T_0}\Big(Y_{1t} - \sum_{j\ge 2} w_j Y_{jt}\Big)^2 \quad \text{s.t.}\quad w_j\ge 0,\ \sum_j w_j = 1

効果の推定量は事後のギャップ：

\hat\tau_t = Y_{1t} - \sum_{j\ge 2} w_j^\* Y_{jt}, \qquad t> T_0

これが $\tau_t$ に一致するための鍵は2つ。(1) 干渉なし：ドナーは処置ユニットの介入に影響されない（ユニット間SUTVA）。(2) 凸包条件＋良好な事前フィット：処置ユニットの事前軌跡が長い事前期間にわたりドナーの凸結合で再現できる（外挿でなく内挿）。凸結合に制約するのは、外挿による過剰適合を防ぎ重みを解釈可能にするためだ。

コード：真の効果を仕込んで重み最適化で回収する

1つの処置ユニットと9つのドナーの擬似パネルを作る。共通ファクター（時間トレンド）に対するローディングを各ユニットに与え、処置ユニットのローディングはドナー0,1,2の平均（凸包の内側）に置く。事後に真の効果 $\tau=5$ を加える。重みは scipy.optimize.minimize（SLSQP、単体上の制約付き最小化）で求める。

import numpy as np
from scipy.optimize import minimize

# === 1処置ユニット＋ドナープールの擬似パネル。事前を凸結合で再現し事後ギャップ＝効果 ===
rng = np.random.default_rng(20)
n_donors = 9
T0, T1 = 25, 10                    # 事前25期・事後10期
T = T0 + T1
tau_true = 5.0                     # 仕込んだ真の介入効果

t = np.arange(T)
F1 = np.sin(t / 3.0)               # 共通の潜在ファクター（時間トレンド）
F2 = t / 10.0
level = rng.uniform(8, 12, n_donors)
mu1 = rng.uniform(0.5, 2.0, n_donors)
mu2 = rng.uniform(0.5, 2.0, n_donors)
donors = np.array([level[j] + mu1[j]*F1 + mu2[j]*F2 + rng.normal(0, 0.15, T)
                   for j in range(n_donors)]).T          # 形状 (T, ドナー数)

# 処置ユニット：ドナー0,1,2 を平均したローディング（ドナーの凸包の内側にある）
lv = level[[0, 1, 2]].mean(); m1 = mu1[[0, 1, 2]].mean(); m2 = mu2[[0, 1, 2]].mean()
treated_Y0 = lv + m1*F1 + m2*F2 + rng.normal(0, 0.15, T)  # 介入が無ければこうなる
treated = treated_Y0 + np.where(t >= T0, tau_true, 0.0)    # 事後に真の効果を加える

# 重み最適化：w>=0, Σw=1 で事前期間のフィットを最小化（凸結合に制約）
def fit_weights(target_pre, pool_pre):
    n = pool_pre.shape[1]
    cons = {"type": "eq", "fun": lambda w: np.sum(w) - 1}
    res = minimize(lambda w: np.sum((target_pre - pool_pre @ w)**2),
                   np.ones(n)/n, method="SLSQP", bounds=[(0, 1)]*n, constraints=cons)
    return res.x

w = fit_weights(treated[:T0], donors[:T0])
synthetic = donors @ w
pre_rmse = np.sqrt(np.mean((treated[:T0] - synthetic[:T0])**2))
post_effect = (treated[T0:] - synthetic[T0:]).mean()

print("重み（ドナー0..8）:", np.round(w, 3))
print(f"事前フィットRMSE  = {pre_rmse:.3f}")
print(f"事後ギャップ=効果 = {post_effect:.3f}（真値 {tau_true:.3f}）")

出力は次の通り。

重み（ドナー0..8）: [0.169 0.382 0.111 0.    0.098 0.072 0.167 0.    0.   ]
事前フィットRMSE  = 0.200
事後ギャップ=効果 = 4.919（真値 5.000）

重みはドナー0・1・2・6 など似た少数のユニットに集中し（ローディングが近い相手を選んだ）、残りは0。事前フィットRMSEは0.200とノイズ水準（0.15）並みに小さく、合成対照が処置ユニットの事前軌跡をよく再現できている。その結果、事後ギャップは4.919で真値 5 をほぼ回収した。凸結合の制約がスパースで解釈可能な重みを生むのもSCMの利点だ。

プラセボ検定：偽陽性でないか

ギャップが本物の効果か、たまたまの当てはめ誤差かを区別するため、各ドナーを「偽の処置」にして同じ手続きを回す（in-space placebo）。介入が無いドナーのギャップは0近傍に留まるはずで、処置ユニットのギャップだけが突出すれば効果は本物だと言える。

import numpy as np
from scipy.optimize import minimize

# === プラセボ検定：各ドナーを「偽の処置」にして当てはめ、ギャップが0近傍か確認 ===
rng = np.random.default_rng(20)
n_donors = 9
T0, T1 = 25, 10
T = T0 + T1
tau_true = 5.0
t = np.arange(T)
F1 = np.sin(t / 3.0); F2 = t / 10.0
level = rng.uniform(8, 12, n_donors)
mu1 = rng.uniform(0.5, 2.0, n_donors); mu2 = rng.uniform(0.5, 2.0, n_donors)
donors = np.array([level[j] + mu1[j]*F1 + mu2[j]*F2 + rng.normal(0, 0.15, T)
                   for j in range(n_donors)]).T
lv = level[[0, 1, 2]].mean(); m1 = mu1[[0, 1, 2]].mean(); m2 = mu2[[0, 1, 2]].mean()
treated_Y0 = lv + m1*F1 + m2*F2 + rng.normal(0, 0.15, T)
treated = treated_Y0 + np.where(t >= T0, tau_true, 0.0)

def fit_weights(target_pre, pool_pre):
    n = pool_pre.shape[1]
    cons = {"type": "eq", "fun": lambda w: np.sum(w) - 1}
    res = minimize(lambda w: np.sum((target_pre - pool_pre @ w)**2),
                   np.ones(n)/n, method="SLSQP", bounds=[(0, 1)]*n, constraints=cons)
    return res.x

# 本物の処置ユニット
w_t = fit_weights(treated[:T0], donors[:T0])
treated_gap = (treated[T0:] - (donors @ w_t)[T0:]).mean()

# プラセボ：各ドナーを偽の処置に、残りをプールにして同じ手続き
placebo_gaps = []
for p in range(n_donors):
    pool = np.delete(donors, p, axis=1)
    w_p = fit_weights(donors[:T0, p], pool[:T0])
    gap = (donors[T0:, p] - (pool @ w_p)[T0:]).mean()
    placebo_gaps.append(gap)
placebo_gaps = np.array(placebo_gaps)

print(f"処置ユニットの事後ギャップ = {treated_gap:.3f}")
print(f"プラセボ事後ギャップ: 平均={placebo_gaps.mean():.3f} "
      f"絶対値の最大={np.abs(placebo_gaps).max():.3f}")

出力は次の通り。

処置ユニットの事後ギャップ = 4.919
プラセボ事後ギャップ: 平均=0.031 絶対値の最大=1.603

プラセボのギャップは平均0.031とほぼ0を中心にばらつき、絶対値の最大でも1.603。処置ユニットの4.919はすべてのプラセボを上回る。9個のプラセボ中ゼロ個しか処置ユニットを超えないので、並べ替え検定的な p 値は約 $1/10=0.1$ ――効果はノイズでは説明しづらい、と読める（ドナー数が増えるほど検出力が上がる）。

コード：軌跡とギャップを可視化する

合成対照が事前を追従し事後で乖離する様子と、ギャップのプラセボ比較を2枚で描く。

import numpy as np
from scipy.optimize import minimize
import matplotlib.pyplot as plt
import japanize_matplotlib

# === 図：処置ユニット vs 合成コントロールの軌跡と、ギャップのプラセボ比較 ===
rng = np.random.default_rng(20)
n_donors = 9
T0, T1 = 25, 10
T = T0 + T1
tau_true = 5.0
t = np.arange(T)
F1 = np.sin(t / 3.0); F2 = t / 10.0
level = rng.uniform(8, 12, n_donors)
mu1 = rng.uniform(0.5, 2.0, n_donors); mu2 = rng.uniform(0.5, 2.0, n_donors)
donors = np.array([level[j] + mu1[j]*F1 + mu2[j]*F2 + rng.normal(0, 0.15, T)
                   for j in range(n_donors)]).T
lv = level[[0, 1, 2]].mean(); m1 = mu1[[0, 1, 2]].mean(); m2 = mu2[[0, 1, 2]].mean()
treated_Y0 = lv + m1*F1 + m2*F2 + rng.normal(0, 0.15, T)
treated = treated_Y0 + np.where(t >= T0, tau_true, 0.0)

def fit_weights(target_pre, pool_pre):
    n = pool_pre.shape[1]
    cons = {"type": "eq", "fun": lambda w: np.sum(w) - 1}
    res = minimize(lambda w: np.sum((target_pre - pool_pre @ w)**2),
                   np.ones(n)/n, method="SLSQP", bounds=[(0, 1)]*n, constraints=cons)
    return res.x

w_t = fit_weights(treated[:T0], donors[:T0])
synthetic = donors @ w_t

fig, axes = plt.subplots(1, 2, figsize=(11, 4.2))
# 左：軌跡
axes[0].plot(t, treated, "o-", color="crimson", label="処置ユニット（実測）")
axes[0].plot(t, synthetic, "s--", color="steelblue", label="合成コントロール")
axes[0].axvline(T0 - 0.5, color="gray", ls=":", label="介入時点")
axes[0].set_title("軌跡：事後にだけ乖離"); axes[0].set_xlabel("期"); axes[0].set_ylabel("結果 Y")
axes[0].legend(fontsize=9)
# 右：ギャップ（処置＋プラセボ）
for p in range(n_donors):
    pool = np.delete(donors, p, axis=1)
    w_p = fit_weights(donors[:T0, p], pool[:T0])
    axes[1].plot(t, donors[:, p] - pool @ w_p, color="lightgray", lw=1)
axes[1].plot(t, treated - synthetic, color="crimson", lw=2.2, label="処置ユニット")
axes[1].axvline(T0 - 0.5, color="gray", ls=":")
axes[1].axhline(0, color="black", lw=0.6)
axes[1].axhline(tau_true, color="green", ls="--", lw=1, label="真の効果=5")
axes[1].set_title("ギャップ：処置だけ事後に跳ねる"); axes[1].set_xlabel("期")
axes[1].set_ylabel("実測 − 合成"); axes[1].legend(fontsize=9)
plt.tight_layout()
plt.show()
print("事後ギャップ平均 =", round((treated[T0:] - synthetic[T0:]).mean(), 3))

左図では処置ユニット（赤）と合成対照（青）が事前は重なり、介入時点で赤だけ上振れする。右図のギャップでは、処置ユニット（赤太線）だけが事後に約5へ跳ね、プラセボ（灰色）は0周辺に留まる。合成対照が反実仮想をよく近似できていることの視覚的証拠だ。

仮定の直観：なぜ凸結合・長い事前期間なのか

合成コントロールが反実仮想を当てるのは、多数の事前期間にわたって処置ユニットを再現できた重みは、見えないファクター・ローディングまで揃えている可能性が高いからだ（短期間の偶然の一致なら過剰適合）。凸結合（重み非負・和1）に縛るのは、ドナーの範囲外への外挿を禁じ、解釈可能で頑健な対照を作るため。逆に処置ユニットがドナーの凸包の外（どのドナーより極端）だと、良いフィットは得られず推定は信用できない。

⚠️ よくある誤解・落とし穴

事前フィットが悪いのに使う：事前で合成が処置ユニットを再現できていないなら、事後ギャップは効果でなく当てはめ誤差。プラセボでも事前RMSEが大きいユニットは比較から外す。
ドナー汚染（干渉）：介入の波及を受けたユニットをドナーに入れると、合成対照自体が汚れる。
過剰適合：事前期間が短く、ノイズに合わせただけの重みは事後で外れる。長い事前期間が要る。
単一の効果しか出ない：処置ユニットへの効果であって母集団平均ではない。
推論はプラセボ／並べ替え：通常の標準誤差は使えない。in-space・in-time プラセボで分布を作る。