Mímisbrunnr知恵の泉

← 因果推論 一覧

🎓 レベル:標準 | 重要度:B(標準)

📎 前提:バックドア基準と識別 | 差分の差分と並行トレンド | 数理:単回帰分析(統計)

要点(BLUF)


概念:閾値が作る「ほぼ実験」

奨学金は GPA 3.0 以上、補助金は所得が基準額以下――現実には連続な指標 XX をある閾値で区切って処置を割り当てる制度が多い。閾値のすぐ上すぐ下の人は、能力も境遇もほとんど変わらないのに、片方だけ処置を受ける。この閾値近傍の偶然を実験の代わりに使うのがRDDだ。

シャープRDでは処置が閾値で確定する。

Ti=1{Xic}T_i = \mathbf{1}\{X_i \ge c\}
flowchart LR
    X["割当変数 X"] --> T["処置 T(X≥c のとき1)"]
    X --> Y["結果 Y(滑らかに依存)"]
    T --> Y

XXYY に滑らかに効く(交絡経路)。だが処置 TTcc不連続に切り替わる。XX 由来の滑らかな部分は cc で連続だから、cc での結果の段差TT の効果だけを映す。

識別の仮定:閾値での連続性

処置が無いときの条件付き期待値が閾値で連続であることを仮定する。

limxcE[Y(0)X=x]=limxcE[Y(0)X=x]\lim_{x\uparrow c} E[Y(0)\mid X=x] = \lim_{x\downarrow c} E[Y(0)\mid X=x]

Y(1)Y(1) も同様に連続)。直観は「閾値ちょうどで人が入れ替わるような操作(駆け込み・改ざん)が無い」こと。これが成り立てば、観測される結果の閾値ジャンプが処置効果に一致する。

τRD=limxcE[YX=x]limxcE[YX=x]=E[Y(1)Y(0)X=c]\tau_{\text{RD}} = \lim_{x\downarrow c} E[Y\mid X=x] - \lim_{x\uparrow c} E[Y\mid X=x] = E[Y(1)-Y(0)\mid X=c]

これは X=cX=c ちょうどの人たちに対する局所的な効果(LATE) であり、閾値から遠い層へは一般化できない。推定は、閾値近傍 Xc<h|X-c|<h で左右別に直線をあて、両者の cc での切片の差を取る。

τ^=α^α^,Y=α+β(Xc)+誤差  (左右別、Xc<h)\hat\tau = \hat\alpha_{\text{右}} - \hat\alpha_{\text{左}}, \qquad Y = \alpha + \beta\,(X-c) + \text{誤差}\ \ (\text{左右別、}|X-c|<h)

コード:真のジャンプを仕込んで局所線形で回収する

閾値 c=0c=0、連続な滑らかトレンド(処置と無関係)に真のジャンプ τ=4\tau=4 を仕込む。閾値の上下で全データを単純比較すると滑らかトレンドのぶん偏る。局所線形なら回収できる。

import numpy as np
import statsmodels.api as sm

# === 閾値で処置が決まる擬似データに真のジャンプを仕込み、局所線形で回収する ===
rng = np.random.default_rng(11)
n = 8000
tau_true = 4.0                         # 仕込んだ真の効果(閾値でのジャンプ)
c = 0.0                                # 閾値(割当変数は中心化済み)

X = rng.uniform(-1, 1, n)              # 割当変数
f = 3 + 1.5*X + 2.0*X**2 + 1.5*X**3    # 連続な滑らかトレンド(処置と無関係)
T = (X >= c).astype(int)              # 閾値以上で処置(シャープRD)
Y = f + tau_true*T + rng.normal(0, 0.3, n)

# 素朴:閾値の上下で全データの平均を比較(トレンドを無視=偏る)
naive = Y[X >= c].mean() - Y[X < c].mean()

# 局所線形回帰:閾値近傍 |X-c|<h で左右別に直線をあて、切片の差=ジャンプ
h = 0.2
right = (X >= c) & (X < c + h)
left = (X < c) & (X > c - h)
intercept_right = sm.OLS(Y[right], sm.add_constant(X[right] - c)).fit().params[0]
intercept_left = sm.OLS(Y[left], sm.add_constant(X[left] - c)).fit().params[0]
jump = intercept_right - intercept_left

print(f"真の効果 tau_true        = {tau_true:.3f}")
print(f"素朴な上下比較           = {naive:.3f}")
print(f"局所線形(h={h})のジャンプ  = {jump:.3f}")

出力は次の通り。

真の効果 tau_true        = 4.000
素朴な上下比較           = 6.230
局所線形(h=0.2)のジャンプ  = 3.968

素朴な上下比較は6.230――閾値の上は XX が大きく、滑らかトレンド f(X)f(X) のぶん YY も高いので、効果に約 +2.2 が上乗せされる。局所線形は3.968で真値 4 をほぼ回収した。閾値の「すぐ近く」だけを見れば f(X)f(X) の差はほとんど無く、残るのはジャンプだけだからだ。

バンド幅依存性:狭いほど偏らない

局所線形の結果はバンド幅 hh に依存する。hh を広げると遠くの曲がった部分まで直線で近似し偏りが増える。狭めると偏りは減るが、使うデータが減り分散が増える。hh を動かしてジャンプ推定を描く。

import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt
import japanize_matplotlib

# === バンド幅 h を変えてジャンプ推定の安定性を見る(バイアス‐分散の綱引き)===
rng = np.random.default_rng(11)
n = 8000
tau_true = 4.0
c = 0.0
X = rng.uniform(-1, 1, n)
f = 3 + 1.5*X + 2.0*X**2 + 1.5*X**3
T = (X >= c).astype(int)
Y = f + tau_true*T + rng.normal(0, 0.3, n)

def rd_jump(h):
    right = (X >= c) & (X < c + h)
    left = (X < c) & (X > c - h)
    br = sm.OLS(Y[right], sm.add_constant(X[right] - c)).fit().params[0]
    bl = sm.OLS(Y[left], sm.add_constant(X[left] - c)).fit().params[0]
    return br - bl

bandwidths = np.arange(0.10, 1.01, 0.05)
jumps = [rd_jump(h) for h in bandwidths]

plt.figure(figsize=(7, 4))
plt.axhline(tau_true, color="green", ls="--", label="真の効果 = 4")
plt.plot(bandwidths, jumps, "o-", color="steelblue", label="局所線形のジャンプ推定")
plt.xlabel("バンド幅 h"); plt.ylabel("ジャンプ推定値")
plt.title("バンド幅依存性:狭いほど偏りが小さい")
plt.legend(); plt.tight_layout()
plt.show()

for h in [0.10, 0.30, 0.50, 1.00]:
    print(f"h={h:.2f}  ジャンプ={rd_jump(h):.3f}")

出力は次の通り。

h=0.10  ジャンプ=3.949
h=0.30  ジャンプ=3.967
h=0.50  ジャンプ=3.923
h=1.00  ジャンプ=3.384

狭い h=0.100.30h=0.10\sim0.30 では 3.95〜3.97 と真値 4 にほぼ一致h=1.00h=1.00 まで広げると 3.384 と下に偏る(左右で曲率が違う部分まで取り込むため)。実務では最適バンド幅の自動選択(Imbens–Kalyanaraman、要最新確認)や複数 hh での頑健性チェックを行う。

コード:閾値の段差を可視化する

ビン平均の散布図に左右の局所線形を重ね、閾値での段差=効果を目で確認する。

import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt
import japanize_matplotlib

# === RDDの可視化:ビン平均の散布図に左右の局所線形をあて、閾値の段差を示す ===
rng = np.random.default_rng(11)
n = 8000
tau_true = 4.0
c = 0.0
X = rng.uniform(-1, 1, n)
f = 3 + 1.5*X + 2.0*X**2 + 1.5*X**3
T = (X >= c).astype(int)
Y = f + tau_true*T + rng.normal(0, 0.3, n)

# ビン平均(散布図の見やすさのため)
bins = np.linspace(-1, 1, 41)
idx = np.digitize(X, bins)
bx = np.array([X[idx == k].mean() for k in range(1, len(bins)) if (idx == k).any()])
by = np.array([Y[idx == k].mean() for k in range(1, len(bins)) if (idx == k).any()])

# 局所線形(h=0.3)の左右フィット
h = 0.3
right = (X >= c) & (X < c + h)
left = (X < c) & (X > c - h)
fr = sm.OLS(Y[right], sm.add_constant(X[right] - c)).fit()
fl = sm.OLS(Y[left], sm.add_constant(X[left] - c)).fit()
xr = np.linspace(0, h, 50); xl = np.linspace(-h, 0, 50)
yr = fr.params[0] + fr.params[1]*xr
yl = fl.params[0] + fl.params[1]*xl
jump = fr.params[0] - fl.params[0]

plt.figure(figsize=(7.5, 4.5))
plt.scatter(bx, by, s=18, color="lightgray", label="ビン平均")
plt.plot(xl, yl, color="steelblue", lw=2, label="左の局所線形")
plt.plot(xr, yr, color="crimson", lw=2, label="右の局所線形")
plt.axvline(c, color="gray", ls="--", lw=1)
plt.annotate(f"ジャンプ ≈ {jump:.2f}", xy=(0, (fr.params[0]+fl.params[0])/2),
             xytext=(0.25, fl.params[0]-0.5), fontsize=11,
             arrowprops=dict(arrowstyle="->"))
plt.xlabel("割当変数 X(閾値=0)"); plt.ylabel("結果 Y")
plt.title("シャープRD:閾値での不連続=因果効果")
plt.legend(loc="upper left"); plt.tight_layout()
plt.show()
print(f"閾値でのジャンプ推定 = {jump:.3f}(真の効果 = {tau_true:.3f})")

出力は 閾値でのジャンプ推定 = 3.967(真の効果 = 4.000)。図では閾値 X=0X=0 で左の直線(約3.0)から右の直線(約7.0)へ約4の段差が現れ、滑らかなビン平均がその不連続を裏づける。

ファジーRD:処置が確率的にしか変わらないとき

閾値で処置確率が 0→1 ではなく 0.3→0.7 のように一部だけ ジャンプする場合はファジーRD。結果のジャンプを処置のジャンプで割る、IV と同じ比の形になる。

τfuzzy=limxcE[YX=x]limxcE[YX=x]limxcE[TX=x]limxcE[TX=x]\tau_{\text{fuzzy}} = \frac{\lim_{x\downarrow c}E[Y\mid X=x]-\lim_{x\uparrow c}E[Y\mid X=x]}{\lim_{x\downarrow c}E[T\mid X=x]-\lim_{x\uparrow c}E[T\mid X=x]}

閾値を操作変数とみなした 操作変数法と2SLS の局所版で、識別されるのは閾値近傍のコンプライアのLATEだ。

仮定の直観:なぜ「閾値の近く」だけ信じるのか

RDDの強みは、識別仮定が連続性ひとつで、しかも閾値前後の密度や共変量の連続性として部分的に検証できることだ(操作の痕跡=閾値直下/直上で人数や属性が跳ねていないか)。代償は外的妥当性:得られるのは X=cX=c ちょうどの効果だけで、閾値から離れた人には何も言えない。狭いバンド幅が偏りを減らすのは、cc の近くでだけ「ほぼランダム」が成り立つからだ。

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

関連ノート