0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

因果推論 Day 19/全30回 合成コントロール法、存在しない対照群を合成する

0
Posted at

この連載について

因果推論を「本を読んだ」で終わらせず、自分の言葉で説明でき、コードで再現できる状態まで落とす30日連載です。直前のDay 18では、処置群と対照群の前後差を2回引くDID(差分の差分法)と、その土台になる平行トレンド仮定、プレトレンドを検査するイベントスタディまでを実データで確かめました。今日はその続き、平行トレンドを任せられる対照群が1つも見つからないときの道具です。データは計量経済学の教科書 Mixtape の texas、テキサス州の刑務所拡張と黒人男性収監者数です。


TL;DR(3行)

  • 合成コントロール=ドナープール50州の凸結合で「処置が無かったテキサス」を作る道具。 重みはフロリダ 37%・ニューヨーク 36%・イリノイ 27% の3州だけに立ち、事前8年の適合誤差はRMSPE 863人(事前平均の4.4%)。
  • 実測と合成の乖離がそのまま効果の推定値。 1993年の刑務所拡張後に実測だけが離陸し、2000年の乖離は +28,901人(実測61,861 vs 合成32,960)、処置後8年の平均で +21,013人
  • 標準誤差は出ません。 代わりに処置ラベルを全州へ付け替えるプラセボ置換検定(Day 5 のFisherの再会)で語ります。テキサスの乖離は 51州中1位、p = 1/51 ≈ 0.020

今日の問い

Day 18 のDIDは、処置群と対照群の「前後差の差」だった。でも処置されたのはテキサス州1つだけで、規模も伸びも似た州は1つも無い。
比較できる対照群が存在しないとき、因果効果はどう測るのか。そして推定値が1つしか出ない設計で、「これは偶然ではない」とどう言うのか?

答えは2段構えです。前半の問いには「対照群が無いなら、複数の州を混ぜて作ればいい」が合成コントロール法の答えです。後半の問いには、標準誤差の代わりに「処置ラベルを他の州に付け替えて同じ分析を回し、本物の順位で珍しさを測る」というFisher流の答えを使います。今日はこの両方を scipy で手組みします。


概念① 合成コントロール:ドナープールの凸結合で対照群を作る

処置を受けていない州の集まりをドナープールと呼びます。合成コントロール法は、ドナープールの加重平均で処置ユニットの「処置が無かった世界」の軌跡を作ります。

$$
\hat{Y}{1t}^{\text{synth}} = \sum{j=2}^{J+1} w_j, Y_{jt}
\qquad\text{s.t.}\quad w_j \ge 0,\quad \sum_j w_j = 1
$$

重み$w$は、処置前の期間で実測の軌跡に最も近づくよう選びます。

$$
\min_{w}\ \sum_{t < T_0} \Bigl( Y_{1t} - \sum_j w_j, Y_{jt} \Bigr)^{2}
$$

ポイントは、重みに「0以上・合計1」という凸結合の制約を課すことです。ただの回帰(制約なしの最小二乗)でも事前期間には適合できますが、その場合「カリフォルニアを$-0.5$州分」のような解が平気で出てきます。凸結合ならそれが起きません。第一に、ドナーの値の範囲の外へ外挿しないので、事前期間だけ辻褄を合わせた不自然な合成ができにくくなります。第二に、解が単体(重みの合計1の三角形の高次元版)の角に落ちるため重みが疎になり、「フロリダ37%+ニューヨーク36%+イリノイ27%」のように合成の中身をそのまま言葉で説明できます。

もうひとつの柱が事前適合(pre-treatment fit)の規律です。処置前の軌跡をなぞれない合成は、処置後の反事実としても信用できません。合成コントロールでは、まず事前適合の誤差(RMSPE)を確認し、なぞれていなければその対象にこの手法を使わない、が運用の第一歩です。

なお Abadie, Diamond & Hainmueller (2010) の原式は、アウトカムのラグに加えて共変量も距離に入れ、変数の重み付け行列$V$まで最適化します。今日は Facure の実装流儀に合わせ、事前期間のアウトカムの軌跡だけに適合させる簡略版でいきます。仕組みの芯は同じです。

概念② プラセボ置換検定:標準誤差の代わりに順位で語る

合成コントロールの推定値は「テキサスの乖離」というたった1つの数です。処置ユニットが1つしかないので、標本を増やして標準誤差を縮める、という漸近論の道具立てが最初から立ちません。

ここで Day 5 のFisherが戻ってきます。鋭い帰無仮説「処置の効果はどの州にも一切無い」が正しいなら、テキサスが処置されたことに意味はなく、どの州に処置ラベルを付け替えても同じ分析ができるはずです。そこで、ドナープールの各州を順に「処置された州」とみなして同じ手続きを回します。これがプラセボ置換検定です。効果が無いはずの州でも、事前適合の誤差やノイズで処置後の乖離はゼロにはなりません。その「偶然でも出る乖離」の束の中で、本物のテキサスの乖離が何位に来るかを数えます。テキサスが51州中1位なら、帰無仮説のもとでたまたま1位を引く確率は $1/51 \approx 0.020$。これがこの設計のp値です。

Day 5 で置いた規律もそのまま効きます。順位を測る統計量(処置後の平均乖離か、RMSPE比か)は事前に1つ決めてコミットします。いくつも試して一番小さいp値を選ぶのは反則です。これは後半で実際に効いてくるので、覚えておいてください。


手を動かす①:ドナープール50州から凸結合の重みを最適化

データは texas です。1980年、テキサス州は受刑者が起こした訴訟(Ruiz v. Estelle)に敗れ、過剰収容の解消を連邦裁判所に監督される立場になりました。当初は仮釈放を増やして収容数を調整していましたが、1993年、Richards 知事のもとで約10億ドルの刑務所建設計画が動き出し、収容能力は3年で概ね2倍になります。処置は「1993年の刑務所拡張」、アウトカムは州ごとの黒人男性の収監者数bmprison )です。データは51ユニット(50州+ワシントンDC)×1985〜2000年の16年分で、処置前の期間は8年あります。

重みの最適化は、概念①の式を scipy.optimize にそのまま渡すだけです。目的関数は事前期間の二乗誤差、制約は「0以上1以下」と「合計1」の2つです。

import numpy as np
import pandas as pd
from causaldata import texas
from scipy.optimize import minimize

d = texas.load_pandas().data
piv = d.pivot(index="year", columns="state", values="bmprison") / 1000  # 千人単位
T0 = 1993                                    # 刑務所拡張の開始年
years = piv.index.to_numpy()
pre = piv.loc[years < T0]                    # 事前期間 1985〜1992 の8年
donors = [s for s in piv.columns if s != "Texas"]

def synth_weights(y_pre, Y0_pre):
    """事前期間の軌跡への二乗誤差を最小にする凸結合の重みを解く"""
    k = Y0_pre.shape[1]
    def loss(w):
        r = y_pre - Y0_pre @ w
        return r @ r
    def grad(w):
        return -2 * Y0_pre.T @ (y_pre - Y0_pre @ w)
    res = minimize(loss, np.full(k, 1 / k), jac=grad,
                   bounds=[(0, 1)] * k,                            # 各重みは0以上1以下
                   constraints={"type": "eq",
                                "fun": lambda w: w.sum() - 1},     # 合計は1
                   method="SLSQP", options={"maxiter": 2000, "ftol": 1e-12})
    return res.x

w = synth_weights(pre["Texas"].to_numpy(), pre[donors].to_numpy())
ws = pd.Series(w, index=donors).sort_values(ascending=False)
print(ws[ws > 0.001].round(3))
print("重みの合計 = %.3f, 0.1%%を超える州 = %d" % (w.sum(), (w > 0.001).sum()))
Florida     0.373
New York    0.355
Illinois    0.272
dtype: float64
重みの合計 = 1.000, 0.1%を超える州 = 3

50州に自由に配ってよい重みが、フロリダ 37.3%・ニューヨーク 35.5%・イリノイ 27.2% の3州だけに立ちました。残り47州はすべて0.1%未満、実質ゼロです。これが概念①で述べた凸結合の疎性です。下の図の右を見ると、選ばれた理屈も読めます。ニューヨークはテキサスより収監者数の水準が高く、イリノイは低く、フロリダは伸び方が近い。水準も伸びも違う3州を混ぜることで、テキサスの事前8年の軌跡を挟み込むように再現しています。

Fig 3: 合成テキサスの中身(左:重みはフロリダ37.3%・ニューヨーク35.5%・イリノイ27.2%の3州だけに立つ、右:事前期間でこの凸結合がテキサスの軌跡をなぞる)


手を動かす②:合成テキサスと実テキサスの乖離=効果

重みが決まれば、合成テキサスの軌跡は行列積1回で出ます。処置後の期間について、実測から合成を引いた乖離が効果の推定値です。

synth = piv[donors].to_numpy() @ w           # 合成テキサスの軌跡
actual = piv["Texas"].to_numpy()
gap = actual - synth                         # 実測との乖離=効果の推定値

pre_rmspe = np.sqrt(np.mean(gap[years < T0] ** 2))
print("事前適合 RMSPE = %.0f 人(事前平均 %.0f 人の %.1f%%)"
      % (pre_rmspe * 1000, actual[years < T0].mean() * 1000,
         100 * pre_rmspe / actual[years < T0].mean()))
for yy in [1993, 1994, 1995, 2000]:
    i = int(np.where(years == yy)[0][0])
    print("%d年: 実測 %5d 人, 合成 %5d 人, 乖離 %+6d 人"
          % (yy, round(actual[i] * 1000), round(synth[i] * 1000), round(gap[i] * 1000)))
print("処置後8年の平均乖離 = %+d 人" % round(gap[years >= T0].mean() * 1000))
事前適合 RMSPE = 863 人(事前平均 19449 人の 4.4%)
1993年: 実測 29260 人, 合成 28584 人, 乖離   +676 人
1994年: 実測 40451 人, 合成 30037 人, 乖離 +10414 人
1995年: 実測 55602 人, 合成 31505 人, 乖離 +24097 人
2000年: 実測 61861 人, 合成 32960 人, 乖離 +28901 人
処置後8年の平均乖離 = +21013 人

まず事前適合を確認します。RMSPEは 863人、事前期間の平均約19,400人に対して4.4%です。8年の軌跡をこの精度でなぞれているので、合成を反事実として使う土俵に乗りました。そのうえで処置後を見ると、乖離は着工年の1993年にはまだ +676人 とほぼゼロで、1994年に +10,414人、収容能力の倍増が完了する時期にあたる1995年に +24,097人 へ跳ね、2000年には +28,901人 に達します。実測61,861人は合成32,960人の約1.9倍です。「拡張が完成するにつれて乖離が積み上がる」という制度の時間割と、乖離の立ち上がり方が噛み合っているのも心強い点です。

Fig 1: 実測テキサス vs 合成テキサス(事前8年をなぞった凸結合から、1993年の拡張後に実測だけが離陸し、2000年の乖離は+28,901人)

Day 18 との対比で言うと、DIDは「選んだ対照群と平行だったはず」という仮定を事後には検証しづらいのに対し、合成コントロールは「対照をなぞれているか」が事前適合という数字で最初から見えます。対照群を1つ選ぶ代わりに、なぞれる対照を最適化で作ってしまう。これがこの手法の発明です。


手を動かす③:全州でプラセボ置換し、効果の順位で有意性を語る

ではこの +21,013人は「偶然でも出る乖離」なのでしょうか。概念②のとおり、全州に処置ラベルを付け替えて同じ手続きを回します。1つだけ注意して、プラセボのドナープールからは本物の処置州テキサスを外します。処置後に爆発するテキサスがプラセボの合成に混ざると、効果が無いはずの州の乖離を歪めるからです。

rows = []
for target in piv.columns:                   # 全51州を順に「処置された州」とみなす
    pool = [s for s in piv.columns if s not in (target, "Texas")]
    wp = synth_weights(pre[target].to_numpy(), pre[pool].to_numpy())
    g = piv[target].to_numpy() - piv[pool].to_numpy() @ wp
    rows.append((target,
                 np.sqrt(np.mean(g[years < T0] ** 2)),    # 事前RMSPE
                 np.sqrt(np.mean(g[years >= T0] ** 2)),   # 事後RMSPE
                 g[years >= T0].mean()))                  # 処置後の平均乖離
r = pd.DataFrame(rows, columns=["state", "pre_rmspe", "post_rmspe", "effect"])
r = r.set_index("state").sort_values("effect", ascending=False)

rank = r.index.get_loc("Texas") + 1
for s, row in r.head(3).iterrows():
    print("%-12s 処置後の平均乖離 %+8.0f 人" % (s, row.effect * 1000))
print("テキサスの順位 = %d位 / 51州, p値 = %d/51 = %.3f" % (rank, rank, rank / 51))
print("事前適合が最悪の州 = %s(RMSPE %.0f 人)"
      % (r.pre_rmspe.idxmax(), r.pre_rmspe.max() * 1000))
ratio_rank = int((r.post_rmspe / r.pre_rmspe).rank(ascending=False)["Texas"])
print("事後/事前RMSPE比で並べ直すと、テキサスは %d位" % ratio_rank)
Texas        処置後の平均乖離   +21013 人
California   処置後の平均乖離    +7384 人
Georgia      処置後の平均乖離    +3764 人
テキサスの順位 = 1位 / 51州, p値 = 1/51 = 0.020
事前適合が最悪の州 = California(RMSPE 4489 人)
事後/事前RMSPE比で並べ直すと、テキサスは 9位

処置後の平均乖離で並べると、テキサスは 51州中1位 です。2位のカリフォルニア +7,384人 に3倍近い差をつけています。「効果はどの州にも無い」が正しければ、本物がたまたま1位を引く確率は $1/51 \approx$ 0.020。下の図が示すとおり、プラセボ50州の乖離は、事前適合が悪いカリフォルニア(破線)を除けば概ね ±6,000人 の帯に収まり、テキサスだけが帯のはるか外へ抜けていきます。

Fig 2: プラセボの束と本物の乖離(灰色:処置ラベルを付け替えた50州、赤:本物のテキサス。テキサスだけが束の外に出て51州中1位、p≈0.020)

出力の下2行は、この検定の運用で大事な2つの現実です。ひとつ目、プラセボ側で事前適合が最悪なのはカリフォルニア(RMSPE 4,489人、テキサスの5倍超)です。事前期間で全米最大の収監者数を持つ州は、凸結合では作れません。混ぜて作れるのはドナーの範囲の内側だけだからです。原論文の流儀にならい「事前適合がテキサスの2倍より悪い州はプラセボから除外」しても、落ちるのはカリフォルニアだけで、テキサスは 50州中1位、p ≈ 0.02 のまま動きません。

ふたつ目は統計量の選択です。原論文が推す統計量は、事前適合の悪さで割り引いた事後RMSPE/事前RMSPE比ですが、この比で並べ直すとテキサスは 9位 に落ちます。順位が落ちた時点で私はまず自分の実装を疑い、比の上位から確認しました。1位はペンシルベニアで、比は339倍。その事前RMSPEは 約6人 です。事前8時点に対してドナーが50州もあると、中規模の州は凸結合でほぼ完璧に補間できてしまい、ほぼゼロの分母が比を爆発させるのです。これは効果の証拠ではなく過剰適合の産物ですが、機械的に比で並べれば順位はこうなります。Day 5 の「統計量は事前に1つ決めてコミットする。並べてから選ばない」という規律は、まさにこの場面のためにあります。ちなみに共変量と入れ子最適化を使う Mixtape の本式では、テキサスは46州中2位(p = 0.04)。どの流儀でも「テキサスは束の外」という図の結論は変わりません。


つまづき・誤解しやすい点

  • 事前適合が悪ければ、その対象には使えません。 合成が処置前をなぞれていないなら、処置後の乖離は効果とノイズの区別がつきません。カリフォルニアのように、ドナーの凸包の外にあるユニット(全州で最大・最小の州など)は原理的に合成できません。事前RMSPEの報告を必須にし、事前期間の後半を検証用に取り置いて「合成が当てられるか」を確かめてから本番に進むのが安全です。
  • 事前適合が良すぎるのも危険です。 事前が8時点しかないのにドナーが50州あれば、効果と無関係にほぼ完璧な補間が偶然できてしまいます(ペンシルベニアの事前誤差6人がその実例)。手を動かす③のRMSPE比の爆発はこの過剰適合の裏返しです。原式が共変量で距離を縛り、実務でドナープールを処置ユニットと似た州へ絞るのは、この「なんでも作れてしまう」自由度を殺すためです。
  • p = 0.020 は「効果の証明」ではありません。 プラセボ検定が語れるのは「効果ゼロの世界でこの乖離の順位は珍しい」ことだけです。1993年前後にテキサスだけを襲った別の要因(州固有の政策変更や経済ショック)があれば、それも乖離に混ざります。ドナーへの処置の波及(テキサスの受刑者移送や犯罪の州間移動)があればSUTVA違反で対照側も汚れます。干渉の扱いは Day 28 でまとめて掘ります。

GISデータ実務での使い方

  • 処置された自治体が1つしかない施策評価の受け皿になります。 ある政令市だけが衛星漏水調査を全域導入した、ある県だけが管路更新の補助制度を変えた。DIDの対照選びで悩む典型場面ですが、残りの自治体をドナープールにして「導入しなかったその市」を合成すれば、月次の漏水修繕件数や管路事故件数の反事実軌跡を作れます。
  • 事前適合が「適用可否の門番」になります。 導入前2〜3年の系列を合成がなぞれるかをまず確認し、なぞれなければその自治体にはこの手法を使わない、と最初に判定できます。仮定が数字で見えるのは、非交絡性を祈るしかない場面(Day 6)より運用がずっと楽です。
  • 順位ベースの推測は非統計家に通ります。 「効果が無ければ、50自治体のプラセボの中で1位の乖離になる確率は2%」という語り口は、標準誤差や信頼区間より説明の通りがよく、自治体向け報告書と相性が良い言い回しです。ただし隣接自治体は広域連携や同じ工事業者経由で処置が波及しがちなので、ドナーから外す判断が要ります。空間の対照選びは Day 28 の主題です。

参考(本棚)

  • Abadie, Diamond & Hainmueller (2010) "Synthetic Control Methods for Comparative Case Studies"。原論文。カリフォルニアのたばこ規制 Proposition 99 を題材に、共変量と重み行列$V$を含む本式とプラセボ検定の運用を定式化
  • Scott Cunningham『Causal Inference: The Mixtape』第10章。texas データの本家。テキサス刑務所拡張の制度背景と、共変量込みの本式で46州中2位(p = 0.04)という結果もこの章
  • Matheus Facure『Causal Inference in Python』第9章。事前アウトカムの軌跡に適合させる今日の実装流儀の下敷き
  • Imbens & Rubin『Causal Inference for Statistics, Social, and Biomedical Sciences』第5章。プラセボ置換検定の土台になるFisher流ランダム化推測(Day 5 の下敷き)
  • データ: texas( causaldata パッケージ経由、MITライセンス、取得日2026-08-15)。1985〜2000年の州別の黒人男性収監者数と州属性。Cunningham の教材データ

次回予告(Day 20)

今日は「対照群が無いなら合成する」という力技を見ました。次回は逆に、制度の側が対照群を用意してくれている状況です。所得スコアが基準以下なら給付、超えたら対象外。こうして閾値で機械的に処置が決まる制度は、閾値のすぐ両側では処置がほぼ偶然で決まる、つまり実験に化けます。回帰不連続デザイン(RDD)で、制度の「際」から因果効果を取り出します。

0
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?