この連載について
因果推論を「本を読んだ」で終わらせず、自分の言葉で説明でき、コードで再現できる状態まで落とす30日連載です。直前のDay 22では、平均効果の裏にある効果の異質性を扱いました。CATE(条件付き平均処置効果)という推定対象、サブグループを事後的に漁る罠、そして Phase 4 の地図(「誰に効くか」の CATE 系と「平均を頑健に」の DML 系)です。今日はその CATE を、手持ちの機械学習でそのまま組みにいきます。データは合成です。真の $\tau(x)$ を私たちが握ることで、推定の答え合わせができるからです。この「正解を知っている側が採点する」やり方を、本稿ではオラクル採点と呼びます(オラクルは「正解を知っている神託役」の意味です)。
TL;DR(3行)
- CATE の教師ラベル $\tau_i$ は誰にも観測できない(因果推論の根本問題の再来)。だから普通の教師あり回帰は組めず、予測モデルの組み合わせで $\hat\tau(x)$ を作るメタラーナーを使います。S/T/X の3つは sklearn だけで手組みできます。
- 弱点も手組みで見えます。S-learner は正則化が $\tau$ を潰す(正則化はモデルが複雑になりすぎないよう締め付ける仕組みで、締めすぎると小さな信号が消えます。木20本で全評価点の $\hat\tau=0$、「効果は存在しない」と言い出す)。T-learner は群サイズの不均衡に弱い(処置15%の今日の設定で、$\hat\tau$ の散らばり 0.98、真の散らばりは 0.71)。
- 採点は PEHE(効果版の RMSE。推定した効果と真の効果のずれを2乗平均して平方根を取った値)。今日の設定では X-learner が 0.239 で最良(S 0.329、T 0.677)です。ただし PEHE は真の $\tau$ を知る合成データでしか計算できず、実データでの評価は別の道具(Day 26)が要ります。
今日の問い
「誰に効くか」を予測モデルとして学びたい。でも1人につき観測できるのは $Y(1)$ か $Y(0)$ の片方だけだから、正解ラベル $\tau_i = Y_i(1) - Y_i(0)$ はデータのどの行にも存在しない。
ラベルの無い量を、手持ちの教師あり学習でどう推定するのか。そしてその推定を、正解ラベル無しでどう採点するのか?
Day 2 で見た因果推論の根本問題が、機械学習の言葉で戻ってきました。教師あり学習は (特徴量, ラベル) のペアを要求しますが、$\tau_i$ のラベル列は原理的に作れません。今日のゴールは、この壁を「予測問題への分解」で迂回するメタラーナーの設計図を、S/T/X の順に1つずつコードにすることです。
概念① メタラーナー:ラベルの無い回帰を、予測問題の組み合わせに分解する
出発点は、Day 6 の非交絡性です。割付を決める共変量 $X$ をすべて観測できているなら、CATE は条件付き期待値の差に書き直せます。
$$
\tau(x) = E[Y \mid W=1, X=x] - E[Y \mid W=0, X=x]
$$
読み下すと「特徴量が $x$ の人たちに絞ったとき、処置を受けた人の $Y$ の平均から、受けなかった人の $Y$ の平均を引いたもの」です。$E[\cdot \mid \cdot]$ は「縦棒の右の条件を満たす人だけで取った平均」と読んでください。左辺はラベルの無い量ですが、右辺に並んだ2つの条件付き期待値は、どちらもただの予測問題です。「処置群での $Y$ の予測」と「対照群での $Y$ の予測」なら、ラベルは観測済みの $Y$ そのもので、勾配ブースティングでもランダムフォレストでも(どちらも決定木をたくさん組み合わせる機械学習の手法です)、好きな学習器で当てにいけます。この分解の仕方には複数の流儀があり、その設計図(レシピ)をメタラーナーと呼びます。中に入れる学習器(base learner)は交換自由で、今日は全部 sklearn の GradientBoostingRegressor に固定します。
以下の記号は、$\hat\mu$(ミューハット)が「$Y$ を予測するモデル」、$\hat\tau$ が「効果の推定値」です。ハット記号は「データから推定した値」の印で、ハットの無い $\mu$ や $\tau$ が真の値を指します。
- S-learner(Single):$W$ をただの特徴量として $X$ に足し、モデル1本 $\hat\mu(x, w)$ を学習。$\hat\tau(x) = \hat\mu(x,1) - \hat\mu(x,0)$。
- T-learner(Two):対照群だけで $\hat\mu_0(x)$、処置群だけで $\hat\mu_1(x)$ の2本を学習。$\hat\tau(x) = \hat\mu_1(x) - \hat\mu_0(x)$。
- X-learner(Cross):T の2本を第1段に使い、第2段で $\tau$ の代理ラベルを作って $\tau$ を直接学習。仕上げに傾向スコアで2つの推定を混ぜます(手を動かす②で1ステップずつ)。
概念② PEHE:効果の予測誤差、ただし現実には計算できない採点表
CATE 推定の標準的な採点は PEHE(precision in estimating heterogeneous effects)です。見た目は RMSE(予測誤差を2乗して平均し、平方根を取った値。機械学習で最もよく使う精度指標の1つ)そのもので、対象が $Y$ ではなく $\tau$ に替わっています(原義は2乗平均で、平方根を取って報告する流儀が多く、本稿でも平方根の値を PEHE と呼びます)。
$$
\mathrm{PEHE} = \sqrt{E\bigl[(\hat\tau(X) - \tau(X))^2\bigr]}
$$
読み下すと「推定した効果 $\hat\tau$ と真の効果 $\tau$ の差を2乗し、全員分の平均を取って平方根に戻した値」です。0 に近いほど効果の予測が真に近く、単位は $Y$ と同じです。
ここで立ち止まる価値があります。なぜ通常の精度指標ではだめなのでしょうか。$Y$ に対する RMSE や決定係数(当てはまりの良さを 0〜1 で表す指標)が測るのは「アウトカムの当てはまり」であって、効果の当てはまりではありません。後で実際に見せますが、$Y$ の予測としてはそこそこのモデルが、効果については「全員ゼロ」という大外れの回答をすることがあります。交差検証(データを分割し、学習に使わなかった部分で精度を測る手順)も同じ理由で頼れません。検証データにあるのは $Y$ のラベルだけなので、交差検証で選ばれるのは「$Y$ をよく当てるモデル」であり、$\tau$ をよく当てるモデルではないのです。
そして本丸の問題があります。$\tau_i$ の正解列が無い以上、実データでは PEHE そのものが計算できません。採点表は定義できるのに、答案と照合する解答用紙が現実には手に入らない。だから今日は合成データで真の $\tau(x)$ を握り、オラクル採点で3つのラーナーの癖を掴みます。実データでの代理評価(順位づけの良さで測る Qini 曲線)は Day 26 に置きます。
手を動かす①:S/T-learner を手組みして、弱点まで見る
データの設計を先に言葉にします。ベースライン $\mu_0(x)$ は $x_2$ から $x_5$ が絡む複雑な関数、効果 $\tau(x)$ は $x_1$ だけで決まる単純なシグモイド(S字カーブ。0 から 2 へなだらかに上がり、ATE は 1.0)、処置確率は $x_2$ に依存させて(交絡)平均15%に絞ります。「背景は複雑、効果は単純、処置は少数派」という、Künzel らが X-learner の得意場面として挙げた構図です。実務でも施策を受けた側が少数派になることは多く、わざとらしい設定ではありません。
import numpy as np
from sklearn.ensemble import GradientBoostingRegressor
rng = np.random.default_rng(23)
N = 6000
# 共変量5個。ベースラインは x2〜x5 の複雑な関数、効果の異質性は x1 だけが決める
X = rng.uniform(-1, 1, size=(N, 5))
x1, x2, x3, x4, x5 = X.T
mu0 = 5 * np.sin(np.pi * x2) + 3 * x3 * x4 + 2 * x5 + 2 * x2 # ベースライン(複雑)
tau = 2 / (1 + np.exp(-4 * x1)) # 真の効果 τ(x):x1 のシグモイド(単純)
e = 1 / (1 + np.exp(2.0 - 1.5 * x2)) # 傾向スコア:x2 が高いほど処置(交絡)
w = rng.binomial(1, e)
y = mu0 + w * tau + rng.normal(0, 1, N)
print("処置群 %d 人 / 対照群 %d 人" % (w.sum(), (1 - w).sum()))
print("真のATE = %.2f(τ(x) は %.2f から %.2f)" % (tau.mean(), tau.min(), tau.max()))
print("素朴な差 = %.2f" % (y[w == 1].mean() - y[w == 0].mean()))
処置群 939 人 / 対照群 5061 人
真のATE = 1.00(τ(x) は 0.04 から 1.96)
素朴な差 = 4.02
出力の読み方は、1行目が群のサイズ、2行目が真の効果の平均と範囲、3行目が何も調整しない単純な平均差です。素朴な引き算は 4.02 で、真の ATE 1.00 の4倍です。$x_2$ が処置確率とベースラインの両方を押し上げているからで、Day 6 で見た交絡そのものです。メタラーナーは共変量 $X$ で条件付けた予測の差を取るので、割付を決める $x_2$ が $X$ に入っている限り、この交絡は調整されます。では S と T を組みます。
def pehe(tau_hat, tau_true):
return np.sqrt(np.mean((tau_hat - tau_true) ** 2))
# 採点用の新データ。合成だから真の τ(x) を知っていて、答え合わせできる
X_test = rng.uniform(-1, 1, size=(2000, 5))
tau_test = 2 / (1 + np.exp(-4 * X_test[:, 0]))
gbm = dict(n_estimators=200, max_depth=3, learning_rate=0.05,
min_samples_leaf=30, random_state=0)
# S-learner:W をただの特徴量として1本のモデルに混ぜる
m_s = GradientBoostingRegressor(**gbm).fit(np.column_stack([X, w]), y)
tau_s = (m_s.predict(np.column_stack([X_test, np.ones(2000)]))
- m_s.predict(np.column_stack([X_test, np.zeros(2000)])))
# T-learner:処置群と対照群で別々のモデル
m0 = GradientBoostingRegressor(**gbm).fit(X[w == 0], y[w == 0])
m1 = GradientBoostingRegressor(**gbm).fit(X[w == 1], y[w == 1])
tau_t = m1.predict(X_test) - m0.predict(X_test)
for name, t in [("S", tau_s), ("T", tau_t)]:
print("%s-learner: PEHE %.3f, 推定ATE %.2f, τ̂の標準偏差 %.2f"
% (name, pehe(t, tau_test), t.mean(), t.std()))
print("真のτ(x)の標準偏差 = %.2f" % tau_test.std())
S-learner: PEHE 0.329, 推定ATE 1.02, τ̂の標準偏差 0.62
T-learner: PEHE 0.677, 推定ATE 1.07, τ̂の標準偏差 0.98
真のτ(x)の標準偏差 = 0.71
各行の3つの数は、左から「効果の予測誤差(PEHE)」「$\hat\tau$ の平均、つまり推定 ATE」「$\hat\tau$ の散らばり(標準偏差)」です。散らばりは「人によって効果がどれだけ違うと推定したか」の目安で、最終行にある真の散らばり 0.71 と見比べます。両者とも推定 ATE は 1.0 付近に来ました。差が出るのは異質性の質です。T-learner の $\hat\tau$ の散らばり 0.98 は真の 0.71 を大きく超えています。異質性を捉えたのではなく、ノイズを異質性と誤認した膨らみです。原因は群サイズの不均衡で、$\hat\mu_1$ は939人だけで複雑なベースラインをまるごと学ぶ羽目になり、その推定誤差が引き算 $\hat\mu_1 - \hat\mu_0$ にそのまま乗ります。一方 S-learner の散らばり 0.62 は真より小さく、こちらは異質性が1割ほど潰れています。この「潰れ」は正則化を強めると劇的に進みます。勾配ブースティングで最も素朴な正則化は木の本数を減らすことなので、木の本数を絞ってみます。
for n_est in [200, 50, 20]:
m = GradientBoostingRegressor(**dict(gbm, n_estimators=n_est)).fit(
np.column_stack([X, w]), y)
t = (m.predict(np.column_stack([X_test, np.ones(2000)]))
- m.predict(np.column_stack([X_test, np.zeros(2000)])))
print("木 %3d 本: PEHE %.3f, 推定ATE %.2f, τ̂=0 の割合 %3.0f%%"
% (n_est, pehe(t, tau_test), t.mean(), 100 * np.mean(np.abs(t) < 0.01)))
木 200 本: PEHE 0.329, 推定ATE 1.02, τ̂=0 の割合 0%
木 50 本: PEHE 0.941, 推定ATE 0.35, τ̂=0 の割合 0%
木 20 本: PEHE 1.244, 推定ATE 0.00, τ̂=0 の割合 100%
木20本の S-learner は、全2,000点で $\hat\tau = 0$、推定 ATE 0.00 です。「処置効果は存在しない」という回答で、しかもこのモデルは $Y$ の予測器としては壊れていません(訓練データへの決定係数は 0.76 あり、ベースラインの大きな構造は20本でも粗く捉えています)。からくりはこうです。木の分岐は、予測の外れ(損失)を最も減らす特徴量から順に採用されますが、$W$ の分岐が効くのは6,000人中939人だけで、しかも効果の大きさは最大2。振幅10近いベースラインを説明する $x_2$ たちの分岐に比べて、利得が小さいのです。容量を絞られたモデルは小さな信号から切り捨てます。そして処置効果はたいてい小さな信号です。これが「正則化が $\tau$ を潰す」の中身で、$Y$ の精度だけ見ていると気づけません(概念②の伏線回収です)。
手を動かす②:X-learner を1ステップずつ
X-learner は、T-learner の2本を土台に「$\tau$ を直接学習する」段を重ねます。鍵は、観測できないはずの $\tau_i$ の代理ラベルを作ってしまうことです。
処置群の人 $i$ を考えます。$Y_i(1)$ は実測済みで、欠けているのは $Y_i(0)$ だけです。その欠けた片方を、対照群5,061人で鍛えた $\hat\mu_0$ の予測で埋めると、$D_i = Y_i - \hat\mu_0(X_i)$ という「半分実測・半分予測」の効果もどきが1人1個できます。対照群には逆向きに $D_j = \hat\mu_1(X_j) - Y_j$ を作ります。これでついに「$\tau$ を目的変数にした回帰」が書けます。ラベルは代理品ですが、処置群側の代理ラベルから $\hat\tau_1(x)$、対照群側から $\hat\tau_0(x)$ を学習します。代理ラベルの作り方を図にすると、処置群と対照群で「実測」と「予測」の役割が入れ替わっているのが見えます。
最後に2つをどう混ぜるか。ここで Day 7 の傾向スコアが再登場します。重み $g(x) = \Pr(W=1 \mid X=x)$(特徴量が $x$ の人が処置を受ける確率)を使い、$\hat\tau(x) = g(x)\hat\tau_0(x) + (1-g(x))\hat\tau_1(x)$ と混ぜます。直感はこうです。$\hat\tau_1$ の材料は「処置群の実測」と「対照群で鍛えた $\hat\mu_0$」でした。処置が少ない場所($g(x)$ が小さい場所)では対照群が豊富なので $\hat\mu_0$ が信頼でき、$\hat\tau_1$ を重く使うべきです。$g$ と $1-g$ の取り合わせは、まさにその配分になっています。今日の設定では処置が15%しかないので、ほとんどの場所で $g(x)$ は小さく、$\hat\tau_1$ の側が重く使われます。コードでは $g(x)$ を Day 7 と同じロジスティック回帰で推定し、C=np.inf は「正則化なし」の指定です。
from sklearn.linear_model import LogisticRegression
# 第1段:T-learner の2本 (m0, m1) をそのまま再利用する
# 第2段:観測できない τ_i の代理ラベルを作り、τ を直接学習する
d1 = y[w == 1] - m0.predict(X[w == 1]) # 処置群:実測Y − 予測Y(0)
d0 = m1.predict(X[w == 0]) - y[w == 0] # 対照群:予測Y(1) − 実測Y
cate = dict(n_estimators=300, max_depth=1, learning_rate=0.05, random_state=0)
tau1 = GradientBoostingRegressor(**cate).fit(X[w == 1], d1)
tau0 = GradientBoostingRegressor(**cate).fit(X[w == 0], d0)
# 第3段:傾向スコア g(x) で2つの推定を混ぜる(Day 7 の再登場)
g = LogisticRegression(C=np.inf).fit(X, w).predict_proba(X_test)[:, 1]
tau_x = g * tau0.predict(X_test) + (1 - g) * tau1.predict(X_test)
print("X-learner: PEHE %.3f, 推定ATE %.2f, τ̂の標準偏差 %.2f"
% (pehe(tau_x, tau_test), tau_x.mean(), tau_x.std()))
X-learner: PEHE 0.239, 推定ATE 1.12, τ̂の標準偏差 0.70
PEHE 0.239 で、S の 0.329 と T の 0.677 を抜いて最良です。$\hat\tau$ の散らばり 0.70 は真の 0.71 とほぼ一致し、潰しても膨らませてもいません(推定 ATE は 1.12 とやや上振れで、この点は S/T に劣ります。完璧な道具ではありません)。
1つ設計上の判断を仕込みました。第2段のモデルを max_depth=1 に絞ったことです(深さ1の木は分岐が1回だけの「切り株」で、1つの特徴量の単純な形しか表せません)。「CATE は単純」という見込みを、$\tau$ を学ぶモデルに正則化として直接注入しています。ここに X-learner の構造的な強みがあります。S/T では正則化がアウトカムモデルに掛かるため、強めると $\tau$ が巻き添えで歪みました。X は $\tau$ を学ぶ場所が独立にあるので、正則化を正しい場所に掛けられるのです。試しに第2段を第1段と同じ深いモデルにすると、代理ラベルに含まれるノイズまで学んでしまいます。
# 第2段に深いモデル(第1段と同じ設定)を使うと、代理ラベルのノイズまで学ぶ
tau1b = GradientBoostingRegressor(**gbm).fit(X[w == 1], d1)
tau0b = GradientBoostingRegressor(**gbm).fit(X[w == 0], d0)
tau_xb = g * tau0b.predict(X_test) + (1 - g) * tau1b.predict(X_test)
print("X-learner(第2段が深い): PEHE %.3f" % pehe(tau_xb, tau_test))
X-learner(第2段が深い): PEHE 0.468
第2段を深くしただけで PEHE は 0.239 から 0.468 に悪化しました。3つのラーナーの推定を真の $\tau(x)$ に重ねたのが下の図です。横軸は効果を決める共変量 $x_1$、縦軸は推定した効果 $\hat\tau$ で、赤い曲線が真の $\tau(x)$、点が評価用2,000点それぞれの推定値です。左の S は形を追えているものの上端がやや潰れ、中央の T は縦に散らばり、右の X が最も赤い真値曲線に沿っています。
手を動かす③:econml で同じものを1行ずつ対応付ける
手組みの中身が分かったので、実務で使うライブラリ実装と突き合わせます。econml は Microsoft が公開している因果推論の Python ライブラリで、そのメタラーナーは、今日書いた部品がそのまま引数名になっています。 SLearner の overall_model が S の1本、 TLearner の models が T の2本(X の第1段も同じ)、 XLearner の cate_models が第2段の $\tau$ モデル、 propensity_model が混合重み $g(x)$ です。
from econml.metalearners import SLearner, TLearner, XLearner
S = SLearner(overall_model=GradientBoostingRegressor(**gbm))
T = TLearner(models=GradientBoostingRegressor(**gbm))
Xl = XLearner(models=GradientBoostingRegressor(**gbm),
cate_models=GradientBoostingRegressor(**cate),
propensity_model=LogisticRegression(C=np.inf))
for name, est, hand in [("S", S, tau_s), ("T", T, tau_t), ("X", Xl, tau_x)]:
est.fit(y, w, X=X)
diff = np.abs(est.effect(X_test) - hand).max()
print("econml %s-learner: PEHE %.3f(手組みとの最大差 %.6f)"
% (name, pehe(est.effect(X_test), tau_test), diff))
econml S-learner: PEHE 0.329(手組みとの最大差 0.000000)
econml T-learner: PEHE 0.677(手組みとの最大差 0.000000)
econml X-learner: PEHE 0.239(手組みとの最大差 0.000000)
3つとも手組みとの最大差は小数第6位まで 0.000000、つまり econml は今日手組みした計算と同じものを実行しています。ライブラリはブラックボックスではなく、レシピの清書だと確認できました。今日の全結果を PEHE の1つの物差しに並べたのが下の図です。木20本の S-learner($\hat\tau$ が全点0)の 1.244 も同じ物差しに載せると、「正則化をどこに掛けるか」で結果がどれだけ動くかが一望できます。
つまづき・誤解しやすい点
- PEHE は実データでは計算できません。 真の $\tau$ が要るからです。今日のオラクル採点は「自分のパイプラインの癖を、答えを知っている環境で検査する」ための手段であって、実データの分析でモデルの良し悪しを PEHE で語ることは原理的に不可能です。$Y$ の交差検証で代用しても $\tau$ を採点したことにはなりません(木20本の S-learner がすり抜けます)。実データでの代理評価は Day 26 の Qini 曲線で、信頼区間の話は Day 25 で扱います。
- S-learner の「効果ゼロ」は正則化の産物かもしれません。 $\hat\tau$ が0に張り付いたら、効果が無いのではなく $W$ の分岐が容量制限で削られただけの可能性があります。正則化の強さを振って $\hat\tau$ の動きを見るのが今日のスイープの使い方です。裏を返すと、効果が本当に無い(または一様に小さい)場面では、ゼロに寄せる S-learner の癖はむしろ保守的で堅実に働きます。どのラーナーが勝つかは場面次第です。
- メタラーナーは交絡を解決しません。 土台は Day 6 の非交絡性のままです。割付を決める変数が $X$ に入っていなければ、S も T も X も同じようにバイアスされます(今日のデータで $x_2$ を特徴量から外すと三者そろってズレます)。機械学習を差し込んだからといって観察データが実験に化けるわけではない、という規律は Day 24 の Double ML でも繰り返します。
GISデータ実務での使い方
- 施策を受けた側はたいてい少数派です。 管路の更新や衛星データ解析の先行導入をした自治体は全体の1〜2割で、今日の処置15%の構図そのままです。この不均衡ではまず T-learner の弱点が出るので、CATE を組むなら X-learner から試すのが手筋になります。
- CATE は予算の優先順位に直結します。 「どの自治体・どの管路セグメントで施策の効果が大きいか」が分かれば、限られた更新予算を効果の高い順に配れます。ただし今日見たとおり CATE 推定は ATE より一段難しく、道具の癖も強いので、その期待値をチームで共有してから始めます。
- オラクルテストを習慣にします。 実データには真の $\tau$ がありません。だからこそ、業務データを模した合成データに真の $\tau(x)$ を仕込み、自分のパイプラインが PEHE で何点取れるかを先に検査します。今日のコードはそのままその雛形です。空間データでは評価の分割にも罠があり(近くの地点同士はよく似ているため、学習に使った情報が検証側へ漏れて見かけの性能が上がる「近傍リーク」)、空間CVとの合わせ技は Day 28 で掘ります。
参考(本棚)
- Matheus Facure『Causal Inference in Python』第7章 Metalearners — S/T/X の構成と、正則化が効果推定を歪める実演。今日の下敷き
- Künzel, Sekhon, Bickel & Yu (2019). Metalearners for estimating heterogeneous treatment effects using machine learning. PNAS 116(10) — X-learner の原論文。「群サイズが不均衡で CATE が単純な場面に強い」という設計意図の出どころ
- Robert Osazuwa Ness『Causal AI』と econml ドキュメント(Metalearners) — メタラーナーの実装視点と API。今日の照合相手
- データ: 合成データ(シード23)。真の $\tau(x)$ を握って PEHE をオラクル採点するための設計で、生成式は本文のコードに全部書いてあります
次回予告(Day 24)
今日は「誰に効くか」を機械学習で組みました。次回は問いをいったん平均に戻します。共変量が高次元でも、機械学習の当てはめ癖に推定を汚させずに ATE を厳密に出す Double Machine Learning(DML) です。今日 S-learner を潰した正則化バイアスが今度は主役として登場し、Day 16 で頭に焼き付けた FWL の残差化が機械学習版に進化します。残差×残差、交差適合、そして「ニュイサンス(推定の目的ではないが途中で必要になる量)を少し間違えても効かない」Neyman 直交性まで、いつもどおり手組みから積み上げます。聞き慣れない言葉が並びましたが、次回で1つずつほどきます。

